diff --git a/lib/SimpleNonlinearSolve/.github/workflows/FormatPR.yml b/lib/SimpleNonlinearSolve/.github/workflows/FormatPR.yml new file mode 100644 index 000000000..b0cb28af5 --- /dev/null +++ b/lib/SimpleNonlinearSolve/.github/workflows/FormatPR.yml @@ -0,0 +1,29 @@ +name: format-pr +on: + schedule: + - cron: '0 0 * * *' +jobs: + build: + runs-on: ubuntu-latest + steps: + - uses: actions/checkout@v3 + - name: Install JuliaFormatter and format + run: | + julia -e 'using Pkg; Pkg.add(PackageSpec(name="JuliaFormatter"))' + julia -e 'using JuliaFormatter; format(".")' + # https://github.com/marketplace/actions/create-pull-request + # https://github.com/peter-evans/create-pull-request#reference-example + - name: Create Pull Request + id: cpr + uses: peter-evans/create-pull-request@v5 + with: + token: ${{ secrets.GITHUB_TOKEN }} + commit-message: Format .jl files + title: 'Automatic JuliaFormatter.jl run' + branch: auto-juliaformatter-pr + delete-branch: true + labels: formatting, automated pr, no changelog + - name: Check outputs + run: | + echo "Pull Request Number - ${{ steps.cpr.outputs.pull-request-number }}" + echo "Pull Request URL - ${{ steps.cpr.outputs.pull-request-url }}" \ No newline at end of file diff --git a/lib/SimpleNonlinearSolve/ext/SimpleBatchedNonlinearSolveExt.jl b/lib/SimpleNonlinearSolve/ext/SimpleBatchedNonlinearSolveExt.jl index e3f394c3b..dd7628888 100644 --- a/lib/SimpleNonlinearSolve/ext/SimpleBatchedNonlinearSolveExt.jl +++ b/lib/SimpleNonlinearSolve/ext/SimpleBatchedNonlinearSolveExt.jl @@ -31,7 +31,7 @@ function _init_J_batched(x::AbstractMatrix{T}) where {T} end function SciMLBase.__solve(prob::NonlinearProblem, alg::Broyden{true}, args...; - abstol = nothing, reltol = nothing, maxiters = 1000, kwargs...) + abstol = nothing, reltol = nothing, maxiters = 1000, kwargs...) tc = alg.termination_condition mode = DiffEqBase.get_termination_mode(tc) f = Base.Fix2(prob.f, prob.p) @@ -74,7 +74,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::Broyden{true}, args...; J⁻¹Δfₙ = _batched_mul(J⁻¹, Δfₙ) J⁻¹ += _batched_mul(((Δxₙ .- J⁻¹Δfₙ) ./ (_batched_mul(_batch_transpose(Δxₙ), J⁻¹Δfₙ) .+ T(1e-5))), - _batched_mul(_batch_transpose(Δxₙ), J⁻¹)) + _batched_mul(_batch_transpose(Δxₙ), J⁻¹)) if termination_condition(fₙ, xₙ, xₙ₋₁, atol, rtol) return SciMLBase.build_solution(prob, alg, xₙ, fₙ; retcode = ReturnCode.Success) diff --git a/lib/SimpleNonlinearSolve/src/SimpleNonlinearSolve.jl b/lib/SimpleNonlinearSolve/src/SimpleNonlinearSolve.jl index 46f1ab3c4..8749aa72b 100644 --- a/lib/SimpleNonlinearSolve/src/SimpleNonlinearSolve.jl +++ b/lib/SimpleNonlinearSolve/src/SimpleNonlinearSolve.jl @@ -16,7 +16,9 @@ end function __init__() @static if !isdefined(Base, :get_extension) - @require NNlib="872c559c-99b0-510c-b3b7-b6c96a88d5cd" begin include("../ext/SimpleBatchedNonlinearSolveExt.jl") end + @require NNlib="872c559c-99b0-510c-b3b7-b6c96a88d5cd" begin + include("../ext/SimpleBatchedNonlinearSolveExt.jl") + end end end @@ -42,31 +44,35 @@ include("alefeld.jl") import PrecompileTools -PrecompileTools.@compile_workload begin for T in (Float32, Float64) - prob_no_brack = NonlinearProblem{false}((u, p) -> u .* u .- p, T(0.1), T(2)) - for alg in (SimpleNewtonRaphson, Halley, Broyden, Klement, SimpleTrustRegion, - SimpleDFSane) - solve(prob_no_brack, alg(), abstol = T(1e-2)) - end +PrecompileTools.@compile_workload begin + for T in (Float32, Float64) + prob_no_brack = NonlinearProblem{false}((u, p) -> u .* u .- p, T(0.1), T(2)) + for alg in (SimpleNewtonRaphson, Halley, Broyden, Klement, SimpleTrustRegion, + SimpleDFSane) + solve(prob_no_brack, alg(), abstol = T(1e-2)) + end - #= - for alg in (SimpleNewtonRaphson,) - for u0 in ([1., 1.], StaticArraysCore.SA[1.0, 1.0]) - u0 = T.(.1) - probN = NonlinearProblem{false}((u,p) -> u .* u .- p, u0, T(2)) - solve(probN, alg(), tol = T(1e-2)) + #= + for alg in (SimpleNewtonRaphson,) + for u0 in ([1., 1.], StaticArraysCore.SA[1.0, 1.0]) + u0 = T.(.1) + probN = NonlinearProblem{false}((u,p) -> u .* u .- p, u0, T(2)) + solve(probN, alg(), tol = T(1e-2)) + end end - end - =# + =# - prob_brack = IntervalNonlinearProblem{false}((u, p) -> u * u - p, T.((0.0, 2.0)), T(2)) - for alg in (Bisection, Falsi, Ridder, Brent, Alefeld) - solve(prob_brack, alg(), abstol = T(1e-2)) + prob_brack = IntervalNonlinearProblem{false}((u, p) -> u * u - p, + T.((0.0, 2.0)), + T(2)) + for alg in (Bisection, Falsi, Ridder, Brent, Alefeld) + solve(prob_brack, alg(), abstol = T(1e-2)) + end end -end end +end # DiffEq styled algorithms export Bisection, Brent, Broyden, LBroyden, SimpleDFSane, Falsi, Halley, Klement, - Ridder, SimpleNewtonRaphson, SimpleTrustRegion, Alefeld + Ridder, SimpleNewtonRaphson, SimpleTrustRegion, Alefeld end # module diff --git a/lib/SimpleNonlinearSolve/src/ad.jl b/lib/SimpleNonlinearSolve/src/ad.jl index dbd90243a..009cd227b 100644 --- a/lib/SimpleNonlinearSolve/src/ad.jl +++ b/lib/SimpleNonlinearSolve/src/ad.jl @@ -29,50 +29,50 @@ function scalar_nlsolve_ad(prob, alg, args...; kwargs...) end function SciMLBase.solve(prob::NonlinearProblem{<:Union{Number, StaticArraysCore.SVector}, - iip, - <:Dual{T, V, P}}, - alg::AbstractSimpleNonlinearSolveAlgorithm, - args...; kwargs...) where {iip, T, V, P} + iip, + <:Dual{T, V, P}}, + alg::AbstractSimpleNonlinearSolveAlgorithm, + args...; kwargs...) where {iip, T, V, P} sol, partials = scalar_nlsolve_ad(prob, alg, args...; kwargs...) return SciMLBase.build_solution(prob, alg, Dual{T, V, P}(sol.u, partials), sol.resid; - retcode = sol.retcode) + retcode = sol.retcode) end function SciMLBase.solve(prob::NonlinearProblem{<:Union{Number, StaticArraysCore.SVector}, - iip, - <:AbstractArray{<:Dual{T, V, P}}}, - alg::AbstractSimpleNonlinearSolveAlgorithm, args...; - kwargs...) where {iip, T, V, P} + iip, + <:AbstractArray{<:Dual{T, V, P}}}, + alg::AbstractSimpleNonlinearSolveAlgorithm, args...; + kwargs...) where {iip, T, V, P} sol, partials = scalar_nlsolve_ad(prob, alg, args...; kwargs...) return SciMLBase.build_solution(prob, alg, Dual{T, V, P}(sol.u, partials), sol.resid; - retcode = sol.retcode) + retcode = sol.retcode) end # avoid ambiguities for Alg in [Bisection] @eval function SciMLBase.solve(prob::IntervalNonlinearProblem{uType, iip, - <:Dual{T, V, P}}, - alg::$Alg, args...; - kwargs...) where {uType, iip, T, V, P} + <:Dual{T, V, P}}, + alg::$Alg, args...; + kwargs...) where {uType, iip, T, V, P} sol, partials = scalar_nlsolve_ad(prob, alg, args...; kwargs...) return SciMLBase.build_solution(prob, alg, Dual{T, V, P}(sol.u, partials), - sol.resid; retcode = sol.retcode, - left = Dual{T, V, P}(sol.left, partials), - right = Dual{T, V, P}(sol.right, partials)) + sol.resid; retcode = sol.retcode, + left = Dual{T, V, P}(sol.left, partials), + right = Dual{T, V, P}(sol.right, partials)) #return BracketingSolution(Dual{T,V,P}(sol.left, partials), Dual{T,V,P}(sol.right, partials), sol.retcode, sol.resid) end @eval function SciMLBase.solve(prob::IntervalNonlinearProblem{uType, iip, - <:AbstractArray{ - <:Dual{T, - V, - P} - }}, - alg::$Alg, args...; - kwargs...) where {uType, iip, T, V, P} + <:AbstractArray{ + <:Dual{T, + V, + P}, + }}, + alg::$Alg, args...; + kwargs...) where {uType, iip, T, V, P} sol, partials = scalar_nlsolve_ad(prob, alg, args...; kwargs...) return SciMLBase.build_solution(prob, alg, Dual{T, V, P}(sol.u, partials), - sol.resid; retcode = sol.retcode, - left = Dual{T, V, P}(sol.left, partials), - right = Dual{T, V, P}(sol.right, partials)) + sol.resid; retcode = sol.retcode, + left = Dual{T, V, P}(sol.left, partials), + right = Dual{T, V, P}(sol.right, partials)) #return BracketingSolution(Dual{T,V,P}(sol.left, partials), Dual{T,V,P}(sol.right, partials), sol.retcode, sol.resid) end end diff --git a/lib/SimpleNonlinearSolve/src/alefeld.jl b/lib/SimpleNonlinearSolve/src/alefeld.jl index 61861a4a7..a2669ca96 100644 --- a/lib/SimpleNonlinearSolve/src/alefeld.jl +++ b/lib/SimpleNonlinearSolve/src/alefeld.jl @@ -9,9 +9,9 @@ algorithm 4.1 because, in certain sense, the second algorithm(4.2) is an optimal struct Alefeld <: AbstractBracketingAlgorithm end function SciMLBase.solve(prob::IntervalNonlinearProblem, - alg::Alefeld, args...; abstol = nothing, - reltol = nothing, - maxiters = 1000, kwargs...) + alg::Alefeld, args...; abstol = nothing, + reltol = nothing, + maxiters = 1000, kwargs...) f = Base.Fix2(prob.f, prob.p) a, b = prob.tspan c = a - (b - a) / (f(b) - f(a)) * f(a) @@ -19,14 +19,14 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, fc = f(c) (a == c || b == c) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.FloatingPointLimit, - left = a, - right = b) + retcode = ReturnCode.FloatingPointLimit, + left = a, + right = b) iszero(fc) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.Success, - left = a, - right = b) + retcode = ReturnCode.Success, + left = a, + right = b) a, b, d = _bracket(f, a, b, c) e = zero(a) # Set e as 0 before iteration to avoid a non-value f(e) @@ -45,14 +45,14 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, ē, fc = d, f(c) (a == c || b == c) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.FloatingPointLimit, - left = a, - right = b) + retcode = ReturnCode.FloatingPointLimit, + left = a, + right = b) iszero(fc) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.Success, - left = a, - right = b) + retcode = ReturnCode.Success, + left = a, + right = b) ā, b̄, d̄ = _bracket(f, a, b, c) # The second bracketing block @@ -68,14 +68,14 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, fc = f(c) (ā == c || b̄ == c) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.FloatingPointLimit, - left = ā, - right = b̄) + retcode = ReturnCode.FloatingPointLimit, + left = ā, + right = b̄) iszero(fc) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.Success, - left = ā, - right = b̄) + retcode = ReturnCode.Success, + left = ā, + right = b̄) ā, b̄, d̄ = _bracket(f, ā, b̄, c) # The third bracketing block @@ -91,14 +91,14 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, fc = f(c) (ā == c || b̄ == c) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.FloatingPointLimit, - left = ā, - right = b̄) + retcode = ReturnCode.FloatingPointLimit, + left = ā, + right = b̄) iszero(fc) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.Success, - left = ā, - right = b̄) + retcode = ReturnCode.Success, + left = ā, + right = b̄) ā, b̄, d = _bracket(f, ā, b̄, c) # The last bracketing block @@ -110,14 +110,14 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, fc = f(c) (ā == c || b̄ == c) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.FloatingPointLimit, - left = ā, - right = b̄) + retcode = ReturnCode.FloatingPointLimit, + left = ā, + right = b̄) iszero(fc) && return SciMLBase.build_solution(prob, alg, c, fc; - retcode = ReturnCode.Success, - left = ā, - right = b̄) + retcode = ReturnCode.Success, + left = ā, + right = b̄) a, b, d = _bracket(f, ā, b̄, c) end end @@ -132,7 +132,7 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, # Reuturn solution when run out of max interation return SciMLBase.build_solution(prob, alg, c, fc; retcode = ReturnCode.MaxIters, - left = a, right = b) + left = a, right = b) end # Define subrotine function bracket, check fc before bracket to return solution diff --git a/lib/SimpleNonlinearSolve/src/bisection.jl b/lib/SimpleNonlinearSolve/src/bisection.jl index 4e2b232e6..673f77067 100644 --- a/lib/SimpleNonlinearSolve/src/bisection.jl +++ b/lib/SimpleNonlinearSolve/src/bisection.jl @@ -20,16 +20,16 @@ function Bisection(; exact_left = false, exact_right = false) end function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Bisection, args...; - maxiters = 1000, - kwargs...) + maxiters = 1000, + kwargs...) f = Base.Fix2(prob.f, prob.p) left, right = prob.tspan fl, fr = f(left), f(right) if iszero(fl) return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.ExactSolutionLeft, left = left, - right = right) + retcode = ReturnCode.ExactSolutionLeft, left = left, + right = right) end i = 1 @@ -38,8 +38,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Bisection, args... mid = (left + right) / 2 (mid == left || mid == right) && return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.FloatingPointLimit, - left = left, right = right) + retcode = ReturnCode.FloatingPointLimit, + left = left, right = right) fm = f(mid) if iszero(fm) right = mid @@ -60,8 +60,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Bisection, args... mid = (left + right) / 2 (mid == left || mid == right) && return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.FloatingPointLimit, - left = left, right = right) + retcode = ReturnCode.FloatingPointLimit, + left = left, right = right) fm = f(mid) if iszero(fm) right = mid @@ -74,5 +74,5 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Bisection, args... end return SciMLBase.build_solution(prob, alg, left, fl; retcode = ReturnCode.MaxIters, - left = left, right = right) + left = left, right = right) end diff --git a/lib/SimpleNonlinearSolve/src/brent.jl b/lib/SimpleNonlinearSolve/src/brent.jl index 99f645f6a..1cedad134 100644 --- a/lib/SimpleNonlinearSolve/src/brent.jl +++ b/lib/SimpleNonlinearSolve/src/brent.jl @@ -7,8 +7,8 @@ A non-allocating Brent method struct Brent <: AbstractBracketingAlgorithm end function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Brent, args...; - maxiters = 1000, - kwargs...) + maxiters = 1000, + kwargs...) f = Base.Fix2(prob.f, prob.p) a, b = prob.tspan fa, fb = f(a), f(b) @@ -16,8 +16,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Brent, args...; if iszero(fa) return SciMLBase.build_solution(prob, alg, a, fa; - retcode = ReturnCode.ExactSolutionLeft, left = a, - right = b) + retcode = ReturnCode.ExactSolutionLeft, left = a, + right = b) end if abs(fa) < abs(fb) c = b @@ -53,8 +53,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Brent, args...; s = (a + b) / 2 (s == a || s == b) && return SciMLBase.build_solution(prob, alg, a, fa; - retcode = ReturnCode.FloatingPointLimit, - left = a, right = b) + retcode = ReturnCode.FloatingPointLimit, + left = a, right = b) cond = true else cond = false @@ -95,8 +95,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Brent, args...; c = (a + b) / 2 if (c == a || c == b) return SciMLBase.build_solution(prob, alg, a, fa; - retcode = ReturnCode.FloatingPointLimit, - left = a, right = b) + retcode = ReturnCode.FloatingPointLimit, + left = a, right = b) end fc = f(c) if iszero(fc) @@ -110,5 +110,5 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Brent, args...; end return SciMLBase.build_solution(prob, alg, a, fa; retcode = ReturnCode.MaxIters, - left = a, right = b) + left = a, right = b) end diff --git a/lib/SimpleNonlinearSolve/src/broyden.jl b/lib/SimpleNonlinearSolve/src/broyden.jl index e46469619..8ce0d66b5 100644 --- a/lib/SimpleNonlinearSolve/src/broyden.jl +++ b/lib/SimpleNonlinearSolve/src/broyden.jl @@ -16,15 +16,15 @@ struct Broyden{batched, TC <: NLSolveTerminationCondition} <: termination_condition::TC function Broyden(; batched = false, - termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; - abstol = nothing, - reltol = nothing)) + termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; + abstol = nothing, + reltol = nothing)) return new{batched, typeof(termination_condition)}(termination_condition) end end function SciMLBase.__solve(prob::NonlinearProblem, alg::Broyden{false}, args...; - abstol = nothing, reltol = nothing, maxiters = 1000, kwargs...) + abstol = nothing, reltol = nothing, maxiters = 1000, kwargs...) tc = alg.termination_condition mode = DiffEqBase.get_termination_mode(tc) f = Base.Fix2(prob.f, prob.p) diff --git a/lib/SimpleNonlinearSolve/src/dfsane.jl b/lib/SimpleNonlinearSolve/src/dfsane.jl index 3d2e9c8aa..e898f060b 100644 --- a/lib/SimpleNonlinearSolve/src/dfsane.jl +++ b/lib/SimpleNonlinearSolve/src/dfsane.jl @@ -52,15 +52,15 @@ struct SimpleDFSane{T} <: AbstractSimpleNonlinearSolveAlgorithm η_strategy::Function function SimpleDFSane(; σ_min::Real = 1e-10, σ_max::Real = 1e10, σ_1::Real = 1.0, - M::Int = 10, γ::Real = 1e-4, τ_min::Real = 0.1, τ_max::Real = 0.5, - nexp::Int = 2, η_strategy::Function = (f_1, k, x, F) -> f_1 / k^2) + M::Int = 10, γ::Real = 1e-4, τ_min::Real = 0.1, τ_max::Real = 0.5, + nexp::Int = 2, η_strategy::Function = (f_1, k, x, F) -> f_1 / k^2) new{typeof(σ_min)}(σ_min, σ_max, σ_1, M, γ, τ_min, τ_max, nexp, η_strategy) end end function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleDFSane, - args...; abstol = nothing, reltol = nothing, maxiters = 1000, - kwargs...) + args...; abstol = nothing, reltol = nothing, maxiters = 1000, + kwargs...) f = Base.Fix2(prob.f, prob.p) x = float(prob.u0) T = eltype(x) @@ -96,7 +96,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleDFSane, for k in 1:maxiters iszero(F_k) && return SciMLBase.build_solution(prob, alg, x, F_k; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) # Spectral parameter range check if abs(σ_k) > σ_max @@ -136,7 +136,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, alg::SimpleDFSane, if isapprox(x_new, x, atol = atol, rtol = rtol) return SciMLBase.build_solution(prob, alg, x_new, F_new; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) end # Update spectral parameter s_k = x_new - x diff --git a/lib/SimpleNonlinearSolve/src/falsi.jl b/lib/SimpleNonlinearSolve/src/falsi.jl index 81ce68ffa..cce11a811 100644 --- a/lib/SimpleNonlinearSolve/src/falsi.jl +++ b/lib/SimpleNonlinearSolve/src/falsi.jl @@ -4,16 +4,16 @@ struct Falsi <: AbstractBracketingAlgorithm end function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Falsi, args...; - maxiters = 1000, - kwargs...) + maxiters = 1000, + kwargs...) f = Base.Fix2(prob.f, prob.p) left, right = prob.tspan fl, fr = f(left), f(right) if iszero(fl) return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.ExactSolutionLeft, left = left, - right = right) + retcode = ReturnCode.ExactSolutionLeft, left = left, + right = right) end i = 1 @@ -21,8 +21,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Falsi, args...; while i < maxiters if nextfloat_tdir(left, prob.tspan...) == right return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.FloatingPointLimit, - left = left, right = right) + retcode = ReturnCode.FloatingPointLimit, + left = left, right = right) end mid = (fr * left - fl * right) / (fr - fl) for i in 1:10 @@ -51,8 +51,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Falsi, args...; mid = (left + right) / 2 (mid == left || mid == right) && return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.FloatingPointLimit, - left = left, right = right) + retcode = ReturnCode.FloatingPointLimit, + left = left, right = right) fm = f(mid) if iszero(fm) right = mid @@ -68,5 +68,5 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Falsi, args...; end return SciMLBase.build_solution(prob, alg, left, fl; retcode = ReturnCode.MaxIters, - left = left, right = right) + left = left, right = right) end diff --git a/lib/SimpleNonlinearSolve/src/halley.jl b/lib/SimpleNonlinearSolve/src/halley.jl index cdda7beda..9e97a57cd 100644 --- a/lib/SimpleNonlinearSolve/src/halley.jl +++ b/lib/SimpleNonlinearSolve/src/halley.jl @@ -30,16 +30,16 @@ and static array problems. """ struct Halley{CS, AD, FDT} <: AbstractNewtonAlgorithm{CS, AD, FDT} function Halley(; chunk_size = Val{0}(), autodiff = Val{true}(), - diff_type = Val{:forward}) + diff_type = Val{:forward}) new{SciMLBase._unwrap_val(chunk_size), SciMLBase._unwrap_val(autodiff), SciMLBase._unwrap_val(diff_type)}() end end function SciMLBase.__solve(prob::NonlinearProblem, - alg::Halley, args...; abstol = nothing, - reltol = nothing, - maxiters = 1000, kwargs...) + alg::Halley, args...; abstol = nothing, + reltol = nothing, + maxiters = 1000, kwargs...) f = Base.Fix2(prob.f, prob.p) x = float(prob.u0) fx = f(x) @@ -81,18 +81,18 @@ function SciMLBase.__solve(prob::NonlinearProblem, if isa(x, Number) fx = f(x) dfx = FiniteDiff.finite_difference_derivative(f, x, diff_type(alg), - eltype(x)) + eltype(x)) d2fx = FiniteDiff.finite_difference_derivative(x -> FiniteDiff.finite_difference_derivative(f, - x), - x, - diff_type(alg), eltype(x)) + x), + x, + diff_type(alg), eltype(x)) else fx = f(x) dfx = FiniteDiff.finite_difference_jacobian(f, x, diff_type(alg), eltype(x)) d2fx = FiniteDiff.finite_difference_jacobian(x -> FiniteDiff.finite_difference_jacobian(f, - x), - x, - diff_type(alg), eltype(x)) + x), + x, + diff_type(alg), eltype(x)) ai = -(dfx \ fx) A = reshape(d2fx * ai, (n, n)) bi = (dfx) \ (A * ai) diff --git a/lib/SimpleNonlinearSolve/src/klement.jl b/lib/SimpleNonlinearSolve/src/klement.jl index 4d59377df..00264d32f 100644 --- a/lib/SimpleNonlinearSolve/src/klement.jl +++ b/lib/SimpleNonlinearSolve/src/klement.jl @@ -9,9 +9,9 @@ This method is non-allocating on scalar problems. struct Klement <: AbstractSimpleNonlinearSolveAlgorithm end function SciMLBase.__solve(prob::NonlinearProblem, - alg::Klement, args...; abstol = nothing, - reltol = nothing, - maxiters = 1000, kwargs...) + alg::Klement, args...; abstol = nothing, + reltol = nothing, + maxiters = 1000, kwargs...) f = Base.Fix2(prob.f, prob.p) x = float(prob.u0) fₙ = f(x) @@ -39,11 +39,11 @@ function SciMLBase.__solve(prob::NonlinearProblem, iszero(fₙ) && return SciMLBase.build_solution(prob, alg, xₙ, fₙ; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) if isapprox(xₙ, xₙ₋₁, atol = atol, rtol = rtol) return SciMLBase.build_solution(prob, alg, xₙ, fₙ; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) end Δxₙ = xₙ - xₙ₋₁ @@ -81,11 +81,11 @@ function SciMLBase.__solve(prob::NonlinearProblem, iszero(fₙ) && return SciMLBase.build_solution(prob, alg, xₙ, fₙ; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) if isapprox(xₙ, xₙ₋₁, atol = atol, rtol = rtol) return SciMLBase.build_solution(prob, alg, xₙ, fₙ; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) end Δxₙ = xₙ - xₙ₋₁ diff --git a/lib/SimpleNonlinearSolve/src/lbroyden.jl b/lib/SimpleNonlinearSolve/src/lbroyden.jl index 6cca7c11e..fc2b51a88 100644 --- a/lib/SimpleNonlinearSolve/src/lbroyden.jl +++ b/lib/SimpleNonlinearSolve/src/lbroyden.jl @@ -18,16 +18,16 @@ struct LBroyden{batched, TC <: NLSolveTerminationCondition} <: threshold::Int function LBroyden(; batched = false, threshold::Int = 27, - termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; - abstol = nothing, - reltol = nothing)) + termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; + abstol = nothing, + reltol = nothing)) return new{batched, typeof(termination_condition)}(termination_condition, threshold) end end @views function SciMLBase.__solve(prob::NonlinearProblem, alg::LBroyden{batched}, args...; - abstol = nothing, reltol = nothing, maxiters = 1000, - kwargs...) where {batched} + abstol = nothing, reltol = nothing, maxiters = 1000, + kwargs...) where {batched} tc = alg.termination_condition mode = DiffEqBase.get_termination_mode(tc) threshold = min(maxiters, alg.threshold) @@ -93,7 +93,7 @@ end selectdim(U, 1, mod1(i, threshold)) .= u update = -_matvec(selectdim(U, 1, 1:min(threshold, i + 1)), - selectdim(Vᵀ, 2, 1:min(threshold, i + 1)), fₙ) + selectdim(Vᵀ, 2, 1:min(threshold, i + 1)), fₙ) xₙ₋₁ = xₙ fₙ₋₁ = fₙ @@ -116,26 +116,26 @@ function _init_lbroyden_state(batched::Bool, x, threshold) end function _rmatvec(U::AbstractMatrix, Vᵀ::AbstractMatrix, - x::Union{<:AbstractVector, <:Number}) + x::Union{<:AbstractVector, <:Number}) length(U) == 0 && return x return -x .+ vec((x' * Vᵀ) * U) end function _rmatvec(U::AbstractArray{T1, 3}, Vᵀ::AbstractArray{T2, 3}, - x::AbstractMatrix) where {T1, T2} + x::AbstractMatrix) where {T1, T2} length(U) == 0 && return x Vᵀx = sum(Vᵀ .* reshape(x, size(x, 1), 1, size(x, 2)); dims = 1) return -x .+ _drdims_sum(U .* permutedims(Vᵀx, (2, 1, 3)); dims = 1) end function _matvec(U::AbstractMatrix, Vᵀ::AbstractMatrix, - x::Union{<:AbstractVector, <:Number}) + x::Union{<:AbstractVector, <:Number}) length(U) == 0 && return x return -x .+ vec(Vᵀ * (U * x)) end function _matvec(U::AbstractArray{T1, 3}, Vᵀ::AbstractArray{T2, 3}, - x::AbstractMatrix) where {T1, T2} + x::AbstractMatrix) where {T1, T2} length(U) == 0 && return x xUᵀ = sum(reshape(x, size(x, 1), 1, size(x, 2)) .* permutedims(U, (2, 1, 3)); dims = 1) return -x .+ _drdims_sum(xUᵀ .* Vᵀ; dims = 2) diff --git a/lib/SimpleNonlinearSolve/src/raphson.jl b/lib/SimpleNonlinearSolve/src/raphson.jl index 8d2d9796c..386c35053 100644 --- a/lib/SimpleNonlinearSolve/src/raphson.jl +++ b/lib/SimpleNonlinearSolve/src/raphson.jl @@ -30,16 +30,16 @@ and static array problems. """ struct SimpleNewtonRaphson{CS, AD, FDT} <: AbstractNewtonAlgorithm{CS, AD, FDT} function SimpleNewtonRaphson(; chunk_size = Val{0}(), autodiff = Val{true}(), - diff_type = Val{:forward}) + diff_type = Val{:forward}) new{SciMLBase._unwrap_val(chunk_size), SciMLBase._unwrap_val(autodiff), SciMLBase._unwrap_val(diff_type)}() end end function SciMLBase.__solve(prob::NonlinearProblem, - alg::SimpleNewtonRaphson, args...; abstol = nothing, - reltol = nothing, - maxiters = 1000, kwargs...) + alg::SimpleNewtonRaphson, args...; abstol = nothing, + reltol = nothing, + maxiters = 1000, kwargs...) f = Base.Fix2(prob.f, prob.p) x = float(prob.u0) fx = float(prob.u0) @@ -71,7 +71,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, else fx = f(x) dfx = FiniteDiff.finite_difference_derivative(f, x, diff_type(alg), eltype(x), - fx) + fx) end iszero(fx) && return SciMLBase.build_solution(prob, alg, x, fx; retcode = ReturnCode.Success) diff --git a/lib/SimpleNonlinearSolve/src/ridder.jl b/lib/SimpleNonlinearSolve/src/ridder.jl index 9fe7bccfe..62b5a931a 100644 --- a/lib/SimpleNonlinearSolve/src/ridder.jl +++ b/lib/SimpleNonlinearSolve/src/ridder.jl @@ -7,16 +7,16 @@ A non-allocating ridder method struct Ridder <: AbstractBracketingAlgorithm end function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Ridder, args...; - maxiters = 1000, - kwargs...) + maxiters = 1000, + kwargs...) f = Base.Fix2(prob.f, prob.p) left, right = prob.tspan fl, fr = f(left), f(right) if iszero(fl) return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.ExactSolutionLeft, left = left, - right = right) + retcode = ReturnCode.ExactSolutionLeft, left = left, + right = right) end xo = oftype(left, Inf) @@ -26,14 +26,14 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Ridder, args...; mid = (left + right) / 2 (mid == left || mid == right) && return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.FloatingPointLimit, - left = left, right = right) + retcode = ReturnCode.FloatingPointLimit, + left = left, right = right) fm = f(mid) s = sqrt(fm^2 - fl * fr) iszero(s) && return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.Failure, - left = left, right = right) + retcode = ReturnCode.Failure, + left = left, right = right) x = mid + (mid - left) * sign(fl - fr) * fm / s fx = f(x) xo = x @@ -63,8 +63,8 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Ridder, args...; mid = (left + right) / 2 (mid == left || mid == right) && return SciMLBase.build_solution(prob, alg, left, fl; - retcode = ReturnCode.FloatingPointLimit, - left = left, right = right) + retcode = ReturnCode.FloatingPointLimit, + left = left, right = right) fm = f(mid) if iszero(fm) right = mid @@ -77,5 +77,5 @@ function SciMLBase.solve(prob::IntervalNonlinearProblem, alg::Ridder, args...; end return SciMLBase.build_solution(prob, alg, left, fl; retcode = ReturnCode.MaxIters, - left = left, right = right) + left = left, right = right) end diff --git a/lib/SimpleNonlinearSolve/src/trustRegion.jl b/lib/SimpleNonlinearSolve/src/trustRegion.jl index aa893ffc1..1423d59b8 100644 --- a/lib/SimpleNonlinearSolve/src/trustRegion.jl +++ b/lib/SimpleNonlinearSolve/src/trustRegion.jl @@ -64,34 +64,34 @@ struct SimpleTrustRegion{T, CS, AD, FDT} <: AbstractNewtonAlgorithm{CS, AD, FDT} expand_factor::T max_shrink_times::Int function SimpleTrustRegion(; chunk_size = Val{0}(), - autodiff = Val{true}(), - diff_type = Val{:forward}, - max_trust_radius::Real = 0.0, - initial_trust_radius::Real = 0.0, - step_threshold::Real = 0.1, - shrink_threshold::Real = 0.25, - expand_threshold::Real = 0.75, - shrink_factor::Real = 0.25, - expand_factor::Real = 2.0, - max_shrink_times::Int = 32) + autodiff = Val{true}(), + diff_type = Val{:forward}, + max_trust_radius::Real = 0.0, + initial_trust_radius::Real = 0.0, + step_threshold::Real = 0.1, + shrink_threshold::Real = 0.25, + expand_threshold::Real = 0.75, + shrink_factor::Real = 0.25, + expand_factor::Real = 2.0, + max_shrink_times::Int = 32) new{typeof(initial_trust_radius), SciMLBase._unwrap_val(chunk_size), SciMLBase._unwrap_val(autodiff), SciMLBase._unwrap_val(diff_type)}(max_trust_radius, - initial_trust_radius, - step_threshold, - shrink_threshold, - expand_threshold, - shrink_factor, - expand_factor, - max_shrink_times) + initial_trust_radius, + step_threshold, + shrink_threshold, + expand_threshold, + shrink_factor, + expand_factor, + max_shrink_times) end end function SciMLBase.__solve(prob::NonlinearProblem, - alg::SimpleTrustRegion, args...; abstol = nothing, - reltol = nothing, - maxiters = 1000, kwargs...) + alg::SimpleTrustRegion, args...; abstol = nothing, + reltol = nothing, + maxiters = 1000, kwargs...) f = Base.Fix2(prob.f, prob.p) x = float(prob.u0) T = typeof(x) @@ -112,7 +112,6 @@ function SciMLBase.__solve(prob::NonlinearProblem, real(oneunit(eltype(T))) * (eps(real(one(eltype(T)))))^(4 // 5) rtol = reltol !== nothing ? reltol : eps(real(one(eltype(T))))^(4 // 5) - if DiffEqBase.has_jac(prob.f) ∇f = prob.f.jac(x, prob.p) F = f(x) @@ -156,7 +155,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, shrink_counter += 1 if shrink_counter > max_shrink_times return SciMLBase.build_solution(prob, alg, x, F; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) end else shrink_counter = 0 @@ -164,7 +163,7 @@ function SciMLBase.__solve(prob::NonlinearProblem, if r > η₁ if isapprox(xₖ₊₁, x, atol = atol, rtol = rtol) return SciMLBase.build_solution(prob, alg, xₖ₊₁, Fₖ₊₁; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) end # Take the step. x = xₖ₊₁ @@ -173,16 +172,16 @@ function SciMLBase.__solve(prob::NonlinearProblem, F, ∇f = value_derivative(f, x) elseif x isa AbstractArray ∇f = FiniteDiff.finite_difference_jacobian(f, x, diff_type(alg), eltype(x), - F) + F) else ∇f = FiniteDiff.finite_difference_derivative(f, x, diff_type(alg), - eltype(x), - F) + eltype(x), + F) end iszero(F) && return SciMLBase.build_solution(prob, alg, x, F; - retcode = ReturnCode.Success) + retcode = ReturnCode.Success) # Update the trust region radius. if r > η₃ && norm(δ) ≈ Δ diff --git a/lib/SimpleNonlinearSolve/test/basictests.jl b/lib/SimpleNonlinearSolve/test/basictests.jl index 1a0249b65..e1d6f9bf8 100644 --- a/lib/SimpleNonlinearSolve/test/basictests.jl +++ b/lib/SimpleNonlinearSolve/test/basictests.jl @@ -18,7 +18,7 @@ for mode in instances(NLSolveTerminationMode.T) end termination_condition = NLSolveTerminationCondition(mode; abstol = nothing, - reltol = nothing) + reltol = nothing) push!(BROYDEN_SOLVERS, Broyden(; batched = false, termination_condition)) push!(BATCHED_BROYDEN_SOLVERS, Broyden(; batched = true, termination_condition)) push!(LBROYDEN_SOLVERS, LBroyden(; batched = false, termination_condition)) @@ -131,7 +131,7 @@ using ForwardDiff f, u0 = (u, p) -> u .* u .- p, @SVector[1.0, 1.0] for alg in (SimpleNewtonRaphson(), LBroyden(), Klement(), SimpleTrustRegion(), - SimpleDFSane(), Halley(), BROYDEN_SOLVERS...) + SimpleDFSane(), Halley(), BROYDEN_SOLVERS...) g = function (p) probN = NonlinearProblem{false}(f, csu0, p) sol = solve(probN, alg, abstol = 1e-9) @@ -154,7 +154,7 @@ end # Scalar f, u0 = (u, p) -> u * u - p, 1.0 for alg in (SimpleNewtonRaphson(), Klement(), SimpleTrustRegion(), - SimpleDFSane(), Halley(), BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) + SimpleDFSane(), Halley(), BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) g = function (p) probN = NonlinearProblem{false}(f, oftype(p, u0), p) sol = solve(probN, alg) @@ -251,7 +251,7 @@ for alg in [Bisection(), Falsi(), Ridder(), Brent()] end for alg in (SimpleNewtonRaphson(), Klement(), SimpleTrustRegion(), - SimpleDFSane(), Halley(), BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) + SimpleDFSane(), Halley(), BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) global g, p g = function (p) probN = NonlinearProblem{false}(f, 0.5, p) @@ -267,10 +267,10 @@ f, u0 = (u, p) -> u .* u .- 2.0, @SVector[1.0, 1.0] probN = NonlinearProblem(f, u0) for alg in (SimpleNewtonRaphson(), SimpleNewtonRaphson(; autodiff = false), - SimpleTrustRegion(), - SimpleTrustRegion(; autodiff = false), Halley(), Halley(; autodiff = false), - Klement(), SimpleDFSane(), - BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) + SimpleTrustRegion(), + SimpleTrustRegion(; autodiff = false), Halley(), Halley(; autodiff = false), + Klement(), SimpleDFSane(), + BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) sol = solve(probN, alg) @test sol.retcode == ReturnCode.Success @@ -284,8 +284,8 @@ for u0 in [1.0, [1, 1.0]] sol = sqrt(2) * u0 for alg in (SimpleNewtonRaphson(), SimpleNewtonRaphson(; autodiff = false), - SimpleTrustRegion(), SimpleTrustRegion(; autodiff = false), Klement(), - SimpleDFSane(), BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) + SimpleTrustRegion(), SimpleTrustRegion(; autodiff = false), Klement(), + SimpleDFSane(), BROYDEN_SOLVERS..., LBROYDEN_SOLVERS...) sol2 = solve(probN, alg) @test sol2.retcode == ReturnCode.Success @@ -399,18 +399,18 @@ expand_factor = [1.5, 2.0, 3.0] max_shrink_times = [10, 20, 30] list_of_options = zip(max_trust_radius, initial_trust_radius, step_threshold, - shrink_threshold, expand_threshold, shrink_factor, - expand_factor, max_shrink_times) + shrink_threshold, expand_threshold, shrink_factor, + expand_factor, max_shrink_times) for options in list_of_options local probN, sol, alg alg = SimpleTrustRegion(max_trust_radius = options[1], - initial_trust_radius = options[2], - step_threshold = options[3], - shrink_threshold = options[4], - expand_threshold = options[5], - shrink_factor = options[6], - expand_factor = options[7], - max_shrink_times = options[8]) + initial_trust_radius = options[2], + step_threshold = options[3], + shrink_threshold = options[4], + expand_threshold = options[5], + shrink_factor = options[6], + expand_factor = options[7], + max_shrink_times = options[8]) probN = NonlinearProblem(f, u0, p) sol = solve(probN, alg) @@ -454,18 +454,18 @@ nexp = [2, 1, 2] ] list_of_options = zip(σ_min, σ_max, σ_1, M, γ, τ_min, τ_max, nexp, - η_strategy) + η_strategy) for options in list_of_options local probN, sol, alg alg = SimpleDFSane(σ_min = options[1], - σ_max = options[2], - σ_1 = options[3], - M = options[4], - γ = options[5], - τ_min = options[6], - τ_max = options[7], - nexp = options[8], - η_strategy = options[9]) + σ_max = options[2], + σ_1 = options[3], + M = options[4], + γ = options[5], + τ_min = options[6], + τ_max = options[7], + nexp = options[8], + η_strategy = options[9]) probN = NonlinearProblem(f, u0, p) sol = solve(probN, alg) @@ -495,7 +495,9 @@ end f, u0 = (u, p) -> u .* u .- p, randn(3) -f_jac(u, p) = begin diagm(2 * u) end +f_jac(u, p) = begin + diagm(2 * u) +end p = [2.0, 1.0, 5.0]; diff --git a/lib/SimpleNonlinearSolve/test/runtests.jl b/lib/SimpleNonlinearSolve/test/runtests.jl index 7269ed6bf..94a0086e4 100644 --- a/lib/SimpleNonlinearSolve/test/runtests.jl +++ b/lib/SimpleNonlinearSolve/test/runtests.jl @@ -6,7 +6,9 @@ const GROUP = get(ENV, "GROUP", "All") const is_APPVEYOR = Sys.iswindows() && haskey(ENV, "APPVEYOR") @time begin - -if GROUP == "All" || GROUP == "Core" - @time @safetestset "Basic Tests + Some AD" begin include("basictests.jl") end -end end + if GROUP == "All" || GROUP == "Core" + @time @safetestset "Basic Tests + Some AD" begin + include("basictests.jl") + end + end +end