refactor: format

SciML · Jan 27, 2025 · ac6c3dc · ac6c3dc
1 parent f66af4e
commit ac6c3dc
Show file tree

Hide file tree

Showing 5 changed files with 88 additions and 54 deletions.
diff --git a/lib/NonlinearSolveHomotopyContinuation/src/NonlinearSolveHomotopyContinuation.jl b/lib/NonlinearSolveHomotopyContinuation/src/NonlinearSolveHomotopyContinuation.jl
@@ -45,7 +45,8 @@ specified using ADTypes.jl.
 HomotopyContinuation.jl requires the taylor series of the polynomial system for the single
 root method. This is automatically computed using TaylorSeries.jl.
 """
-@concrete struct HomotopyContinuationJL{AllRoots} <: NonlinearSolveBase.AbstractNonlinearSolveAlgorithm
+@concrete struct HomotopyContinuationJL{AllRoots} <:
+                 NonlinearSolveBase.AbstractNonlinearSolveAlgorithm
     autodiff
     kwargs
 end

diff --git a/lib/NonlinearSolveHomotopyContinuation/src/interface_types.jl b/lib/NonlinearSolveHomotopyContinuation/src/interface_types.jl
@@ -18,23 +18,26 @@ regardless of the signature of ``f``.
     f
 end
 
-function (cjw::ComplexJacobianWrapper{Inplace})(u::AbstractVector{T}, x::AbstractVector{T}, p) where {T}
+function (cjw::ComplexJacobianWrapper{Inplace})(
+        u::AbstractVector{T}, x::AbstractVector{T}, p) where {T}
     x = reinterpret(Complex{T}, x)
     u = reinterpret(Complex{T}, u)
     cjw.f(u, x, p)
     u = parent(u)
     return u
 end
 
-function (cjw::ComplexJacobianWrapper{OutOfPlace})(u::AbstractVector{T}, x::AbstractVector{T}, p) where {T}
+function (cjw::ComplexJacobianWrapper{OutOfPlace})(
+        u::AbstractVector{T}, x::AbstractVector{T}, p) where {T}
     x = reinterpret(Complex{T}, x)
     u_tmp = cjw.f(x, p)
     u_tmp = reinterpret(T, u_tmp)
     copyto!(u, u_tmp)
     return u
 end
 
-function (cjw::ComplexJacobianWrapper{Scalar})(u::AbstractVector{T}, x::AbstractVector{T}, p) where {T}
+function (cjw::ComplexJacobianWrapper{Scalar})(
+        u::AbstractVector{T}, x::AbstractVector{T}, p) where {T}
     x = reinterpret(Complex{T}, x)
     u_tmp = cjw.f(x[1], p)
     u[1] = real(u_tmp)
@@ -52,7 +55,8 @@ polynomial systems specified using `NonlinearProblem`.
 
 $(FIELDS)
 """
-@concrete struct HomotopySystemWrapper{variant <: HomotopySystemVariant} <: HC.AbstractSystem
+@concrete struct HomotopySystemWrapper{variant <: HomotopySystemVariant} <:
+                 HC.AbstractSystem
     """
     The wrapped polynomial function.
     """
@@ -107,14 +111,16 @@ function HC.ModelKit.evaluate!(u, sys::HomotopySystemWrapper{Scalar}, x, p = not
     return u
 end
 
-function HC.ModelKit.evaluate_and_jacobian!(u, U, sys::HomotopySystemWrapper{Inplace}, x, p = nothing)
+function HC.ModelKit.evaluate_and_jacobian!(
+        u, U, sys::HomotopySystemWrapper{Inplace}, x, p = nothing)
     p = sys.p
     sys.f(u, x, p)
     sys.jac(U, x, p)
     return u, U
 end
 
-function HC.ModelKit.evaluate_and_jacobian!(u, U, sys::HomotopySystemWrapper{OutOfPlace}, x, p = nothing)
+function HC.ModelKit.evaluate_and_jacobian!(
+        u, U, sys::HomotopySystemWrapper{OutOfPlace}, x, p = nothing)
     p = sys.p
     u_tmp = sys.f(x, p)
     copyto!(u, u_tmp)
@@ -123,15 +129,18 @@ function HC.ModelKit.evaluate_and_jacobian!(u, U, sys::HomotopySystemWrapper{Out
     return u, U
 end
 
-function HC.ModelKit.evaluate_and_jacobian!(u, U, sys::HomotopySystemWrapper{Scalar}, x, p = nothing)
+function HC.ModelKit.evaluate_and_jacobian!(
+        u, U, sys::HomotopySystemWrapper{Scalar}, x, p = nothing)
     p = sys.p
     u[1] = sys.f(x[1], p)
     U[1] = sys.jac(x[1], p)
     return u, U
 end
 
 for V in (Inplace, OutOfPlace, Scalar)
-    @eval function HC.ModelKit.evaluate_and_jacobian!(u, U, sys::HomotopySystemWrapper{$V, F, J}, x, p = nothing) where {F, J <: ComplexJacobianWrapper}
+    @eval function HC.ModelKit.evaluate_and_jacobian!(
+            u, U, sys::HomotopySystemWrapper{$V, F, J}, x,
+            p = nothing) where {F, J <: ComplexJacobianWrapper}
         p = sys.p
         U_tmp = sys.jacobian_buffers
         x = reinterpret(Float64, x)
@@ -151,9 +160,10 @@ for V in (Inplace, OutOfPlace, Scalar)
     end
 end
 
-function update_taylorvars_from_taylorvector!(vars, x::HC.ModelKit.TaylorVector{M}) where {M}
+function update_taylorvars_from_taylorvector!(
+        vars, x::HC.ModelKit.TaylorVector{M}) where {M}
     for i in eachindex(vars)
-        for j in 0:M-1
+        for j in 0:(M - 1)
             vars[i][j] = x[i, j + 1]
         end
         for j in M:4
@@ -171,19 +181,22 @@ function update_taylorvars_from_taylorvector!(vars, x::AbstractVector)
     end
 end
 
-function update_maybe_taylorvector_from_taylorvars!(u::Vector, vars, buffer, ::Val{N}) where {N}
+function update_maybe_taylorvector_from_taylorvars!(
+        u::Vector, vars, buffer, ::Val{N}) where {N}
     for i in eachindex(vars)
         u[i] = buffer[i][N]
     end
 end
 
-function update_maybe_taylorvector_from_taylorvars!(u::HC.ModelKit.TaylorVector, vars, buffer, ::Val{N}) where {N}
+function update_maybe_taylorvector_from_taylorvars!(
+        u::HC.ModelKit.TaylorVector, vars, buffer, ::Val{N}) where {N}
     for i in eachindex(vars)
         u[i] = ntuple(j -> buffer[i][j - 1], Val(N + 1))
     end
 end
 
-function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N}, sys::HomotopySystemWrapper{Inplace}, x, p = nothing) where {N}
+function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N},
+        sys::HomotopySystemWrapper{Inplace}, x, p = nothing) where {N}
     f = sys.f
     p = sys.p
     buffer, vars = sys.taylorvars
@@ -193,7 +206,8 @@ function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N}, sys::HomotopySystemWra
     return u
 end
 
-function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N}, sys::HomotopySystemWrapper{OutOfPlace}, x, p = nothing) where {N}
+function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N},
+        sys::HomotopySystemWrapper{OutOfPlace}, x, p = nothing) where {N}
     f = sys.f
     p = sys.p
     vars = sys.taylorvars
@@ -203,7 +217,8 @@ function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N}, sys::HomotopySystemWra
     return u
 end
 
-function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N}, sys::HomotopySystemWrapper{Scalar}, x, p = nothing) where {N}
+function HC.ModelKit.taylor!(u::AbstractVector, ::Val{N},
+        sys::HomotopySystemWrapper{Scalar}, x, p = nothing) where {N}
     f = sys.f
     p = sys.p
     var = sys.taylorvars
@@ -270,8 +285,10 @@ function HC.ModelKit.evaluate_and_jacobian!(u, U, h::GuessHomotopy, x, t, p = no
     return u, U
 end
 
-HC.ModelKit.taylor!(u, v::Val{N}, H::GuessHomotopy, tx, t, incremental::Bool) where {N} =
+function HC.ModelKit.taylor!(
+        u, v::Val{N}, H::GuessHomotopy, tx, t, incremental::Bool) where {N}
     HC.ModelKit.taylor!(u, v, H, tx, t)
+end
 
 function HC.ModelKit.taylor!(u, ::Val{N}, h::GuessHomotopy, x, t, p = nothing) where {N}
     HC.ModelKit.taylor!(h.taylorbuffer, Val(N), h.sys, x, p)

diff --git a/lib/NonlinearSolveHomotopyContinuation/src/solve.jl b/lib/NonlinearSolveHomotopyContinuation/src/solve.jl
@@ -4,7 +4,8 @@
 Create and return the appropriate `HomotopySystemWrapper` to use for solving the given
 `prob` with `alg`.
 """
-function homotopy_continuation_preprocessing(prob::NonlinearProblem, alg::HomotopyContinuationJL)
+function homotopy_continuation_preprocessing(
+        prob::NonlinearProblem, alg::HomotopyContinuationJL)
     # cast to a `HomotopyNonlinearFunction`
     f = if prob.f isa HomotopyNonlinearFunction
         prob.f
@@ -45,7 +46,8 @@ function homotopy_continuation_preprocessing(prob::NonlinearProblem, alg::Homoto
     taylorvars = if isscalar
         Taylor1(zeros(ComplexF64, 5), 4)
     elseif iip
-        ([Taylor1(zeros(ComplexF64, 5), 4) for _ in u0], [Taylor1(zeros(ComplexF64, 5), 4) for _ in u0])
+        ([Taylor1(zeros(ComplexF64, 5), 4) for _ in u0],
+            [Taylor1(zeros(ComplexF64, 5), 4) for _ in u0])
     else
         [Taylor1(zeros(ComplexF64, 5), 4) for _ in u0]
     end
@@ -57,12 +59,14 @@ function homotopy_continuation_preprocessing(prob::NonlinearProblem, alg::Homoto
     end
 
     # HC-compatible system
-    hcsys = HomotopySystemWrapper{variant}(f.f.f, jac, p, alg.autodiff, prep, vars, taylorvars, jacobian_buffers)
+    hcsys = HomotopySystemWrapper{variant}(
+        f.f.f, jac, p, alg.autodiff, prep, vars, taylorvars, jacobian_buffers)
 
     return f, hcsys
 end
 
-function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{true}; denominator_abstol = 1e-7, kwargs...)
+function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{true};
+        denominator_abstol = 1e-7, kwargs...)
     f, hcsys = homotopy_continuation_preprocessing(prob, alg)
 
     u0 = state_values(prob)
@@ -110,7 +114,8 @@ function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{t
     return SciMLBase.EnsembleSolution(nlsols, 0.0, true, nothing)
 end
 
-function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{false}; denominator_abstol = 1e-7, kwargs...)
+function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{false};
+        denominator_abstol = 1e-7, kwargs...)
     f, hcsys = homotopy_continuation_preprocessing(prob, alg)
 
     u0 = state_values(prob)
@@ -120,14 +125,16 @@ function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{f
     fu0 = NonlinearSolveBase.Utils.evaluate_f(prob, u0_p)
 
     homotopy = GuessHomotopy(hcsys, fu0)
-    orig_sol = HC.solve(homotopy, u0_p isa Number ? [[u0_p]] : [u0_p]; alg.kwargs..., kwargs...)
+    orig_sol = HC.solve(
+        homotopy, u0_p isa Number ? [[u0_p]] : [u0_p]; alg.kwargs..., kwargs...)
     realsols = map(res -> res.solution, HC.results(orig_sol; only_real = true))
     if u0 isa Number
         realsols = map(only, realsols)
     end
 
     # no real solutions or infeasible solution
-    if isempty(realsols) || any(<=(denominator_abstol), map(abs, f.denominator(real.(only(realsols)), p)))
+    if isempty(realsols) ||
+       any(<=(denominator_abstol), map(abs, f.denominator(real.(only(realsols)), p)))
         retcode = if isempty(realsols)
             SciMLBase.ReturnCode.ConvergenceFailure
         else
@@ -153,4 +160,3 @@ function CommonSolve.solve(prob::NonlinearProblem, alg::HomotopyContinuationJL{f
     retcode = SciMLBase.ReturnCode.Success
     return SciMLBase.build_solution(prob, alg, u, resid; retcode, original = orig_sol)
 end
-
diff --git a/lib/NonlinearSolveHomotopyContinuation/test/allroots.jl b/lib/NonlinearSolveHomotopyContinuation/test/allroots.jl
@@ -11,7 +11,8 @@ alg = HomotopyContinuationJL{true}(; threading = false)
     jac = function (u, p)
         return 2u - p[1]
     end
-    @testset "`NonlinearProblem` - $name" for (jac, name) in [(nothing, "no jac"), (jac, "jac")]
+    @testset "`NonlinearProblem` - $name" for (jac, name) in [
+        (nothing, "no jac"), (jac, "jac")]
         fn = NonlinearFunction(rhs; jac)
         prob = NonlinearProblem(fn, 1.0, [5.0, 6.0])
         sol = solve(prob, alg)
@@ -21,14 +22,14 @@ alg = HomotopyContinuationJL{true}(; threading = false)
         sort!(sol.u; by = x -> x.u)
         @test sol.u[1] isa NonlinearSolution
         @test SciMLBase.successful_retcode(sol.u[1])
-        @test sol.u[1].u ≈ 2.0 atol = 1e-10
+        @test sol.u[1].u≈2.0 atol=1e-10
         @test sol.u[2] isa NonlinearSolution
         @test SciMLBase.successful_retcode(sol.u[2])
-        @test sol.u[2].u ≈ 3.0 atol = 1e-10
+        @test sol.u[2].u≈3.0 atol=1e-10
 
         @testset "no real solutions" begin
             prob = NonlinearProblem(1.0, 0.5) do u, p
-                return u * u - 2p * u  + p
+                return u * u - 2p * u + p
             end
             sol = solve(prob, alg)
             @test length(sol) == 1
@@ -47,7 +48,8 @@ alg = HomotopyContinuationJL{true}(; threading = false)
         unpolynomialize = function (u, p)
             return [asin(u)]
         end
-        fn = HomotopyNonlinearFunction(; denominator, polynomialize, unpolynomialize) do u, p
+        fn = HomotopyNonlinearFunction(;
+            denominator, polynomialize, unpolynomialize) do u, p
             return (u - p[1]) * (u - p[2])
         end
         prob = NonlinearProblem(fn, 0.0, [0.5, 0.7])
@@ -74,12 +76,12 @@ alg = HomotopyContinuationJL{true}(; threading = false)
 end
 
 f! = function (du, u, p)
-    du[1] = u[1] * u[1] - p[1] * u[2] + u[2] ^ 3 + 1
-    du[2] = u[2] ^ 3 + 2 * p[2] * u[1] * u[2] + u[2]
+    du[1] = u[1] * u[1] - p[1] * u[2] + u[2]^3 + 1
+    du[2] = u[2]^3 + 2 * p[2] * u[1] * u[2] + u[2]
 end
 
 f = function (u, p)
-    [u[1] * u[1] - p[1] * u[2] + u[2] ^ 3 + 1, u[2] ^ 3 + 2 * p[2] * u[1] * u[2] + u[2]]
+    [u[1] * u[1] - p[1] * u[2] + u[2]^3 + 1, u[2]^3 + 2 * p[2] * u[1] * u[2] + u[2]]
 end
 
 jac! = function (j, u, p)
@@ -89,11 +91,13 @@ jac! = function (j, u, p)
     j[2, 2] = 3 * u[2]^2 + 2 * p[2] * u[1] + 1
 end
 jac = function (u, p)
-    [2u[1] -p[1] + 3 * u[2]^2;
-     2 * p[2] * u[2] 3 * u[2]^2 + 2 * p[2] * u[1] + 1]
+    [2u[1] -p[1]+3 * u[2]^2;
+     2*p[2]*u[2] 3*u[2]^2+2*p[2]*u[1]+1]
 end
 
-@testset "vector u - $name" for (rhs, jac, name) in [(f, nothing, "oop"), (f, jac, "oop + jac"), (f!, nothing, "iip"), (f!, jac!, "iip + jac")]
+@testset "vector u - $name" for (rhs, jac, name) in [
+    (f, nothing, "oop"), (f, jac, "oop + jac"),
+    (f!, nothing, "iip"), (f!, jac!, "iip + jac")]
     sol = nothing
     @testset "`NonlinearProblem`" begin
         fn = NonlinearFunction(rhs; jac)
@@ -103,7 +107,7 @@ end
         @test sol.converged
         for nlsol in sol.u
             @test SciMLBase.successful_retcode(nlsol)
-            @test f(nlsol.u, prob.p) ≈ [0.0, 0.0] atol = 1e-10
+            @test f(nlsol.u, prob.p)≈[0.0, 0.0] atol=1e-10
         end
 
         @testset "no real solutions" begin
@@ -123,7 +127,7 @@ end
             return [[cbrt(u[1]), sin(u[2] / 40)]]
         end
         polynomialize = function (u, p)
-            return [u[1] ^ 3, 40asin(u[2])]
+            return [u[1]^3, 40asin(u[2])]
         end
         nlfn = NonlinearFunction(rhs; jac)
         fn = HomotopyNonlinearFunction(nlfn; denominator, polynomialize, unpolynomialize)
@@ -133,7 +137,9 @@ end
         @test sol2.converged
         @test length(sol.u) == length(sol2.u)
         for nlsol2 in sol2.u
-            @test any(nlsol -> isapprox(polynomialize(nlsol2.u, prob.p), nlsol.u; rtol = 1e-8), sol.u)
+            @test any(
+                nlsol -> isapprox(polynomialize(nlsol2.u, prob.p), nlsol.u; rtol = 1e-8),
+                sol.u)
         end
 
         @testset "some invalid solutions" begin