Skip to content

Commit

Permalink
feat: add safe_HomotopyContinuationProblem
Browse files Browse the repository at this point in the history
  • Loading branch information
AayushSabharwal committed Nov 30, 2024
1 parent 0f1e80e commit fecca74
Show file tree
Hide file tree
Showing 3 changed files with 49 additions and 2 deletions.
10 changes: 8 additions & 2 deletions ext/MTKHomotopyContinuationExt.jl
Original file line number Diff line number Diff line change
Expand Up @@ -96,16 +96,22 @@ All other keyword arguments are forwarded to `HomotopyContinuation.solver_starts
"""
function MTK.HomotopyContinuationProblem(
sys::NonlinearSystem, u0map, parammap = nothing; kwargs...)
prob = MTK._safe_HomotopyContinuationProblem(sys, u0map, parammap; kwargs...)
prob isa MTK.HomotopyContinuationProblem || throw(prob)
return prob
end

function MTK._safe_HomotopyContinuationProblem(sys, u0map, parammap = nothing; kwargs...)
if !iscomplete(sys)
error("A completed `NonlinearSystem` is required. Call `complete` or `structural_simplify` on the system before creating a `HomotopyContinuationProblem`")
end
transformation = MTK.PolynomialTransformation(sys)
if transformation isa MTK.NotPolynomialError
throw(transformation)
return transformation
end
result = MTK.transform_system(sys, transformation)
if result isa MTK.NotPolynomialError
throw(result)
return result
end
MTK.HomotopyContinuationProblem(sys, transformation, result, u0map, parammap; kwargs...)
end
Expand Down
21 changes: 21 additions & 0 deletions src/systems/nonlinear/homotopy_continuation.jl
Original file line number Diff line number Diff line change
Expand Up @@ -48,6 +48,27 @@ function HomotopyContinuationProblem(::AbstractSystem, _u0, _p; kwargs...)
error("HomotopyContinuation.jl is required to create and solve `HomotopyContinuationProblem`s. Please run `Pkg.add(\"HomotopyContinuation\")` to continue.")
end

"""
$(TYPEDSIGNATURES)
Utility function for `safe_HomotopyContinuationProblem`, implemented in the extension.
"""
function _safe_HomotopyContinuationProblem end

"""
$(TYPEDSIGNATURES)
Return a `HomotopyContinuationProblem` if the extension is loaded and the system is
polynomial. If the extension is not loaded, return `nothing`. If the system is not
polynomial, return the appropriate `NotPolynomialError`.
"""
function safe_HomotopyContinuationProblem(sys::NonlinearSystem, args...; kwargs...)
if Base.get_extension(ModelingToolkit, :MTKHomotopyContinuationExt) === nothing
return nothing
end
return _safe_HomotopyContinuationProblem(sys, args...; kwargs...)
end

SymbolicIndexingInterface.symbolic_container(p::HomotopyContinuationProblem) = p.sys
SymbolicIndexingInterface.state_values(p::HomotopyContinuationProblem) = p.u0
function SymbolicIndexingInterface.set_state!(p::HomotopyContinuationProblem, args...)
Expand Down
20 changes: 20 additions & 0 deletions test/extensions/homotopy_continuation.jl
Original file line number Diff line number Diff line change
@@ -1,6 +1,19 @@
using ModelingToolkit, NonlinearSolve, SymbolicIndexingInterface
import ModelingToolkit as MTK
using LinearAlgebra
using Test

@testset "Safe HCProblem" begin
@variables x y z
eqs = [0 ~ x^2 + y^2 + 2x * y
0 ~ x^2 + 4x + 4
0 ~ y * z + 4x^2]
@mtkbuild sys = NonlinearSystem(eqs)
prob = MTK.safe_HomotopyContinuationProblem(sys, [x => 1.0, y => 1.0, z => 1.0], [])
@test prob === nothing
end


import HomotopyContinuation

@testset "No parameters" begin
Expand Down Expand Up @@ -78,30 +91,37 @@ end
@test_throws ["Cannot convert", "Unable", "symbolically solve",
"Exponent", "not an integer", "not a polynomial"] HomotopyContinuationProblem(
sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError
@mtkbuild sys = NonlinearSystem([x^x - x ~ 0])
@test_throws ["Cannot convert", "Unable", "symbolically solve",
"Exponent", "unknowns", "not a polynomial"] HomotopyContinuationProblem(
sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError
@mtkbuild sys = NonlinearSystem([((x^2) / sin(x))^2 + x ~ 0])
@test_throws ["Cannot convert", "both polynomial", "non-polynomial",
"recognized", "sin", "not a polynomial"] HomotopyContinuationProblem(
sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError

@variables y = 2.0
@mtkbuild sys = NonlinearSystem([x^2 + y^2 + 2 ~ 0, y ~ sin(x)])
@test_throws ["Cannot convert", "recognized", "sin", "not a polynomial"] HomotopyContinuationProblem(
sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError

@mtkbuild sys = NonlinearSystem([x^2 + y^2 - 2 ~ 0, sin(x + y) ~ 0])
@test_throws ["Cannot convert", "function of multiple unknowns"] HomotopyContinuationProblem(
sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError

@mtkbuild sys = NonlinearSystem([sin(x)^2 + 1 ~ 0, cos(y) - cos(x) - 1 ~ 0])
@test_throws ["Cannot convert", "multiple non-polynomial terms", "same unknown"] HomotopyContinuationProblem(
sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError

@mtkbuild sys = NonlinearSystem([sin(x^2)^2 + sin(x^2) - 1 ~ 0])
@test_throws ["import Nemo"] HomotopyContinuationProblem(sys, [])
@test MTK.safe_HomotopyContinuationProblem(sys, []) isa MTK.NotPolynomialError
end

import Nemo
Expand Down

0 comments on commit fecca74

Please sign in to comment.