From 88856385e9d86354f851ea9667dd869d05d7cb65 Mon Sep 17 00:00:00 2001 From: Utkarsh Date: Mon, 25 Sep 2023 21:18:59 -0400 Subject: [PATCH] Update all algorithms to use termination condition --- src/NonlinearSolve.jl | 2 +- src/levenberg.jl | 59 +++++++++++++++++++++++------ src/raphson.jl | 40 ++++++++++++-------- src/trustRegion.jl | 86 ++++++++++++++++++++++++++++++++----------- test/basictests.jl | 42 +++++++++++++++++++++ 5 files changed, 179 insertions(+), 50 deletions(-) diff --git a/src/NonlinearSolve.jl b/src/NonlinearSolve.jl index f8ae69939..e3ee2af81 100644 --- a/src/NonlinearSolve.jl +++ b/src/NonlinearSolve.jl @@ -26,7 +26,7 @@ const AbstractSparseADType = Union{ADTypes.AbstractSparseFiniteDifferences, ADTypes.AbstractSparseForwardMode, ADTypes.AbstractSparseReverseMode} abstract type AbstractNonlinearSolveAlgorithm <: AbstractNonlinearAlgorithm end -abstract type AbstractNewtonAlgorithm{CJ, AD} <: AbstractNonlinearSolveAlgorithm end +abstract type AbstractNewtonAlgorithm{CJ, AD, TC} <: AbstractNonlinearSolveAlgorithm end abstract type AbstractNonlinearSolveCache{iip} end diff --git a/src/levenberg.jl b/src/levenberg.jl index 1a0343010..b2cd6c6d0 100644 --- a/src/levenberg.jl +++ b/src/levenberg.jl @@ -73,7 +73,8 @@ numerically-difficult nonlinear systems. [this paper](https://arxiv.org/abs/1201.5885) to use a minimum value of the elements in `DᵀD` to prevent the damping from being too small. Defaults to `1e-8`. """ -@concrete struct LevenbergMarquardt{CJ, AD, T} <: AbstractNewtonAlgorithm{CJ, AD} +@concrete struct LevenbergMarquardt{CJ, AD, T, TC <: NLSolveTerminationCondition} <: + AbstractNewtonAlgorithm{CJ, AD, TC} ad::AD linsolve precs @@ -84,23 +85,29 @@ numerically-difficult nonlinear systems. α_geodesic::T b_uphill::T min_damping_D::T + termination_condition::TC end function LevenbergMarquardt(; concrete_jac = nothing, linsolve = nothing, precs = DEFAULT_PRECS, damping_initial::Real = 1.0, damping_increase_factor::Real = 2.0, damping_decrease_factor::Real = 3.0, finite_diff_step_geodesic::Real = 0.1, α_geodesic::Real = 0.75, b_uphill::Real = 1.0, min_damping_D::AbstractFloat = 1e-8, + termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; + abstol = nothing, + reltol = nothing), adkwargs...) ad = default_adargs_to_adtype(; adkwargs...) return LevenbergMarquardt{_unwrap_val(concrete_jac)}(ad, linsolve, precs, damping_initial, damping_increase_factor, damping_decrease_factor, - finite_diff_step_geodesic, α_geodesic, b_uphill, min_damping_D) + finite_diff_step_geodesic, α_geodesic, b_uphill, min_damping_D, + termination_condition) end @concrete mutable struct LevenbergMarquardtCache{iip} <: AbstractNonlinearSolveCache{iip} f alg u + u_prev fu1 fu2 du @@ -114,6 +121,7 @@ end internalnorm retcode::ReturnCode.T abstol + reltol prob DᵀD JᵀJ @@ -138,11 +146,12 @@ end Jv mat_tmp stats::NLStats + tc_storage end -function SciMLBase.__init(prob::Union{NonlinearProblem{uType, iip}, - NonlinearLeastSquaresProblem{uType, iip}}, alg::LevenbergMarquardt, - args...; alias_u0 = false, maxiters = 1000, abstol = 1e-6, internalnorm = DEFAULT_NORM, +function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::LevenbergMarquardt, + args...; alias_u0 = false, maxiters = 1000, abstol = nothing, reltol = nothing, + internalnorm = DEFAULT_NORM, linsolve_kwargs = (;), kwargs...) where {uType, iip} @unpack f, u0, p = prob u = alias_u0 ? u0 : deepcopy(u0) @@ -176,15 +185,30 @@ function SciMLBase.__init(prob::Union{NonlinearProblem{uType, iip}, fu_tmp = zero(fu1) mat_tmp = zero(JᵀJ) - return LevenbergMarquardtCache{iip}(f, alg, u, fu1, fu2, du, p, uf, linsolve, J, - jac_cache, false, maxiters, internalnorm, ReturnCode.Default, abstol, prob, DᵀD, + tc = alg.termination_condition + mode = DiffEqBase.get_termination_mode(tc) + + atol = _get_tolerance(abstol, tc.abstol, eltype(u)) + rtol = _get_tolerance(reltol, tc.reltol, eltype(u)) + + storage = mode ∈ DiffEqBase.SAFE_TERMINATION_MODES ? NLSolveSafeTerminationResult() : + nothing + + return LevenbergMarquardtCache{iip}(f, alg, u, copy(u), fu1, fu2, du, p, uf, linsolve, + J, + jac_cache, false, maxiters, internalnorm, ReturnCode.Default, atol, rtol, prob, DᵀD, JᵀJ, λ, λ_factor, damping_increase_factor, damping_decrease_factor, h, α_geodesic, b_uphill, min_damping_D, v, a, tmp_vec, v_old, loss, δ, loss, make_new_J, fu_tmp, - zero(u), zero(fu1), mat_tmp, NLStats(1, 0, 0, 0, 0)) + zero(u), zero(fu1), mat_tmp, NLStats(1, 0, 0, 0, 0), storage) + end function perform_step!(cache::LevenbergMarquardtCache{true}) @unpack fu1, f, make_new_J = cache + + tc_storage = cache.tc_storage + termination_condition = cache.alg.termination_condition(tc_storage) + if iszero(fu1) cache.force_stop = true return nothing @@ -197,7 +221,7 @@ function perform_step!(cache::LevenbergMarquardtCache{true}) cache.make_new_J = false cache.stats.njacs += 1 end - @unpack u, p, λ, JᵀJ, DᵀD, J, alg, linsolve = cache + @unpack u, u_prev, p, λ, JᵀJ, DᵀD, J, alg, linsolve = cache # Usual Levenberg-Marquardt step ("velocity"). # The following lines do: cache.v = -cache.mat_tmp \ cache.u_tmp @@ -237,7 +261,11 @@ function perform_step!(cache::LevenbergMarquardtCache{true}) if (1 - β)^b_uphill * loss ≤ loss_old # Accept step. cache.u .+= δ - if loss < cache.abstol + if termination_condition(cache.fu_tmp, + cache.u, + u_prev, + cache.abstol, + cache.reltol) cache.force_stop = true return nothing end @@ -249,6 +277,7 @@ function perform_step!(cache::LevenbergMarquardtCache{true}) cache.make_new_J = true end end + @. u_prev = u cache.λ *= cache.λ_factor cache.λ_factor = cache.damping_increase_factor return nothing @@ -256,6 +285,10 @@ end function perform_step!(cache::LevenbergMarquardtCache{false}) @unpack fu1, f, make_new_J = cache + + tc_storage = cache.tc_storage + termination_condition = cache.alg.termination_condition(tc_storage) + if iszero(fu1) cache.force_stop = true return nothing @@ -272,7 +305,8 @@ function perform_step!(cache::LevenbergMarquardtCache{false}) cache.make_new_J = false cache.stats.njacs += 1 end - @unpack u, p, λ, JᵀJ, DᵀD, J, linsolve, alg = cache + + @unpack u, u_prev, p, λ, JᵀJ, DᵀD, J, linsolve, alg = cache cache.mat_tmp = JᵀJ + λ * DᵀD # Usual Levenberg-Marquardt step ("velocity"). @@ -313,7 +347,7 @@ function perform_step!(cache::LevenbergMarquardtCache{false}) if (1 - β)^b_uphill * loss ≤ loss_old # Accept step. cache.u += δ - if loss < cache.abstol + if termination_condition(fu_new, cache.u, u_prev, cache.abstol, cache.reltol) cache.force_stop = true return nothing end @@ -325,6 +359,7 @@ function perform_step!(cache::LevenbergMarquardtCache{false}) cache.make_new_J = true end end + cache.u_prev = @. cache.u cache.λ *= cache.λ_factor cache.λ_factor = cache.damping_increase_factor return nothing diff --git a/src/raphson.jl b/src/raphson.jl index ad64cab41..60340da28 100644 --- a/src/raphson.jl +++ b/src/raphson.jl @@ -29,7 +29,8 @@ for large-scale and numerically-difficult nonlinear systems. which means that no line search is performed. Algorithms from `LineSearches.jl` can be used here directly, and they will be converted to the correct `LineSearch`. """ -@concrete struct NewtonRaphson{CJ, AD, TC <: NLSolveTerminationCondition} <: AbstractNewtonAlgorithm{CJ, AD} +@concrete struct NewtonRaphson{CJ, AD, TC <: NLSolveTerminationCondition} <: + AbstractNewtonAlgorithm{CJ, AD, TC} ad::AD linsolve precs @@ -38,19 +39,24 @@ for large-scale and numerically-difficult nonlinear systems. end function NewtonRaphson(; concrete_jac = nothing, linsolve = nothing, - linesearch = LineSearch(), precs = DEFAULT_PRECS, termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; - abstol = nothing, - reltol = nothing), adkwargs...) + linesearch = LineSearch(), precs = DEFAULT_PRECS, + termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; + abstol = nothing, + reltol = nothing), adkwargs...) ad = default_adargs_to_adtype(; adkwargs...) linesearch = linesearch isa LineSearch ? linesearch : LineSearch(; method = linesearch) - return NewtonRaphson{_unwrap_val(concrete_jac)}(ad, linsolve, precs, linesearch, termination_condition) + return NewtonRaphson{_unwrap_val(concrete_jac)}(ad, + linsolve, + precs, + linesearch, + termination_condition) end @concrete mutable struct NewtonRaphsonCache{iip} <: AbstractNonlinearSolveCache{iip} f alg u - uprev + u_prev fu1 fu2 du @@ -72,7 +78,8 @@ end end function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::NewtonRaphson, args...; - alias_u0 = false, maxiters = 1000, abstol = 1e-6, internalnorm = DEFAULT_NORM, + alias_u0 = false, maxiters = 1000, abstol = nothing, reltol = nothing, + internalnorm = DEFAULT_NORM, linsolve_kwargs = (;), kwargs...) where {uType, iip} @unpack f, u0, p = prob u = alias_u0 ? u0 : deepcopy(u0) @@ -80,7 +87,6 @@ function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::NewtonRaphson uf, linsolve, J, fu2, jac_cache, du = jacobian_caches(alg, f, u, p, Val(iip); linsolve_kwargs) - tc = alg.termination_condition mode = DiffEqBase.get_termination_mode(tc) @@ -92,11 +98,12 @@ function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::NewtonRaphson return NewtonRaphsonCache{iip}(f, alg, u, copy(u), fu1, fu2, du, p, uf, linsolve, J, jac_cache, false, maxiters, internalnorm, ReturnCode.Default, atol, rtol, prob, - NLStats(1, 0, 0, 0, 0), LineSearchCache(alg.linesearch, f, u, p, fu1, Val(iip)), storage) + NLStats(1, 0, 0, 0, 0), LineSearchCache(alg.linesearch, f, u, p, fu1, Val(iip)), + storage) end function perform_step!(cache::NewtonRaphsonCache{true}) - @unpack u, uprev, fu1, f, p, alg, J, linsolve, du = cache + @unpack u, u_prev, fu1, f, p, alg, J, linsolve, du = cache jacobian!!(J, cache) tc_storage = cache.tc_storage @@ -112,9 +119,10 @@ function perform_step!(cache::NewtonRaphsonCache{true}) @. u = u - α * du f(cache.fu1, u, p) - termination_condition(cache.fu1, u, uprev, cache.abstol, cache.reltol) && (cache.force_stop = true) + termination_condition(cache.fu1, u, u_prev, cache.abstol, cache.reltol) && + (cache.force_stop = true) - @. uprev = u + @. u_prev = u cache.stats.nf += 1 cache.stats.njacs += 1 cache.stats.nsolve += 1 @@ -123,12 +131,11 @@ function perform_step!(cache::NewtonRaphsonCache{true}) end function perform_step!(cache::NewtonRaphsonCache{false}) - @unpack u, uprev, fu1, f, p, alg, linsolve = cache + @unpack u, u_prev, fu1, f, p, alg, linsolve = cache tc_storage = cache.tc_storage termination_condition = cache.alg.termination_condition(tc_storage) - cache.J = jacobian!!(cache.J, cache) # u = u - J \ fu if linsolve === nothing @@ -144,9 +151,10 @@ function perform_step!(cache::NewtonRaphsonCache{false}) cache.u = @. u - α * cache.du # `u` might not support mutation cache.fu1 = f(cache.u, p) - termination_condition(cache.fu1, cache.u, uprev, cache.abstol, cache.reltol) && (cache.force_stop = true) + termination_condition(cache.fu1, cache.u, u_prev, cache.abstol, cache.reltol) && + (cache.force_stop = true) - cache.uprev = cache.u + cache.u_prev = @. cache.u cache.stats.nf += 1 cache.stats.njacs += 1 cache.stats.nsolve += 1 diff --git a/src/trustRegion.jl b/src/trustRegion.jl index 50d1078f3..f0d1295e3 100644 --- a/src/trustRegion.jl +++ b/src/trustRegion.jl @@ -28,14 +28,14 @@ states as `RadiusUpdateSchemes.T`. Simply put the desired scheme as follows: """ `RadiusUpdateSchemes.NLsolve` - The same updating scheme as in NLsolve's (https://github.com/JuliaNLSolvers/NLsolve.jl) trust region dogleg implementation. + The same updating scheme as in NLsolve's (https://github.com/JuliaNLSolvers/NLsolve.jl) trust region dogleg implementation. """ NLsolve """ `RadiusUpdateSchemes.NocedalWright` - Trust region updating scheme as in Nocedal and Wright [see Alg 11.5, page 291]. + Trust region updating scheme as in Nocedal and Wright [see Alg 11.5, page 291]. """ NocedalWright @@ -147,7 +147,8 @@ for large-scale and numerically-difficult nonlinear systems. `linsolve` and `precs` are used exclusively for the inplace version of the algorithm. Support for the OOP version is planned! """ -@concrete struct TrustRegion{CJ, AD, MTR} <: AbstractNewtonAlgorithm{CJ, AD} +@concrete struct TrustRegion{CJ, AD, MTR, TC <: NLSolveTerminationCondition} <: + AbstractNewtonAlgorithm{CJ, AD, TC} ad::AD linsolve precs @@ -160,6 +161,7 @@ for large-scale and numerically-difficult nonlinear systems. shrink_factor::MTR expand_factor::MTR max_shrink_times::Int + termination_condition::TC end function TrustRegion(; concrete_jac = nothing, linsolve = nothing, precs = DEFAULT_PRECS, @@ -167,11 +169,15 @@ function TrustRegion(; concrete_jac = nothing, linsolve = nothing, precs = DEFAU max_trust_radius::Real = 0 // 1, initial_trust_radius::Real = 0 // 1, step_threshold::Real = 1 // 10000, shrink_threshold::Real = 1 // 4, expand_threshold::Real = 3 // 4, shrink_factor::Real = 1 // 4, - expand_factor::Real = 2 // 1, max_shrink_times::Int = 32, adkwargs...) + expand_factor::Real = 2 // 1, max_shrink_times::Int = 32, + termination_condition = NLSolveTerminationCondition(NLSolveTerminationMode.NLSolveDefault; + abstol = nothing, + reltol = nothing), adkwargs...) ad = default_adargs_to_adtype(; adkwargs...) return TrustRegion{_unwrap_val(concrete_jac)}(ad, linsolve, precs, radius_update_scheme, max_trust_radius, initial_trust_radius, step_threshold, shrink_threshold, - expand_threshold, shrink_factor, expand_factor, max_shrink_times) + expand_threshold, shrink_factor, expand_factor, max_shrink_times, + termination_condition) end @concrete mutable struct TrustRegionCache{iip, trustType, floatType} <: @@ -193,6 +199,7 @@ end internalnorm retcode::ReturnCode.T abstol + reltol prob radius_update_scheme::RadiusUpdateSchemes.T trust_r::trustType @@ -220,10 +227,12 @@ end p4::floatType ϵ::floatType stats::NLStats + tc_storage end function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::TrustRegion, args...; - alias_u0 = false, maxiters = 1000, abstol = 1e-8, internalnorm = DEFAULT_NORM, + alias_u0 = false, maxiters = 1000, abstol = nothing, reltol = nothing, + internalnorm = DEFAULT_NORM, linsolve_kwargs = (;), kwargs...) where {uType, iip} @unpack f, u0, p = prob u = alias_u0 ? u0 : deepcopy(u0) @@ -324,13 +333,22 @@ function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::TrustRegion, initial_trust_radius = convert(trustType, 1.0) end + tc = alg.termination_condition + mode = DiffEqBase.get_termination_mode(tc) + + atol = _get_tolerance(abstol, tc.abstol, eltype(u)) + rtol = _get_tolerance(reltol, tc.reltol, eltype(u)) + + storage = mode ∈ DiffEqBase.SAFE_TERMINATION_MODES ? NLSolveSafeTerminationResult() : + nothing + return TrustRegionCache{iip}(f, alg, u_prev, u, fu_prev, fu1, fu2, p, uf, linsolve, J, - jac_cache, false, maxiters, internalnorm, ReturnCode.Default, abstol, prob, + jac_cache, false, maxiters, internalnorm, ReturnCode.Default, atol, rtol, prob, radius_update_scheme, initial_trust_radius, max_trust_radius, step_threshold, shrink_threshold, expand_threshold, shrink_factor, expand_factor, loss, loss_new, H, g, shrink_counter, du, u_tmp, u_gauss_newton, u_cauchy, fu_new, make_new_J, r, p1, p2, p3, p4, ϵ, - NLStats(1, 0, 0, 0, 0)) + NLStats(1, 0, 0, 0, 0), storage) end function perform_step!(cache::TrustRegionCache{true}) @@ -341,7 +359,7 @@ function perform_step!(cache::TrustRegionCache{true}) mul!(cache.g, J', fu) cache.stats.njacs += 1 - # do not use A = cache.H, b = _vec(cache.g) since it is equivalent + # do not use A = cache.H, b = _vec(cache.g) since it is equivalent # to A = cache.J, b = _vec(fu) as long as the Jacobian is non-singular linres = dolinsolve(alg.precs, linsolve, A = J, b = _vec(fu), linu = _vec(u_gauss_newton), @@ -407,6 +425,10 @@ end function trust_region_step!(cache::TrustRegionCache) @unpack fu_new, du, g, H, loss, max_trust_r, radius_update_scheme = cache + + tc_storage = cache.tc_storage + termination_condition = cache.alg.termination_condition(tc_storage) + cache.loss_new = get_loss(fu_new) # Compute the ratio of the actual reduction to the predicted reduction. @@ -435,8 +457,11 @@ function trust_region_step!(cache::TrustRegionCache) # No need to make a new J, no step was taken, so we try again with a smaller trust_r cache.make_new_J = false end - - if iszero(cache.fu) || cache.internalnorm(cache.fu) < cache.abstol + if iszero(cache.fu) || termination_condition(cache.fu, + cache.u, + cache.u_prev, + cache.abstol, + cache.reltol) cache.force_stop = true end @@ -450,12 +475,12 @@ function trust_region_step!(cache::TrustRegionCache) cache.make_new_J = false end - # trust region update - if r < 1 // 10 # cache.shrink_threshold - cache.trust_r *= 1 // 2 # cache.shrink_factor - elseif r >= 9 // 10 # cache.expand_threshold - cache.trust_r = 2 * norm(cache.du) # cache.expand_factor * norm(cache.du) - elseif r >= 1 // 2 # cache.p1 + # trust region update + if r < 1 // 10 # cache.shrink_threshold + cache.trust_r *= 1 // 2 # cache.shrink_factor + elseif r >= 9 // 10 # cache.expand_threshold + cache.trust_r = 2 * norm(cache.du) # cache.expand_factor * norm(cache.du) + elseif r >= 1 // 2 # cache.p1 cache.trust_r = max(cache.trust_r, 2 * norm(cache.du)) # cache.expand_factor * norm(cache.du)) end @@ -504,7 +529,12 @@ function trust_region_step!(cache::TrustRegionCache) cache.trust_r = rfunc(r, shrink_threshold, p1, p3, p4, p2) * cache.internalnorm(du) - if iszero(cache.fu) || cache.internalnorm(cache.fu) < cache.abstol || + if iszero(cache.fu) || + termination_condition(cache.fu, + cache.u, + cache.u_prev, + cache.abstol, + cache.reltol) || cache.internalnorm(g) < cache.ϵ cache.force_stop = true end @@ -529,7 +559,12 @@ function trust_region_step!(cache::TrustRegionCache) @unpack p1 = cache cache.trust_r = p1 * cache.internalnorm(jvp!(cache)) - if iszero(cache.fu) || cache.internalnorm(cache.fu) < cache.abstol || + if iszero(cache.fu) || + termination_condition(cache.fu, + cache.u, + cache.u_prev, + cache.abstol, + cache.reltol) || cache.internalnorm(g) < cache.ϵ cache.force_stop = true end @@ -553,7 +588,12 @@ function trust_region_step!(cache::TrustRegionCache) @unpack p1 = cache cache.trust_r = p1 * (cache.internalnorm(cache.fu)^0.99) - if iszero(cache.fu) || cache.internalnorm(cache.fu) < cache.abstol || + if iszero(cache.fu) || + termination_condition(cache.fu, + cache.u, + cache.u_prev, + cache.abstol, + cache.reltol) || cache.internalnorm(g) < cache.ϵ cache.force_stop = true end @@ -571,7 +611,11 @@ function trust_region_step!(cache::TrustRegionCache) cache.trust_r *= cache.p2 cache.shrink_counter += 1 end - if iszero(cache.fu) || cache.internalnorm(cache.fu) < cache.abstol + if iszero(cache.fu) || termination_condition(cache.fu, + cache.u, + cache.u_prev, + cache.abstol, + cache.reltol) cache.force_stop = true end end diff --git a/test/basictests.jl b/test/basictests.jl index 1e1ded563..5a50b5760 100644 --- a/test/basictests.jl +++ b/test/basictests.jl @@ -123,6 +123,20 @@ end probN = NonlinearProblem(quadratic_f, u0, 2.0) @test all(solve(probN, NewtonRaphson(; autodiff)).u .≈ sqrt(2.0)) end + + @testset "Termination condition: $(mode) u0: $(_nameof(u0))" for mode in instances(NLSolveTerminationMode.T), + u0 in (1.0, [1.0, 1.0]) + + if mode ∈ + (NLSolveTerminationMode.SteadyStateDefault, NLSolveTerminationMode.RelSafeBest, + NLSolveTerminationMode.AbsSafeBest) + continue + end + termination_condition = NLSolveTerminationCondition(mode; abstol = nothing, + reltol = nothing) + probN = NonlinearProblem(quadratic_f, u0, 2.0) + @test all(solve(probN, NewtonRaphson(; termination_condition)).u .≈ sqrt(2.0)) + end end # --- TrustRegion tests --- @@ -280,6 +294,20 @@ end @test sol_iip.u ≈ sol_oop.u end end + + @testset "Termination condition: $(mode) u0: $(_nameof(u0))" for mode in instances(NLSolveTerminationMode.T), + u0 in (1.0, [1.0, 1.0]) + + if mode ∈ + (NLSolveTerminationMode.SteadyStateDefault, NLSolveTerminationMode.RelSafeBest, + NLSolveTerminationMode.AbsSafeBest) + continue + end + termination_condition = NLSolveTerminationCondition(mode; abstol = nothing, + reltol = nothing) + probN = NonlinearProblem(quadratic_f, u0, 2.0) + @test all(solve(probN, TrustRegion(; termination_condition)).u .≈ sqrt(2.0)) + end end # --- LevenbergMarquardt tests --- @@ -389,4 +417,18 @@ end @test all(abs.(quadratic_f(sol.u, 2.0)) .< 1e-10) end end + + @testset "Termination condition: $(mode) u0: $(_nameof(u0))" for mode in instances(NLSolveTerminationMode.T), + u0 in (1.0, [1.0, 1.0]) + + if mode ∈ + (NLSolveTerminationMode.SteadyStateDefault, NLSolveTerminationMode.RelSafeBest, + NLSolveTerminationMode.AbsSafeBest) + continue + end + termination_condition = NLSolveTerminationCondition(mode; abstol = nothing, + reltol = nothing) + probN2 = NonlinearProblem(quadratic_f, u0, 2.0) + @test all(solve(probN, LevenbergMarquardt(; termination_condition)).u .≈ sqrt(2.0)) + end end