Julia: Unreachable reached in Atoms.jl yet again

[jagot]

This works on Julia v1.1 but crashes on 20834c3176. Will try bisect.

using Atoms
using SCF
using AtomicLevels
using AtomicPotentials
using FEDVRQuasi

rₘₐₓ, ρ = 10, 0.1
fedvr_order = 10
N = max(ceil(Int, rₘₐₓ/(ρ*fedvr_order)),2)
t = range(0.0, stop=rₘₐₓ, length=N)
R = FEDVR(t, fedvr_order)[:,2:end-1]

nucleus = pc"He"
atom = Atom(R, [spin_configurations(ground_state(nucleus))[1]],
            nucleus, eltype(R))

fock = Fock(atom)

Expected output:

Fock operator with
- quantum system: Atom{Float64}(R=FEDVR{Float64} basis with 9 elements on 0.0..10.0, restricted to basis functions 2..81 ⊂ 1..82; Z = 2 [helium]; 2 e⁻ ⇒ Q = 0) with 1 Configuration{SpinOrbital}: 1s²
- SCF equations:  
  - Hartree–Fock equation: E|1s₀α⟩ = OrbitalEquation(1s₀α): 
  [1, 1]  =  + 1ĥ₀|1s₀α⟩ + 1r⁻¹×Y⁰(1s₀β,1s₀β)|1s₀α⟩

⟨1s₀α| 𝓗 |1s₀α⟩ = -0.7500000012522258 Ha = -20.408250034074314 eV

  - Hartree–Fock equation: E|1s₀β⟩ = OrbitalEquation(1s₀β): 
  [1, 1]  =  + 1ĥ₀|1s₀β⟩ + 1r⁻¹×Y⁰(1s₀α,1s₀α)|1s₀β⟩

⟨1s₀β| 𝓗 |1s₀β⟩ = -0.7500000012522258 Ha = -20.408250034074314 eV

Actual output.txt

Environment:
Project.toml.txt
Manifest.toml.txt

[JeffBezanson]

Looks like the initial cause is some sort of memory corruption:

Internal error: encountered unexpected error in runtime:
TypeError(func=:TypeVar, context="upper bound", expected=Type, got=<?#0x7fd5fbe60150::<?#0x7fd5fbe60160::<?#0x7fd5fbea72d0::<?#0x7fd5fbea7330::(nil)>>>>)
rec_backtrace at /home/jagot/prog/julia/julia/src/stackwalk.c:94
jl_throw at /home/jagot/prog/julia/julia/src/task.c:210
jl_type_error_rt at /home/jagot/prog/julia/julia/src/rtutils.c:118
jl_new_typevar at /home/jagot/prog/julia/julia/src/builtins.c:996
finish_unionall at /home/jagot/prog/julia/julia/src/subtype.c:1771
intersect_unionall_ at /home/jagot/prog/julia/julia/src/subtype.c:1934

[JeffBezanson]

Progress in reducing this:

F = LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}},Diagonal{Float64,Array{Float64,1}},ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{Adjoint{Int64,BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}},ContinuumArrays.QuasiArrays.QuasiAdjoint{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}}}}}}}}

struct H{T,RO<:(LazyArrays.Applied{<:Any,typeof(*),s1495} where s1495<:(Tuple{s1494,Diagonal{T,Array{T,1}},<:Any} where s1494<:(Union{B, s1496} where s1496<:(ContinuumArrays.QuasiArrays.ApplyQuasiArray{<:Any,<:Any,s1496} where s1496<:(LazyArrays.Applied{<:Any,typeof(*),s1495} where s1495<:(Tuple{s1493,<:BandedMatrix})) where s1493<:B))) where B)}
end

typeintersect(H, H{<:Any,F}) # nonempty
typeintersect(H{<:Any,F}, H) # empty

[jagot]

Thanks @JeffBezanson.

CC @dlfivefifty

[JeffBezanson]

Reduced further:

julia> X = Pair{<:Any, Ref{Tuple{Ref{Ref{Tuple{Int}}},Ref{Float64}}}};

julia> Y = Pair{T, S} where S<:(Ref{A} where A<:(Tuple{C,Ref{T}} where C<:(Ref{D} where D<:(Ref{E} where E<:Tuple{F}) where F<:B)) where B) where T;

julia> typeintersect(X,Y)
Union{}

julia> typeintersect(Y,X)
Pair{Float64,Ref{Tuple{C,Ref{Float64}} where C<:Ref{Ref{Tuple{Int64}}}}}

[ViralBShah]

Best issue title.

[JeffBezanson]

Julia: we reach the unreachable.

[JeffBezanson]

I managed to fix my reduced case, but it doesn't fix the package. New reduction:

T = Atoms.HFPotential{_A,SpinOrbital{Orbital{Int64}},Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}},SubArray{Float64,1,Array{Float64,2},Tuple{Base.Slice{Base.OneTo{Int64}},Int64},true}}},LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}},Diagonal{Float64,Array{Float64,1}},ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{Adjoint{Int64,BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}},ContinuumArrays.QuasiArrays.QuasiAdjoint{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}}}}}}}},_B} where _B where _A

typeintersect(Atoms.HFPotential, T)  # nonempty
typeintersect(T, Atoms.HFPotential)  # empty

As a workaround, I suggest reducing the complexity of type parameter constraints. For example, changing HFPotential to this fixes it:

mutable struct HFPotential{kind,O,T,B<:Basis,OV,
                           RO<:HFPotentialOperator,P<:PoissonProblem} <: OrbitalIntegral{1,O,T,B,OV}

(There I just removed the "triangular" part of the constraint RO<:HFPotentialOperator{T,B}.)

[jagot]

Thanks @JeffBezanson! I should in general try to simplify my types, hopefully I could lower my compile times as well.

I finally finished the bisect:

57a3b087a5a9e0146ea0b516e8f83343f5f762b8 is the first bad commit
commit 57a3b087a5a9e0146ea0b516e8f83343f5f762b8
Author: Jeff Bezanson <[email protected]>
Date:   Fri Mar 29 16:25:33 2019 -0400

    fix #31439, type intersection issue with lower bounds (#31512)

:040000 040000 c724b1d15f9c42029f381add39f4516089d4f0a3 86a69e7a96d8869a7a15f2f2b2ee87cc0a4d99ff M  src
:040000 040000 b5d482333dd74cecc9004890f42bcb7da8aa3a21 bc037112f46419ab42b71bd079ef3231fd08c181 M  test

[JeffBezanson]

Ok! I may have a fix for the type intersection issue. Unfortunately, it seems to bring along a new method ambiguity. (What happens is, type intersection incorrectly returning Union{} causes these crashes, but also removes ambiguities since we think two signatures don't overlap when in fact they do. That leads to an arbitrary method being called, which usually works.) This is what I get:

MethodError: BlockBandedMatrices.BlockSkylineMatrix{Float64,DATA,BS} where BS<:BlockArrays.AbstractBlockSizes{2} where DATA<:AbstractArray{Float64,1}(::BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Array{Int64,1},Array{Int64,1},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}, ::BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}) is ambiguous. Candidates:
  (::Type{BlockBandedMatrices.BlockSkylineMatrix{T,Array{T,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}})(A::AbstractArray{T,2} where T, block_sizes::BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}) where T in BlockBandedMatrices at /home/jeff/.julia/packages/BlockBandedMatrices/meBlM/src/BlockSkylineMatrix.jl:176
  (::Type{BlockBandedMatrices.BlockSkylineMatrix{T,DATA,BS} where BS<:BlockArrays.AbstractBlockSizes{2} where DATA<:AbstractArray{T,1}})(A::BlockBandedMatrices.AbstractBlockBandedMatrix, block_sizes::BlockBandedMatrices.BlockSkylineSizes) where T in BlockBandedMatrices at /home/jeff/.julia/packages/BlockBandedMatrices/meBlM/src/BlockSkylineMatrix.jl:168
Possible fix, define
  (::Type{BlockBandedMatrices.BlockSkylineMatrix{T,Array{T,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}})(::BlockBandedMatrices.AbstractBlockBandedMatrix, ::BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}})
Stacktrace:
 [1] BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},Fill{Int64,1,Tuple{Base.OneTo{Int64}}},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}(::BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Array{Int64,1},Array{Int64,1},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}, ::Tuple{Array{Int64,1},Array{Int64,1}}, ::Tuple{Int64,Int64}) at /home/jeff/.julia/packages/BlockBandedMatrices/meBlM/src/BlockSkylineMatrix.jl:217
 [2] Type at /home/jeff/.julia/packages/BlockBandedMatrices/meBlM/src/BlockSkylineMatrix.jl:225 [inlined]
 [3] to_block_banded_matrix(::BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Array{Int64,1},Array{Int64,1},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}) at /home/jeff/.julia/packages/IterativeFactorizations/negfr/src/IterativeFactorizations.jl:135
 [4] factorization(::BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Array{Int64,1},Array{Int64,1},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}}) at /home/jeff/.julia/packages/IterativeFactorizations/negfr/src/IterativeFactorizations.jl:140
 [5] #ShiftInvert#10 at /home/jeff/.julia/packages/IterativeFactorizations/negfr/src/IterativeFactorizations.jl:141 [inlined]
 [6] Type at /home/jeff/.julia/packages/SCF/ajthS/src/shift_invert.jl:34 [inlined]
 [7] #diagonalize_one_body#21(::Symbol, ::Float64, ::Float64, ::Int64, ::Base.GenericIOBuffer{Array{UInt8,1}}, ::typeof(Atoms.diagonalize_one_body), ::LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}},BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Array{Int64,1},Array{Int64,1},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}},ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{Adjoint{Int64,BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}},ContinuumArrays.QuasiArrays.QuasiAdjoint{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}}}}}}}}, ::Int64) at /home/jeff/.julia/packages/Atoms/WYM3H/src/one_body.jl:157
 [8] (::getfield(Atoms, Symbol("#kw##diagonalize_one_body")))(::NamedTuple{(:σ, :io, :verbosity),Tuple{Float64,Base.GenericIOBuffer{Array{UInt8,1}},Int64}}, ::typeof(Atoms.diagonalize_one_body), ::LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}},BlockBandedMatrices.BlockSkylineMatrix{Float64,Array{Float64,1},BlockBandedMatrices.BlockSkylineSizes{BlockArrays.BlockSizes{2,Array{Int64,1}},Array{Int64,1},Array{Int64,1},BandedMatrix{Int64,Array{Int64,2},Base.OneTo{Int64}},Array{Int64,1}}},ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.QuasiArrays.LazyQuasiArrayApplyStyle,typeof(*),Tuple{Adjoint{Int64,BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}},ContinuumArrays.QuasiArrays.QuasiAdjoint{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}}}}}}}}, ::Int64) at ./none:0
 [9] (::getfield(Atoms, Symbol("##25#31")){Array{Float64,2},Dict{Int64,Array{SpinOrbital,1}},getfield(Atoms, Symbol("##24#29")){Int64,Int64},Formatting.FormatSpec,Formatting.FormatSpec,getfield(Atoms, Symbol("#fmt_eng_vec#30")),Int64,Int64,Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},Atom{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},SpinOrbital,Configuration{SpinOrbital},Float64,PointCharge{Int64}}})(::Base.GenericIOBuffer{Array{UInt8,1}}) at /home/jeff/.julia/packages/Atoms/WYM3H/src/hydrogenic.jl:67
 [10] print_block(::getfield(Atoms, Symbol("##25#31")){Array{Float64,2},Dict{Int64,Array{SpinOrbital,1}},getfield(Atoms, Symbol("##24#29")){Int64,Int64},Formatting.FormatSpec,Formatting.FormatSpec,getfield(Atoms, Symbol("#fmt_eng_vec#30")),Int64,Int64,Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},Atom{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},SpinOrbital,Configuration{SpinOrbital},Float64,PointCharge{Int64}}}, ::Base.GenericIOBuffer{Array{UInt8,1}}) at /home/jeff/.julia/packages/SCF/ajthS/src/utils.jl:3
 [11] (::getfield(Atoms, Symbol("##23#28")){SpinOrbital,Float64,Int64,Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},Atom{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},SpinOrbital,Configuration{SpinOrbital},Float64,PointCharge{Int64}}})(::Base.GenericIOBuffer{Array{UInt8,1}}) at /home/jeff/.julia/packages/Atoms/WYM3H/src/hydrogenic.jl:47
 [12] print_block(::getfield(Atoms, Symbol("##23#28")){SpinOrbital,Float64,Int64,Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}},Atom{Float64,FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},SpinOrbital,Configuration{SpinOrbital},Float64,PointCharge{Int64}}}, ::Base.TTY) at /home/jeff/.julia/packages/SCF/ajthS/src/utils.jl:3
 [13] print_block(::Function) at /home/jeff/.julia/packages/SCF/ajthS/src/utils.jl:2
 [14] #hydrogenic!#22 at /home/jeff/.julia/packages/Atoms/WYM3H/src/hydrogenic.jl:17 [inlined]
 [15] hydrogenic! at /home/jeff/.julia/packages/Atoms/WYM3H/src/hydrogenic.jl:14 [inlined]
 [16] #Atom#11(::Base.Iterators.Pairs{Union{},Union{},Tuple{},NamedTuple{(),Tuple{}}}, ::Type{Atom}, ::Symbol, ::Type{Float64}, ::ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}}, ::Array{Configuration{SpinOrbital},1}, ::PointCharge{Int64}, ::Type{Float64}) at /home/jeff/.julia/packages/Atoms/WYM3H/src/atom_types.jl:108
 [17] Type at /home/jeff/.julia/packages/Atoms/WYM3H/src/atom_types.jl:106 [inlined]
 [18] #Atom#12 at /home/jeff/.julia/packages/Atoms/WYM3H/src/atom_types.jl:134 [inlined]
 [19] Type at /home/jeff/.julia/packages/Atoms/WYM3H/src/atom_types.jl:134 [inlined]
 [20] #Atom#14 at /home/jeff/.julia/packages/Atoms/WYM3H/src/atom_types.jl:149 [inlined]
 [21] Atom(::ContinuumArrays.QuasiArrays.ApplyQuasiArray{Float64,2,LazyArrays.Applied{ContinuumArrays.BasisStyle,typeof(*),Tuple{FEDVR{Float64,Float64,Fill{Int64,1,Tuple{Base.OneTo{Int64}}}},BandedMatrix{Int64,Ones{Int64,2,Tuple{Base.OneTo{Int64},Base.OneTo{Int64}}},Base.OneTo{Int64}}}}}, ::Array{Configuration{SpinOrbital},1}, ::PointCharge{Int64}, ::Type{Float64}) at /home/jeff/.julia/packages/Atoms/WYM3H/src/atom_types.jl:149
 [22] top-level scope at /home/jeff/src/julia/i31703/i31703.jl:14

[jagot]

Interesting! I will have a look tomorrow if I am calling factorization incorrectly/ambiguously.

Again CC @dlfivefifty who's written most of BlockBandedMatrices.

[JeffBezanson]

And by sheer luck, I think I now have a version that fixes it without adding the ambiguity. It would still be worth thinking about whether the methods are indeed ambiguous though; there might be another potential problem lurking there.