Skip to content

Commit 9df099e

Browse files
author
oscarddssmith
committed
prepare for switching to Linsolve Interface
1 parent 44014a7 commit 9df099e

File tree

5 files changed

+39
-146
lines changed

5 files changed

+39
-146
lines changed

‎lib/OrdinaryDiffEqNonlinearSolve/src/newton.jl‎

Lines changed: 3 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -220,14 +220,9 @@ end
220220
reltol = eps(eltype(dz))
221221
end
222222

223-
if is_always_new(nlsolver) || (iter == 1 && new_W)
224-
linres = dolinsolve(integrator, linsolve; A = W, b = _vec(b), linu = _vec(dz),
225-
reltol = reltol)
226-
else
227-
linres = dolinsolve(
228-
integrator, linsolve; A = nothing, b = _vec(b), linu = _vec(dz),
229-
reltol = reltol)
230-
end
223+
make_new_W = is_always_new(nlsolver) || (iter == 1 && new_W)
224+
linres = dolinsolve(integrator, linsolve; A = make_new_W ? W : nothing, b = _vec(b),
225+
linu = _vec(dz), reltol)
231226

232227
if !SciMLBase.successful_retcode(linres.retcode) && linres.retcode != SciMLBase.ReturnCode.Default
233228
return convert(eltype(atmp,),Inf)

‎lib/OrdinaryDiffEqNonlinearSolve/src/utils.jl‎

Lines changed: 4 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -191,14 +191,10 @@ function build_nlsolver(
191191
end
192192
jac_config = build_jac_config(alg, nf, uf, du1, uprev, u, ztmp, dz)
193193
end
194-
linprob = LinearProblem(W, _vec(k); u0 = _vec(dz))
195-
Pl, Pr = wrapprecs(
196-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
197-
nothing)...,
198-
weight, dz)
199-
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
200-
Pl = Pl, Pr = Pr,
201-
assumptions = LinearSolve.OperatorAssumptions(true))
194+
linprob = LinearProblem(W, _vec(k), (isdae ? du1 : nothing,u,p,t); u0 = _vec(dz))
195+
linsolve = init(linprob,
196+
wrapprecs(alg.linsolve, W, weight),
197+
(isdae ? du1 : nothing,u,p,t); alias_A = true, alias_b = true)
202198

203199
tType = typeof(t)
204200
invγdt = inv(oneunit(t) * one(uTolType))

‎lib/OrdinaryDiffEqRosenbrock/src/generic_rosenbrock.jl‎

Lines changed: 6 additions & 8 deletions
Original file line numberDiff line numberDiff line change
@@ -247,10 +247,8 @@ function gen_algcache(cacheexpr::Expr,constcachename::Symbol,algname::Symbol,tab
247247
tf = TimeGradientWrapper(f,uprev,p)
248248
uf = UJacobianWrapper(f,t,p)
249249
linsolve_tmp = zero(rate_prototype)
250-
linprob = LinearProblem(W,_vec(linsolve_tmp); u0=_vec(tmp))
251-
linsolve = init(linprob,alg.linsolve,alias_A=true,alias_b=true,
252-
Pl = LinearSolve.InvPreconditioner(Diagonal(_vec(weight))),
253-
Pr = Diagonal(_vec(weight)))
250+
linprob = LinearProblem(W,_vec(linsolve_tmp), (nothing, u, p, t); u0=_vec(tmp))
251+
linsolve = init(linprob,alg.linsolve,alias_A=true,alias_b=true)
254252
grad_config = build_grad_config(alg,f,tf,du1,t)
255253
jac_config = build_jac_config(alg,f,uf,du1,uprev,u,tmp,du2)
256254
$cachename($(valsyms...))
@@ -1036,7 +1034,7 @@ references = """
10361034
""",
10371035
"Rodas3",
10381036
references = """
1039-
- Sandu, Verwer, Van Loon, Carmichael, Potra, Dabdub, Seinfeld, Benchmarking stiff ode solvers for atmospheric chemistry problems-I.
1037+
- Sandu, Verwer, Van Loon, Carmichael, Potra, Dabdub, Seinfeld, Benchmarking stiff ode solvers for atmospheric chemistry problems-I.
10401038
implicit vs explicit, Atmospheric Environment, 31(19), 3151-3166, 1997.
10411039
""",
10421040
with_step_limiter=true) Rodas3
@@ -1096,9 +1094,9 @@ lower if not corrected).
10961094
""",
10971095
"Rodas4P",
10981096
references = """
1099-
- Steinebach, G., Rentrop, P., An adaptive method of lines approach for modelling flow and transport in rivers.
1097+
- Steinebach, G., Rentrop, P., An adaptive method of lines approach for modelling flow and transport in rivers.
11001098
Adaptive method of lines , Wouver, A. Vande, Sauces, Ph., Schiesser, W.E. (ed.),S. 181-205,Chapman & Hall/CRC, 2001,
1101-
- Steinebach, G., Oder-reduction of ROW-methods for DAEs and method of lines applications.
1099+
- Steinebach, G., Oder-reduction of ROW-methods for DAEs and method of lines applications.
11021100
Preprint-Nr. 1741, FB Mathematik, TH Darmstadt, 1995.
11031101
""",
11041102
with_step_limiter=true) Rodas4P
@@ -1111,7 +1109,7 @@ of Roadas4P and in case of inexact Jacobians a second order W method.
11111109
""",
11121110
"Rodas4P2",
11131111
references = """
1114-
- Steinebach G., Improvement of Rosenbrock-Wanner Method RODASP, In: Reis T., Grundel S., Schöps S. (eds)
1112+
- Steinebach G., Improvement of Rosenbrock-Wanner Method RODASP, In: Reis T., Grundel S., Schöps S. (eds)
11151113
Progress in Differential-Algebraic Equations II. Differential-Algebraic Equations Forum. Springer, Cham., 165-184, 2020.
11161114
""",
11171115
with_step_limiter=true) Rodas4P2

‎lib/OrdinaryDiffEqRosenbrock/src/rosenbrock_caches.jl‎

Lines changed: 7 additions & 38 deletions
Original file line numberDiff line numberDiff line change
@@ -150,12 +150,8 @@ function alg_cache(alg::Rosenbrock23, u, rate_prototype, ::Type{uEltypeNoUnits},
150150
uf = UJacobianWrapper(f, t, p)
151151
linsolve_tmp = zero(rate_prototype)
152152

153-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
154-
Pl, Pr = wrapprecs(
155-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
156-
nothing)..., weight, tmp)
153+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing,u,p,t); u0 = _vec(tmp))
157154
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
158-
Pl = Pl, Pr = Pr,
159155
assumptions = LinearSolve.OperatorAssumptions(true))
160156

161157
grad_config = build_grad_config(alg, f, tf, du1, t)
@@ -195,13 +191,8 @@ function alg_cache(alg::Rosenbrock32, u, rate_prototype, ::Type{uEltypeNoUnits},
195191
tf = TimeGradientWrapper(f, uprev, p)
196192
uf = UJacobianWrapper(f, t, p)
197193
linsolve_tmp = zero(rate_prototype)
198-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
199-
200-
Pl, Pr = wrapprecs(
201-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
202-
nothing)..., weight, tmp)
194+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing,u,p,t); u0 = _vec(tmp))
203195
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
204-
Pl = Pl, Pr = Pr,
205196
assumptions = LinearSolve.OperatorAssumptions(true))
206197
grad_config = build_grad_config(alg, f, tf, du1, t)
207198
jac_config = build_jac_config(alg, f, uf, du1, uprev, u, tmp, du2)
@@ -344,12 +335,8 @@ function alg_cache(alg::ROS3P, u, rate_prototype, ::Type{uEltypeNoUnits},
344335
tf = TimeGradientWrapper(f, uprev, p)
345336
uf = UJacobianWrapper(f, t, p)
346337
linsolve_tmp = zero(rate_prototype)
347-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
348-
Pl, Pr = wrapprecs(
349-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
350-
nothing)..., weight, tmp)
338+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing, u, p, t); u0 = _vec(tmp))
351339
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
352-
Pl = Pl, Pr = Pr,
353340
assumptions = LinearSolve.OperatorAssumptions(true))
354341
grad_config = build_grad_config(alg, f, tf, du1, t)
355342
jac_config = build_jac_config(alg, f, uf, du1, uprev, u, tmp, du2)
@@ -430,12 +417,8 @@ function alg_cache(alg::Rodas3, u, rate_prototype, ::Type{uEltypeNoUnits},
430417
tf = TimeGradientWrapper(f, uprev, p)
431418
uf = UJacobianWrapper(f, t, p)
432419
linsolve_tmp = zero(rate_prototype)
433-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
434-
Pl, Pr = wrapprecs(
435-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
436-
nothing)..., weight, tmp)
420+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing, u, p, t); u0 = _vec(tmp))
437421
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
438-
Pl = Pl, Pr = Pr,
439422
assumptions = LinearSolve.OperatorAssumptions(true))
440423
grad_config = build_grad_config(alg, f, tf, du1, t)
441424
jac_config = build_jac_config(alg, f, uf, du1, uprev, u, tmp, du2)
@@ -623,12 +606,8 @@ function alg_cache(alg::Rodas23W, u, rate_prototype, ::Type{uEltypeNoUnits},
623606
tf = TimeGradientWrapper(f, uprev, p)
624607
uf = UJacobianWrapper(f, t, p)
625608
linsolve_tmp = zero(rate_prototype)
626-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
627-
Pl, Pr = wrapprecs(
628-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
629-
nothing)..., weight, tmp)
609+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing, u, p, t); u0 = _vec(tmp))
630610
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
631-
Pl = Pl, Pr = Pr,
632611
assumptions = LinearSolve.OperatorAssumptions(true))
633612
grad_config = build_grad_config(alg, f, tf, du1, t)
634613
jac_config = build_jac_config(alg, f, uf, du1, uprev, u, tmp, du2)
@@ -667,12 +646,8 @@ function alg_cache(alg::Rodas3P, u, rate_prototype, ::Type{uEltypeNoUnits},
667646
tf = TimeGradientWrapper(f, uprev, p)
668647
uf = UJacobianWrapper(f, t, p)
669648
linsolve_tmp = zero(rate_prototype)
670-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
671-
Pl, Pr = wrapprecs(
672-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
673-
nothing)..., weight, tmp)
649+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing, u, p, t); u0 = _vec(tmp))
674650
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
675-
Pl = Pl, Pr = Pr,
676651
assumptions = LinearSolve.OperatorAssumptions(true))
677652
grad_config = build_grad_config(alg, f, tf, du1, t)
678653
jac_config = build_jac_config(alg, f, uf, du1, uprev, u, tmp, du2)
@@ -771,14 +746,8 @@ function alg_cache(alg::Union{Rodas4, Rodas42, Rodas4P, Rodas4P2, Rodas5, Rodas5
771746
tf = TimeGradientWrapper(f, uprev, p)
772747
uf = UJacobianWrapper(f, t, p)
773748
linsolve_tmp = zero(rate_prototype)
774-
linprob = LinearProblem(W, _vec(linsolve_tmp); u0 = _vec(tmp))
775-
776-
Pl, Pr = wrapprecs(
777-
alg.precs(W, nothing, u, p, t, nothing, nothing, nothing,
778-
nothing)..., weight, tmp)
779-
749+
linprob = LinearProblem(W, _vec(linsolve_tmp), (nothing, u, p, t); u0 = _vec(tmp))
780750
linsolve = init(linprob, alg.linsolve, alias_A = true, alias_b = true,
781-
Pl = Pl, Pr = Pr,
782751
assumptions = LinearSolve.OperatorAssumptions(true))
783752

784753
grad_config = build_grad_config(alg, f, tf, du1, t)

0 commit comments

Comments
 (0)