Skip to content

Commit 3565824

Browse files
ChrisRackauckasChrisRackauckas
authored
Merge pull request #264 from avik-pal/ap/patch
Don't change the default termination condition
2 parents 8fe131a + c428bd9 commit 3565824

File tree

11 files changed

+95
-161
lines changed

11 files changed

+95
-161
lines changed

‎src/broyden.jl‎

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -185,9 +185,8 @@ function SciMLBase.reinit!(cache::GeneralBroydenCache{iip}, u0 = cache.u; p = ca
185185
cache.u = u0
186186
cache.fu = cache.f(cache.u, p)
187187
end
188-
termination_condition = _get_reinit_termination_condition(cache,
189-
abstol,
190-
reltol,
188+
189+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
191190
termination_condition)
192191

193192
cache.abstol = abstol

‎src/dfsane.jl‎

Lines changed: 9 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -114,23 +114,18 @@ function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg::DFSane, args.
114114
𝒹, uₙ₋₁, fuₙ, fuₙ₋₁ = copy(uₙ), copy(uₙ), copy(uₙ), copy(uₙ)
115115

116116
if iip
117-
# f = (dx, x) -> prob.f(dx, x, p)
118-
# f(fuₙ₋₁, uₙ₋₁)
119117
prob.f(fuₙ₋₁, uₙ₋₁, p)
120118
else
121-
# f = (x) -> prob.f(x, p)
122-
fuₙ₋₁ = prob.f(uₙ₋₁, p) # f(uₙ₋₁)
119+
fuₙ₋₁ = prob.f(uₙ₋₁, p)
123120
end
124121

125122
f₍ₙₒᵣₘ₎ₙ₋₁ = norm(fuₙ₋₁)^nₑₓₚ
126123
f₍ₙₒᵣₘ₎₀ = f₍ₙₒᵣₘ₎ₙ₋₁
127124

128125
ℋ = fill(f₍ₙₒᵣₘ₎ₙ₋₁, M)
129126

130-
abstol, reltol, termination_condition = _init_termination_elements(abstol,
131-
reltol,
132-
termination_condition,
133-
T)
127+
abstol, reltol, termination_condition = _init_termination_elements(abstol, reltol,
128+
termination_condition, T)
134129

135130
mode = DiffEqBase.get_termination_mode(termination_condition)
136131

@@ -167,14 +162,13 @@ function perform_step!(cache::DFSaneCache{true})
167162

168163
f(cache.fuₙ, cache.uₙ)
169164
f₍ₙₒᵣₘ₎ₙ = norm(cache.fuₙ)^nₑₓₚ
170-
for _ in 1:(cache.alg.max_inner_iterations)
165+
for jjj in 1:(cache.alg.max_inner_iterations)
171166
𝒸 = f̄ + η - γ * α₊^2 * f₍ₙₒᵣₘ₎ₙ₋₁
172167

173168
f₍ₙₒᵣₘ₎ₙ ≤ 𝒸 && break
174169

175170
α₊ = clamp(α₊^2 * f₍ₙₒᵣₘ₎ₙ₋₁ / (f₍ₙₒᵣₘ₎ₙ + (T(2) * α₊ - T(1)) * f₍ₙₒᵣₘ₎ₙ₋₁),
176-
τₘᵢₙ * α₊,
177-
τₘₐₓ * α₊)
171+
τₘᵢₙ * α₊, τₘₐₓ * α₊)
178172
@. cache.uₙ = cache.uₙ₋₁ - α₋ * cache.𝒹
179173

180174
f(cache.fuₙ, cache.uₙ)
@@ -183,8 +177,7 @@ function perform_step!(cache::DFSaneCache{true})
183177
f₍ₙₒᵣₘ₎ₙ .≤ 𝒸 && break
184178

185179
α₋ = clamp(α₋^2 * f₍ₙₒᵣₘ₎ₙ₋₁ / (f₍ₙₒᵣₘ₎ₙ + (T(2) * α₋ - T(1)) * f₍ₙₒᵣₘ₎ₙ₋₁),
186-
τₘᵢₙ * α₋,
187-
τₘₐₓ * α₋)
180+
τₘᵢₙ * α₋, τₘₐₓ * α₋)
188181

189182
@. cache.uₙ = cache.uₙ₋₁ + α₊ * cache.𝒹
190183
f(cache.fuₙ, cache.uₙ)
@@ -207,7 +200,7 @@ function perform_step!(cache::DFSaneCache{true})
207200
# Spectral parameter bounds check
208201
if abs(cache.σₙ) > σₘₐₓ || abs(cache.σₙ) < σₘᵢₙ
209202
test_norm = sqrt(sum(abs2, cache.fuₙ₋₁))
210-
cache.σₙ = clamp(1.0 / test_norm, 1, 1e5)
203+
cache.σₙ = clamp(T(1) / test_norm, T(1), T(1e5))
211204
end
212205

213206
# Take step
@@ -283,7 +276,7 @@ function perform_step!(cache::DFSaneCache{false})
283276
# Spectral parameter bounds check
284277
if abs(cache.σₙ) > σₘₐₓ || abs(cache.σₙ) < σₘᵢₙ
285278
test_norm = sqrt(sum(abs2, cache.fuₙ₋₁))
286-
cache.σₙ = clamp(1.0 / test_norm, 1, 1e5)
279+
cache.σₙ = clamp(T(1) / test_norm, T(1), T(1e5))
287280
end
288281

289282
# Take step
@@ -337,9 +330,7 @@ function SciMLBase.reinit!(cache::DFSaneCache{iip}, u0 = cache.uₙ; p = cache.p
337330
T = eltype(cache.uₙ)
338331
cache.σₙ = T(cache.alg.σ_1)
339332

340-
termination_condition = _get_reinit_termination_condition(cache,
341-
abstol,
342-
reltol,
333+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
343334
termination_condition)
344335

345336
cache.abstol = abstol

‎src/gaussnewton.jl‎

Lines changed: 7 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -109,21 +109,17 @@ function SciMLBase.__init(prob::NonlinearLeastSquaresProblem{uType, iip}, alg_::
109109
JᵀJ, Jᵀf = nothing, nothing
110110
end
111111

112-
abstol, reltol, termination_condition = _init_termination_elements(abstol,
113-
reltol,
114-
termination_condition,
115-
eltype(u); mode = NLSolveTerminationMode.AbsNorm)
112+
abstol, reltol, termination_condition = _init_termination_elements(abstol, reltol,
113+
termination_condition, eltype(u); mode = NLSolveTerminationMode.AbsNorm)
116114

117115
mode = DiffEqBase.get_termination_mode(termination_condition)
118116

119117
storage = mode ∈ DiffEqBase.SAFE_TERMINATION_MODES ? NLSolveSafeTerminationResult() :
120118
nothing
121119

122120
return GaussNewtonCache{iip}(f, alg, u, copy(u), fu1, fu2, zero(fu1), du, p, uf,
123-
linsolve, J,
124-
JᵀJ, Jᵀf, jac_cache, false, maxiters, internalnorm, ReturnCode.Default, abstol,
125-
reltol,
126-
prob, NLStats(1, 0, 0, 0, 0), storage, termination_condition)
121+
linsolve, J, JᵀJ, Jᵀf, jac_cache, false, maxiters, internalnorm, ReturnCode.Default,
122+
abstol, reltol, prob, NLStats(1, 0, 0, 0, 0), storage, termination_condition)
127123
end
128124

129125
function perform_step!(cache::GaussNewtonCache{true})
@@ -149,10 +145,7 @@ function perform_step!(cache::GaussNewtonCache{true})
149145
@. u = u - du
150146
f(cache.fu_new, u, p)
151147

152-
(termination_condition(cache.fu_new .- cache.fu1,
153-
cache.u,
154-
u_prev,
155-
cache.abstol,
148+
(termination_condition(cache.fu_new .- cache.fu1, cache.u, u_prev, cache.abstol,
156149
cache.reltol) ||
157150
termination_condition(cache.fu_new, cache.u, u_prev, cache.abstol, cache.reltol)) &&
158151
(cache.force_stop = true)
@@ -219,9 +212,8 @@ function SciMLBase.reinit!(cache::GaussNewtonCache{iip}, u0 = cache.u; p = cache
219212
cache.u = u0
220213
cache.fu1 = cache.f(cache.u, p)
221214
end
222-
termination_condition = _get_reinit_termination_condition(cache,
223-
abstol,
224-
reltol,
215+
216+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
225217
termination_condition)
226218

227219
cache.abstol = abstol

‎src/klement.jl‎

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -238,10 +238,9 @@ function SciMLBase.reinit!(cache::GeneralKlementCache{iip}, u0 = cache.u; p = ca
238238
cache.fu = cache.f(cache.u, p)
239239
end
240240

241-
termination_condition = _get_reinit_termination_condition(cache,
242-
abstol,
243-
reltol,
241+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
244242
termination_condition)
243+
245244
cache.abstol = abstol
246245
cache.reltol = reltol
247246
cache.termination_condition = termination_condition

‎src/levenberg.jl‎

Lines changed: 2 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -185,10 +185,8 @@ function SciMLBase.__init(prob::Union{NonlinearProblem{uType, iip},
185185
v = similar(du)
186186
end
187187

188-
abstol, reltol, termination_condition = _init_termination_elements(abstol,
189-
reltol,
190-
termination_condition,
191-
eltype(u); mode = NLSolveTerminationMode.AbsNorm)
188+
abstol, reltol, termination_condition = _init_termination_elements(abstol, reltol,
189+
termination_condition, eltype(u); mode = NLSolveTerminationMode.AbsNorm)
192190

193191
λ = convert(eltype(u), alg.damping_initial)
194192
λ_factor = convert(eltype(u), alg.damping_increase_factor)

‎src/pseudotransient.jl‎

Lines changed: 13 additions & 35 deletions
Original file line numberDiff line numberDiff line change
@@ -2,10 +2,13 @@
22
PseudoTransient(; concrete_jac = nothing, linsolve = nothing,
33
precs = DEFAULT_PRECS, alpha_initial = 1e-3, adkwargs...)
44
5-
An implementation of PseudoTransient method that is used to solve steady state problems in an accelerated manner. It uses an adaptive time-stepping to
6-
integrate an initial value of nonlinear problem until sufficient accuracy in the desired steady-state is achieved to switch over to Newton's method and
7-
gain a rapid convergence. This implementation specifically uses "switched evolution relaxation" SER method. For detail information about the time-stepping and algorithm,
8-
please see the paper: [Coffey, Todd S. and Kelley, C. T. and Keyes, David E. (2003), Pseudotransient Continuation and Differential-Algebraic Equations,
5+
An implementation of PseudoTransient method that is used to solve steady state problems in
6+
an accelerated manner. It uses an adaptive time-stepping to integrate an initial value of
7+
nonlinear problem until sufficient accuracy in the desired steady-state is achieved to
8+
switch over to Newton's method and gain a rapid convergence. This implementation
9+
specifically uses "switched evolution relaxation" SER method. For detail information about
10+
the time-stepping and algorithm, please see the paper:
11+
[Coffey, Todd S. and Kelley, C. T. and Keyes, David E. (2003), Pseudotransient Continuation and Differential-Algebraic Equations,
912
SIAM Journal on Scientific Computing,25, 553-569.](https://doi.org/10.1137/S106482750241044X)
1013
1114
### Keyword Arguments
@@ -48,7 +51,7 @@ function PseudoTransient(; concrete_jac = nothing, linsolve = nothing,
4851
return PseudoTransient{_unwrap_val(concrete_jac)}(ad, linsolve, precs, alpha_initial)
4952
end
5053

51-
@concrete mutable struct PseudoTransientCache{iip}
54+
@concrete mutable struct PseudoTransientCache{iip} <: AbstractNonlinearSolveCache{iip}
5255
f
5356
alg
5457
u
@@ -75,14 +78,10 @@ end
7578
tc_storage
7679
end
7780

78-
isinplace(::PseudoTransientCache{iip}) where {iip} = iip
79-
8081
function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg_::PseudoTransient,
81-
args...;
82-
alias_u0 = false, maxiters = 1000, abstol = nothing, reltol = nothing,
82+
args...; alias_u0 = false, maxiters = 1000, abstol = nothing, reltol = nothing,
8383
termination_condition = nothing, internalnorm = DEFAULT_NORM,
84-
linsolve_kwargs = (;),
85-
kwargs...) where {uType, iip}
84+
linsolve_kwargs = (;), kwargs...) where {uType, iip}
8685
alg = get_concrete_algorithm(alg_, prob)
8786

8887
@unpack f, u0, p = prob
@@ -99,9 +98,7 @@ function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg_::PseudoTransi
9998
res_norm = internalnorm(fu1)
10099

101100
abstol, reltol, termination_condition = _init_termination_elements(abstol,
102-
reltol,
103-
termination_condition,
104-
eltype(u))
101+
reltol, termination_condition, eltype(u))
105102

106103
mode = DiffEqBase.get_termination_mode(termination_condition)
107104

@@ -111,8 +108,7 @@ function SciMLBase.__init(prob::NonlinearProblem{uType, iip}, alg_::PseudoTransi
111108
return PseudoTransientCache{iip}(f, alg, u, copy(u), fu1, fu2, du, p, alpha, res_norm,
112109
uf,
113110
linsolve, J, jac_cache, false, maxiters, internalnorm, ReturnCode.Default, abstol,
114-
reltol,
115-
prob, NLStats(1, 0, 0, 0, 0), termination_condition, storage)
111+
reltol, prob, NLStats(1, 0, 0, 0, 0), termination_condition, storage)
116112
end
117113

118114
function perform_step!(cache::PseudoTransientCache{true})
@@ -176,22 +172,6 @@ function perform_step!(cache::PseudoTransientCache{false})
176172
return nothing
177173
end
178174

179-
function SciMLBase.solve!(cache::PseudoTransientCache)
180-
while !cache.force_stop && cache.stats.nsteps < cache.maxiters
181-
perform_step!(cache)
182-
cache.stats.nsteps += 1
183-
end
184-
185-
if cache.stats.nsteps == cache.maxiters
186-
cache.retcode = ReturnCode.MaxIters
187-
else
188-
cache.retcode = ReturnCode.Success
189-
end
190-
191-
return SciMLBase.build_solution(cache.prob, cache.alg, cache.u, cache.fu1;
192-
cache.retcode, cache.stats)
193-
end
194-
195175
function SciMLBase.reinit!(cache::PseudoTransientCache{iip}, u0 = cache.u; p = cache.p,
196176
alpha_new,
197177
abstol = cache.abstol, reltol = cache.reltol,
@@ -207,9 +187,7 @@ function SciMLBase.reinit!(cache::PseudoTransientCache{iip}, u0 = cache.u; p = c
207187
cache.fu1 = cache.f(cache.u, p)
208188
end
209189

210-
termination_condition = _get_reinit_termination_condition(cache,
211-
abstol,
212-
reltol,
190+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
213191
termination_condition)
214192

215193
cache.alpha = convert(eltype(cache.u), alpha_new)

‎src/raphson.jl‎

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -179,10 +179,9 @@ function SciMLBase.reinit!(cache::NewtonRaphsonCache{iip}, u0 = cache.u; p = cac
179179
cache.fu1 = cache.f(cache.u, p)
180180
end
181181

182-
termination_condition = _get_reinit_termination_condition(cache,
183-
abstol,
184-
reltol,
182+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
185183
termination_condition)
184+
186185
cache.abstol = abstol
187186
cache.reltol = reltol
188187
cache.termination_condition = termination_condition

‎src/trustRegion.jl‎

Lines changed: 2 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -736,9 +736,8 @@ function SciMLBase.reinit!(cache::TrustRegionCache{iip}, u0 = cache.u; p = cache
736736
cache.u = u0
737737
cache.fu = cache.f(cache.u, p)
738738
end
739-
termination_condition = _get_reinit_termination_condition(cache,
740-
abstol,
741-
reltol,
739+
740+
termination_condition = _get_reinit_termination_condition(cache, abstol, reltol,
742741
termination_condition)
743742

744743
cache.abstol = abstol

‎src/utils.jl‎

Lines changed: 22 additions & 24 deletions
Original file line numberDiff line numberDiff line change
@@ -215,46 +215,44 @@ function _get_tolerance(η, tc_η, ::Type{T}) where {T}
215215
return T(ifelse(η !== nothing, η, ifelse(tc_η !== nothing, tc_η, fallback_η)))
216216
end
217217

218-
function _init_termination_elements(abstol,
219-
reltol,
220-
termination_condition,
221-
::Type{T}; mode = NLSolveTerminationMode.NLSolveDefault) where {T}
218+
function _init_termination_elements(abstol, reltol, termination_condition,
219+
::Type{T}; mode = NLSolveTerminationMode.AbsNorm) where {T}
222220
if termination_condition !== nothing
223-
abstol !== nothing ?
224-
(abstol != termination_condition.abstol ?
225-
error("Incompatible absolute tolerances found. The tolerances supplied as the keyword argument and the one supplied in the termination condition should be same.") :
226-
nothing) : nothing
227-
reltol !== nothing ?
228-
(reltol != termination_condition.abstol ?
229-
error("Incompatible relative tolerances found. The tolerances supplied as the keyword argument and the one supplied in the termination condition should be same.") :
230-
nothing) : nothing
221+
if abstol !== nothing && abstol != termination_condition.abstol
222+
error("Incompatible absolute tolerances found. The tolerances supplied as the \
223+
keyword argument and the one supplied in the termination condition should \
224+
be same.")
225+
end
226+
if reltol !== nothing && reltol != termination_condition.reltol
227+
error("Incompatible relative tolerances found. The tolerances supplied as the \
228+
keyword argument and the one supplied in the termination condition should \
229+
be same.")
230+
end
231231
abstol = _get_tolerance(abstol, termination_condition.abstol, T)
232232
reltol = _get_tolerance(reltol, termination_condition.reltol, T)
233233
return abstol, reltol, termination_condition
234234
else
235235
abstol = _get_tolerance(abstol, nothing, T)
236236
reltol = _get_tolerance(reltol, nothing, T)
237-
termination_condition = NLSolveTerminationCondition(mode;
238-
abstol,
239-
reltol)
237+
termination_condition = NLSolveTerminationCondition(mode; abstol, reltol)
240238
return abstol, reltol, termination_condition
241239
end
242240
end
243241

244242
function _get_reinit_termination_condition(cache, abstol, reltol, termination_condition)
245243
if termination_condition != cache.termination_condition
246-
if abstol != cache.abstol
247-
if abstol != termination_condition.abstol
248-
error("Incompatible absolute tolerances found. The tolerances supplied as the keyword argument and the one supplied in the termination condition should be same.")
249-
end
244+
if abstol != cache.abstol && abstol != termination_condition.abstol
245+
error("Incompatible absolute tolerances found. The tolerances supplied as the \
246+
keyword argument and the one supplied in the termination condition \
247+
should be same.")
250248
end
251249

252-
if reltol != cache.reltol
253-
if reltol != termination_condition.reltol
254-
error("Incompatible absolute tolerances found. The tolerances supplied as the keyword argument and the one supplied in the termination condition should be same.")
255-
end
250+
if reltol != cache.reltol && reltol != termination_condition.reltol
251+
error("Incompatible absolute tolerances found. The tolerances supplied as the \
252+
keyword argument and the one supplied in the termination condition \
253+
should be same.")
256254
end
257-
termination_condition
255+
return termination_condition
258256
else
259257
# Build the termination_condition with new abstol and reltol
260258
return NLSolveTerminationCondition{

‎test/23_test_problems.jl‎

Lines changed: 2 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -11,8 +11,9 @@ function test_on_library(problems, dicts, alg_ops, broken_tests, ϵ = 1e-4)
1111
@testset "$idx: $(dict["title"])" begin
1212
for alg in alg_ops
1313
try
14-
sol = solve(nlprob, alg, abstol = 1e-18, reltol = 1e-18)
14+
sol = solve(nlprob, alg)
1515
problem(res, sol.u, nothing)
16+
1617
broken = idx in broken_tests[alg] ? true : false
1718
@test norm(res)≤ϵ broken=broken
1819
catch

0 commit comments

Comments
 (0)