JuliaHomotopyContinuation / HomotopyContinuation.jl

A Julia package for solving systems of polynomials via homotopy continuation.
https://www.JuliaHomotopyContinuation.org
MIT License
178 stars 30 forks source link

solve errors for real target subspace #512

Closed jamblejoe closed 1 year ago

jamblejoe commented 1 year ago

Consider

using HomotopyContinuation, LinearAlgebra
@var x[1:3] y[1:3]
ψ = [x..., y...]

H = - (x[1]*x[2]+y[1]*y[2]) - (x[2]*x[3]+y[2]*y[3]) + (x[1]^2+y[1]^2)^2 + (x[2]^2+y[2]^2)^2 + (x[3]^2+y[3]^2)^2
S = x[1]^2 + x[2]^2 + x[3]^2 + y[1]^2 + y[2]^2 + y[3]^2 - 1
F = [H; S]
J = differentiate(F, ψ)

L₀ = rand_subspace(6; dim = 2)
start = solve(F, target_subspace = L₀)
start_sols = solutions(start)

solve(
    F,
    start_sols;
    start_subspace =  L₀,
    target_subspaces = [rand_subspace(6; dim = 2)],
)

All fine, then trying to target a real subspace

solve(
    F,
    start_sols;
    start_subspace =  L₀,
    target_subspaces = [rand_subspace(6; dim = 2, real=true)],
)

I get the following error

ERROR: MethodError: no method matching LRUCache.LinkedNode{Tuple{LinearSubspace{ComplexF64}, LinearSubspace{ComplexF64}}}(::Tuple{LinearSubspace{ComplexF64}, LinearSubspace{Float64}})
Closest candidates are:
  LRUCache.LinkedNode{T}(::T) where T at ~/.julia/packages/LRUCache/44dJX/src/cyclicorderedset.jl:9
Stacktrace:
  [1] _unsafe_addindex!(lru::LRUCache.LRU{Tuple{LinearSubspace{ComplexF64}, LinearSubspace{ComplexF64}}, HomotopyContinuation.LinearSubspaceGeodesicInfo}, v::HomotopyContinuation.LinearSubspaceGeodesicInfo, key::Tuple{LinearSubspace{ComplexF64}, LinearSubspace{Float64}})
    @ LRUCache ~/.julia/packages/LRUCache/44dJX/src/LRUCache.jl:136
  [2] get!(default::HomotopyContinuation.var"#100#101"{LinearSubspace{ComplexF64}, LinearSubspace{Float64}}, lru::LRUCache.LRU{Tuple{LinearSubspace{ComplexF64}, LinearSubspace{ComplexF64}}, HomotopyContinuation.LinearSubspaceGeodesicInfo}, key::Tuple{LinearSubspace{ComplexF64}, LinearSubspace{Float64}})
    @ LRUCache ~/.julia/packages/LRUCache/44dJX/src/LRUCache.jl:112
  [3] set_subspaces!(H::IntrinsicSubspaceHomotopy{MixedSystem{Int32, (0x84ec11465e833cef, 1)}}, start::LinearSubspace{ComplexF64}, target::LinearSubspace{Float64})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/homotopies/intrinsic_subspace_homotopy.jl:134
  [4] target_parameters!(H::IntrinsicSubspaceHomotopy{MixedSystem{Int32, (0x84ec11465e833cef, 1)}}, q::LinearSubspace{Float64})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/homotopies/intrinsic_subspace_homotopy.jl:141
  [5] target_parameters!(T::Tracker{IntrinsicSubspaceHomotopy{MixedSystem{Int32, (0x84ec11465e833cef, 1)}}, Matrix{ComplexF64}}, p::LinearSubspace{Float64})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/tracker.jl:1028
  [6] target_parameters!(T::EndgameTracker{IntrinsicSubspaceHomotopy{MixedSystem{Int32, (0x84ec11465e833cef, 1)}}, Matrix{ComplexF64}}, p::LinearSubspace{Float64})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/endgame_tracker.jl:826
  [7] target_parameters!
    @ ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:647 [inlined]
  [8] many_solve(solver::Solver{EndgameTracker{IntrinsicSubspaceHomotopy{MixedSystem{Int32, (0x84ec11465e833cef, 1)}}, Matrix{ComplexF64}}}, starts::Vector{Vector{ComplexF64}}, many_target_parameters::Vector{LinearSubspace{Float64}}, progress::ProgressMeter.Progress, transform_result::typeof(tuple), transform_parameters::typeof(identity), ::Val{false}; threading::Bool, catch_interrupt::Bool)
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:753
  [9] #solve#291
    @ ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:698 [inlined]
 [10] solve(::Vector{Expression}, ::Vararg{Any}; show_progress::Bool, threading::Bool, catch_interrupt::Bool, target_parameters::Nothing, stop_early_cb::Function, transform_result::Nothing, transform_parameters::typeof(identity), flatten::Nothing, target_subspaces::Vector{LinearSubspace{Float64}}, kwargs::Base.Pairs{Symbol, LinearSubspace{Float64}, Tuple{Symbol}, NamedTuple{(:start_subspace,), Tuple{LinearSubspace{Float64}}}})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:475
 [11] top-level scope
    @ REPL[16]:1
julia> versioninfo()
Julia Version 1.8.0
Commit 5544a0fab76 (2022-08-17 13:38 UTC)
Platform Info:
  OS: Linux (x86_64-linux-gnu)
  CPU: 48 × AMD Ryzen Threadripper 3960X 24-Core Processor
  WORD_SIZE: 64
  LIBM: libopenlibm
  LLVM: libLLVM-13.0.1 (ORCJIT, znver2)
  Threads: 1 on 48 virtual cores
[f213a82b] HomotopyContinuation v2.6.4
PBrdng commented 1 year ago

Hi, this seems to be a type error.

When I do this, it works:

L = convert(LinearSubspace{ComplexF64}, rand_subspace(6; dim = 2, real = true))

solve(
    F,
    start_sols;
    start_subspace =  L₀,
    target_subspaces = [L],
)

I'll look into the code to find the origin of this problem and fix it in the next release.

jamblejoe commented 1 year ago

Hi, this seems to be a type error.

When I do this, it works:

L = convert(LinearSubspace{ComplexF64}, rand_subspace(6; dim = 2, real = true))

solve(
    F,
    start_sols;
    start_subspace =  L₀,
    target_subspaces = [L],
)

I'll look into the code to find the origin of this problem and fix it in the next release.

It works for me as well.

jamblejoe commented 1 year ago

I realized something weird. After this

solve(
    F,
    start_sols;
    start_subspace =  L₀,
    target_subspaces = [rand_subspace(6; dim = 2, real=true)],
)

throwing the above error once, calling solve on F again never finishes. E.g. running the same code again and interrupting the process after a while I get

ERROR: InterruptException:
Stacktrace:
  [1] lock(l::Base.Threads.SpinLock)
    @ Base.Threads ./locks-mt.jl:46
  [2] get!(default::HomotopyContinuation.var"#96#98"{LinearSubspace{ComplexF64}, LinearSubspace{ComplexF64}}, lru::LRUCache.LRU{Tuple{LinearSubspace{ComplexF64}, LinearSubspace{ComplexF64}}, HomotopyContinuation.LinearSubspaceGeodesicInfo}, key::Tuple{LinearSubspace{ComplexF64}, LinearSubspace{ComplexF64}})
    @ LRUCache ~/.julia/packages/LRUCache/44dJX/src/LRUCache.jl:98
  [3] IntrinsicSubspaceHomotopy(system::MixedSystem{Int32, (0x84ec11465e833cef, 1)}, start::LinearSubspace{ComplexF64}, target::LinearSubspace{ComplexF64})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/homotopies/intrinsic_subspace_homotopy.jl:97
  [4] IntrinsicSubspaceHomotopy(system::MixedSystem{Int32, (0x84ec11465e833cef, 1)}, start::LinearSubspace{ComplexF64}, target::LinearSubspace{Float64}; compile::Symbol)
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/homotopies/intrinsic_subspace_homotopy.jl:85
  [5] IntrinsicSubspaceHomotopy(system::MixedSystem{Int32, (0x84ec11465e833cef, 1)}, start::LinearSubspace{ComplexF64}, target::LinearSubspace{Float64})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/homotopies/intrinsic_subspace_homotopy.jl:79
  [6] IntrinsicSubspaceHomotopy(F::System, start::LinearSubspace{ComplexF64}, target::LinearSubspace{Float64}; compile::Symbol)
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/homotopies/intrinsic_subspace_homotopy.jl:76
  [7] linear_subspace_homotopy(F::System, V::LinearSubspace{ComplexF64}, W::LinearSubspace{Float64}; compile::Symbol, intrinsic::Nothing)
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:277
  [8] solver_startsolutions(F::System, starts::Vector{Vector{ComplexF64}}; seed::UInt32, start_system::Symbol, generic_parameters::Nothing, p₁::Nothing, start_parameters::Nothing, p₀::Nothing, target_parameters::Nothing, compile::Symbol, start_subspace::LinearSubspace{ComplexF64}, target_subspace::LinearSubspace{Float64}, intrinsic::Nothing, kwargs::Base.Pairs{Symbol, Union{}, Tuple{}, NamedTuple{(), Tuple{}}})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:138
  [9] #solver_startsolutions#267
    @ ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:84 [inlined]
 [10] solve(::Vector{Expression}, ::Vararg{Any}; show_progress::Bool, threading::Bool, catch_interrupt::Bool, target_parameters::Nothing, stop_early_cb::Function, transform_result::Nothing, transform_parameters::typeof(identity), flatten::Nothing, target_subspaces::Vector{LinearSubspace{Float64}}, kwargs::Base.Pairs{Symbol, LinearSubspace{ComplexF64}, Tuple{Symbol}, NamedTuple{(:start_subspace,), Tuple{LinearSubspace{ComplexF64}}}})
    @ HomotopyContinuation ~/.julia/packages/HomotopyContinuation/I1faM/src/solve.jl:449
 [11] top-level scope
    @ REPL[19]:1
PBrdng commented 1 year ago

I also observed this, and I don't know yet the reason for this.

Thanks for finding all these bugs, btw 👍

jamblejoe commented 1 year ago

Sure! I like your work and would like to use this package for sampling points from manifolds defined by ${H(x) = constant}$ for Hamiltonians $H$. So thanks for putting together this package!

saschatimme commented 1 year ago

This looks like that the method error doesn't play well with the multi-threading. If I remember correctly (and this is a while) we only have proper error catching inside the path tracking but not in the "assemble the necessary pieces for tracking" part. I would expect that solving the type error in #513 should fix the problem.

PBrdng commented 1 year ago

The fix is now on the master.