diff --git a/docs/src/stiff/rosenbrock.md b/docs/src/stiff/rosenbrock.md index d592d93e67..a72026da5a 100644 --- a/docs/src/stiff/rosenbrock.md +++ b/docs/src/stiff/rosenbrock.md @@ -3,6 +3,13 @@ ## Standard Rosenbrock Methods ```@docs +ROS2 +ROS3 +ROS2PR +ROS3PR +Scholz47 +ROS3PRL +ROS3PRL2 ROS3P Rodas3 Rodas3P @@ -31,5 +38,7 @@ ROS34PW1a ROS34PW1b ROS34PW2 ROS34PW3 +ROS34PRw +ROS2S RosenbrockW6S4OS ``` diff --git a/src/OrdinaryDiffEq.jl b/src/OrdinaryDiffEq.jl index 6973d2b75a..eb4cc169b9 100644 --- a/src/OrdinaryDiffEq.jl +++ b/src/OrdinaryDiffEq.jl @@ -392,7 +392,7 @@ export Rosenbrock23, Rosenbrock32, RosShamp4, Veldd4, Velds4, GRK4T, GRK4A, Rodas5, Rodas5P, Rodas5Pe, Rodas5Pr, RosenbrockW6S4OS, ROS34PW1a, ROS34PW1b, ROS34PW2, ROS34PW3, ROS34PRw, ROS3PRL, ROS3PRL2, - ROS2PR, ROS2S, ROS3, ROS3PR, Scholz4_7 + ROS2, ROS2PR, ROS2S, ROS3, ROS3PR, Scholz4_7 export LawsonEuler, NorsettEuler, ETD1, ETDRK2, ETDRK3, ETDRK4, HochOst4, Exp4, EPIRK4s3A, EPIRK4s3B, diff --git a/src/alg_utils.jl b/src/alg_utils.jl index 6c93dbc671..62ad86e2af 100644 --- a/src/alg_utils.jl +++ b/src/alg_utils.jl @@ -634,6 +634,7 @@ alg_order(alg::Feagin12) = 12 alg_order(alg::Feagin14) = 14 alg_order(alg::PFRK87) = 8 +alg_order(alg::ROS2) = 2 alg_order(alg::ROS2PR) = 2 alg_order(alg::ROS2S) = 2 alg_order(alg::ROS3) = 3 diff --git a/src/algorithms.jl b/src/algorithms.jl index fc367c6427..b8ff73dee7 100644 --- a/src/algorithms.jl +++ b/src/algorithms.jl @@ -2898,6 +2898,11 @@ end - Shampine L.F. and Reichelt M., (1997) The MATLAB ODE Suite, SIAM Journal of Scientific Computing, 18 (1), pp. 1-22. +#### ROS2 + +- J. G. Verwer et al. (1999): A second-order Rosenbrock method applied to photochemical dispersion problems + https://doi.org/10.1137/S1064827597326651 + #### ROS3P - Lang, J. & Verwer, ROS3Pβ€”An Accurate Third-Order Rosenbrock Solver Designed for @@ -2972,6 +2977,7 @@ University of Geneva, Switzerland. for Alg in [ :Rosenbrock23, :Rosenbrock32, + :ROS2, :ROS2PR, :ROS2S, :ROS3, diff --git a/src/caches/rosenbrock_caches.jl b/src/caches/rosenbrock_caches.jl index ec38dc6739..7aeac95ebe 100644 --- a/src/caches/rosenbrock_caches.jl +++ b/src/caches/rosenbrock_caches.jl @@ -407,6 +407,12 @@ end ################################################################################ +### ROS2 methods + +@ROS2(:cache) + +################################################################################ + ### ROS23 methods @ROS23(:cache) diff --git a/src/generic_rosenbrock.jl b/src/generic_rosenbrock.jl index a484d16f66..44ecf3c2ae 100644 --- a/src/generic_rosenbrock.jl +++ b/src/generic_rosenbrock.jl @@ -805,7 +805,7 @@ macro Rosenbrock4(part) end end -#ROS23 and ROS34PW methods (Rang and Angermann, 2005) +#ROS2, ROS23 and ROS34PW methods (Rang and Angermann, 2005) """ Ros34dummyTableau() @@ -849,6 +849,25 @@ function Ros23dummyTableau() RosenbrockAdaptiveTableau(a,C,b,btilde,gamma,d,c) end +""" + Ros2dummyTableau() + +Generate a dummy tableau for ROS2 methods. This type of methods has 2 steps. +""" +function Ros2dummyTableau() + a=[false false; + true false] + C=[false false; + true false] + b=[true,true] + btilde=[true,true] + gamma=true + c=[false,true] + d=[true,true] + RosenbrockAdaptiveTableau(a,C,b,btilde,gamma,d,c) +end + + """ _transformtab(Alpha,Gamma,B,Bhat) @@ -872,6 +891,78 @@ end +# 2 step ROS Methods +""" + ROS2Tableau() +2nd order stiffly accurate Rosenbrock-Wanner method with 2 internal stages with (Rinf=0). +The embedded method is taken from Kinetic PreProcessor (KPP). +J. G. Verwer et al. (1999): A second-order Rosenbrock method applied to photochemical dispersion problems +https://doi.org/10.1137/S1064827597326651 +""" +function ROS2Tableau() # 2nd order + gamma=1.7071067811865475 # 1+1/sqrt(2) + Alpha=[0 0; + 1. 0] + Gamma=[gamma 0; + -3.414213562373095 gamma] + B=[0.5, 0.5] + Bhat=[1, 0] + a,C,b,btilde,d,c=_transformtab(Alpha,Gamma,B,Bhat) + RosenbrockAdaptiveTableau(a,C,b,btilde,gamma,d,c) +end + +@doc """ + ROS2() + +2nd order stiffly accurate Rosenbrock-Wanner method with 2 internal stages with (Rinf=0). +The embedded method is taken from Kinetic PreProcessor (KPP). +J. G. Verwer et al. (1999): A second-order Rosenbrock method applied to photochemical dispersion problems +More Information at https://doi.org/10.1137/S1064827597326651 +""" ROS2 + + + +""" + @ROS2(part) + +Generate code for the 2 step ROS methods: ROS2 +`part` should be one of `:tableau`, `:cache`, `:init`, `:performstep`. +`@ROS2(:tableau)` should be placed in `tableaus/rosenbrock_tableaus.jl`. +`@ROS2(:cache)` should be placed in `caches/rosenbrock_caches.jl`. +`@ROS2(:init)` and `@ROS2(:performstep)` should be placed in +`perform_step/rosenbrock_perform_step.jl`. +""" +macro ROS2(part) + tabmask=Ros2dummyTableau() + cachename=:ROS2Cache + constcachename=:ROS2ConstantCache + ROS2tabname=:ROS2Tableau + n_normalstep=length(tabmask.b)-1 + if part.value==:tableau + tabstructexpr=gen_tableau_struct(tabmask,:Ros2Tableau) + tabexprs=Array{Expr,1}([tabstructexpr]) + push!(tabexprs,gen_tableau(ROS2Tableau(),tabstructexpr,ROS2tabname)) + return esc(quote $(tabexprs...) end) + elseif part.value==:cache + constcacheexpr,cacheexpr=gen_cache_struct(tabmask,cachename,constcachename) + cacheexprs=Array{Expr,1}([constcacheexpr,cacheexpr]) + push!(cacheexprs,gen_algcache(cacheexpr,constcachename,:ROS2,ROS2tabname)) + return esc(quote $(cacheexprs...) end) + elseif part.value==:init + return esc(gen_initialize(cachename,constcachename)) + elseif part.value==:performstep + performstepexprs=Array{Expr,1}() + push!(performstepexprs,gen_constant_perform_step(tabmask,constcachename,n_normalstep)) + push!(performstepexprs,gen_perform_step(tabmask,cachename,n_normalstep)) + return esc(quote $(performstepexprs...) end) + else + throw(ArgumentError("Unknown parameter!")) + nothing + end +end + + + # 3 step ROS Methods """ ROS2PRTableau() diff --git a/src/perform_step/rosenbrock_perform_step.jl b/src/perform_step/rosenbrock_perform_step.jl index b803e0612d..4937407918 100644 --- a/src/perform_step/rosenbrock_perform_step.jl +++ b/src/perform_step/rosenbrock_perform_step.jl @@ -933,6 +933,13 @@ end ################################################################################ +#### ROS2 type method + +@ROS2(:init) +@ROS2(:performstep) + +################################################################################ + #### ROS23 type method @ROS23(:init) diff --git a/src/tableaus/rosenbrock_tableaus.jl b/src/tableaus/rosenbrock_tableaus.jl index 616799d1a8..493c67111a 100644 --- a/src/tableaus/rosenbrock_tableaus.jl +++ b/src/tableaus/rosenbrock_tableaus.jl @@ -214,6 +214,8 @@ function Rodas3PTableau(T, T2) h21, h22, h23, h24, h25, h31, h32, h33, h34, h35, h2_21, h2_22, h2_23, h2_24, h2_25) end +@ROS2(:tableau) + @ROS23(:tableau) @ROS34PW(:tableau) diff --git a/test/algconvergence/ode_rosenbrock_tests.jl b/test/algconvergence/ode_rosenbrock_tests.jl index 6c4841fa4c..257d424af9 100644 --- a/test/algconvergence/ode_rosenbrock_tests.jl +++ b/test/algconvergence/ode_rosenbrock_tests.jl @@ -90,6 +90,23 @@ import LinearSolve sol = solve(prob, Rodas3()) @test length(sol) < 20 + ### ROS2 + prob = prob_ode_linear + + sim = test_convergence(dts, prob, ROS2()) + @test sim.π’ͺest[:final]β‰ˆ2 atol=testTol + + sol = solve(prob, ROS2()) + @test length(sol) < 61 + + prob = prob_ode_2Dlinear + + sim = test_convergence(dts, prob, ROS2()) + @test sim.π’ͺest[:final]β‰ˆ2 atol=testTol + + sol = solve(prob, ROS2PR()) + @test length(sol) < 60 + ### ROS2PR prob = prob_ode_linear diff --git a/test/interface/mass_matrix_tests.jl b/test/interface/mass_matrix_tests.jl index 08d7eb2c4d..ff177d5f96 100644 --- a/test/interface/mass_matrix_tests.jl +++ b/test/interface/mass_matrix_tests.jl @@ -107,6 +107,7 @@ dependent_M2 = MatrixOperator(ones(3, 3), update_func = update_func2, @test _norm_dsol(Rosenbrock32(), prob, prob2)β‰ˆ0 atol=1e-11 @test _norm_dsol(ROS3P(), prob, prob2)β‰ˆ0 atol=1e-11 @test _norm_dsol(Rodas3(), prob, prob2)β‰ˆ0 atol=1e-11 + @test _norm_dsol(ROS2(), prob, prob2)β‰ˆ0 atol=1e-11 @test _norm_dsol(ROS2PR(), prob, prob2)β‰ˆ0 atol=1e-11 @test _norm_dsol(ROS2S(), prob, prob2)β‰ˆ0 atol=1e-11 @test _norm_dsol(ROS3(), prob, prob2)β‰ˆ0 atol=1e-11 diff --git a/test/regression/iipvsoop_tests.jl b/test/regression/iipvsoop_tests.jl index c978289c9a..df0fb77696 100644 --- a/test/regression/iipvsoop_tests.jl +++ b/test/regression/iipvsoop_tests.jl @@ -116,7 +116,10 @@ end working_rosenbrock_algs = [Rosenbrock23(), ROS3P(), Rodas3(), RosShamp4(), Veldd4(), Velds4(), GRK4T(), GRK4A(), Ros4LStab(), Rodas4(), Rodas42(), Rodas4P(), Rodas5(), - Rodas23W(), Rodas3P(), Rodas5Pe(), Rodas5P()] + Rodas23W(), Rodas3P(), Rodas5Pe(), Rodas5P(), + ROS2(), ROS2PR(), ROS2S(), ROS3(), ROS3PR(), Scholz4_7(), + ROS34PW1a(), ROS34PW1b(), ROS34PW2(), ROS34PW3(), + ROS34PRw(), ROS3PRL(), ROS3PRL2()] rosenbrock_algs = [Rosenbrock32() ]