diff --git a/src/OrdinaryDiffEq.jl b/src/OrdinaryDiffEq.jl index 4521e6a117..b7ba5f8021 100644 --- a/src/OrdinaryDiffEq.jl +++ b/src/OrdinaryDiffEq.jl @@ -176,7 +176,7 @@ module OrdinaryDiffEq export FunctionMap, Euler, Heun, Ralston, Midpoint, RK4, ExplicitRK, OwrenZen3, OwrenZen4, OwrenZen5, BS3, BS5, RK46NL, DP5, Tsit5, DP8, Vern6, Vern7, Vern8, TanYam7, TsitPap8, - Vern9, Feagin10, Feagin12, Feagin14, CompositeAlgorithm, Anas5, RKO65, FRK65 + Vern9, Feagin10, Feagin12, Feagin14, CompositeAlgorithm, Anas5, RKO65, FRK65, PFRK87 export SSPRK22, SSPRK33, SSPRK53, SSPRK53_2N1, SSPRK53_2N2, SSPRK63, SSPRK73, SSPRK83, SSPRK432, SSPRKMSVS32, SSPRKMSVS43, SSPRK932, SSPRK54, SSPRK104 diff --git a/src/alg_utils.jl b/src/alg_utils.jl index 52da7dfb7e..26a6e0cf06 100644 --- a/src/alg_utils.jl +++ b/src/alg_utils.jl @@ -290,6 +290,7 @@ alg_order(alg::Hairer42) = 4 alg_order(alg::Feagin10) = 10 alg_order(alg::Feagin12) = 12 alg_order(alg::Feagin14) = 14 +alg_order(alg::PFRK87) = 8 alg_order(alg::Rosenbrock23) = 2 alg_order(alg::Rosenbrock32) = 3 diff --git a/src/algorithms.jl b/src/algorithms.jl index 64f958c0d4..c92f07514c 100644 --- a/src/algorithms.jl +++ b/src/algorithms.jl @@ -403,6 +403,11 @@ struct FRK65{T} <: OrdinaryDiffEqAdaptiveAlgorithm FRK65(omega=0.0) = new{typeof(omega)}(omega) end +struct PFRK87{T} <: OrdinaryDiffEqAdaptiveAlgorithm + omega::T + PFRK87(omega=0.0) = new{typeof(omega)}(omega) +end + ################################################################################ diff --git a/src/caches/high_order_rk_caches.jl b/src/caches/high_order_rk_caches.jl index 1cad1ba6bf..0910d5d8c1 100644 --- a/src/caches/high_order_rk_caches.jl +++ b/src/caches/high_order_rk_caches.jl @@ -123,3 +123,50 @@ function alg_cache(alg::TsitPap8,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNo end alg_cache(alg::TsitPap8,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Val{false}) = TsitPap8ConstantCache(constvalue(uBottomEltypeNoUnits),constvalue(tTypeNoUnits)) + +@cache struct PFRK87Cache{uType,rateType,uNoUnitsType,TabType} <: OrdinaryDiffEqMutableCache + u::uType + uprev::uType + fsalfirst::rateType + k2::rateType + k3::rateType + k4::rateType + k5::rateType + k6::rateType + k7::rateType + k8::rateType + k9::rateType + k10::rateType + k11::rateType + k12::rateType + k13::rateType + utilde::uType + tmp::uType + atmp::uNoUnitsType + k::rateType + tab::TabType +end + +function alg_cache(alg::PFRK87,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Val{true}) + tab = PFRK87ConstantCache(constvalue(uBottomEltypeNoUnits),constvalue(tTypeNoUnits)) + k1 = zero(rate_prototype) + k2 = zero(rate_prototype) + k3 = zero(rate_prototype) + k4 = zero(rate_prototype) + k5 = zero(rate_prototype) + k6 = zero(rate_prototype) + k7 = zero(rate_prototype) + k8 = zero(rate_prototype) + k9 = zero(rate_prototype) + k10 = zero(rate_prototype) + k11 = zero(rate_prototype) + k12 = zero(rate_prototype) + k13 = zero(rate_prototype) + utilde = similar(u) + k = zero(rate_prototype) + tmp = similar(u) + atmp = similar(u,uEltypeNoUnits) + PFRK87Cache(u,uprev,k1,k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,utilde,tmp,atmp,k,tab) +end + +alg_cache(alg::PFRK87,u,rate_prototype,uEltypeNoUnits,uBottomEltypeNoUnits,tTypeNoUnits,uprev,uprev2,f,t,dt,reltol,p,calck,::Val{false}) = PFRK87ConstantCache(constvalue(uBottomEltypeNoUnits),constvalue(tTypeNoUnits)) diff --git a/src/perform_step/high_order_rk_perform_step.jl b/src/perform_step/high_order_rk_perform_step.jl index 094707b9d5..1d0f7bf9a2 100644 --- a/src/perform_step/high_order_rk_perform_step.jl +++ b/src/perform_step/high_order_rk_perform_step.jl @@ -601,3 +601,176 @@ end integrator.destats.nf += 1 end =# + +function initialize!(integrator, cache::PFRK87ConstantCache) + integrator.fsalfirst = integrator.f(integrator.uprev, integrator.p, integrator.t) # Pre-start fsal + integrator.destats.nf += 1 + integrator.kshortsize = 2 + integrator.k = typeof(integrator.k)(undef, integrator.kshortsize) + + # Avoid undefined entries if k is an array of arrays + integrator.fsallast = zero(integrator.fsalfirst) + integrator.k[1] = integrator.fsalfirst + integrator.k[2] = integrator.fsallast +end + +@muladd function perform_step!(integrator, cache::PFRK87ConstantCache, repeat_step=false) + @unpack t,dt,uprev,u,f,p = integrator + @unpack α0201, α0301, α0401, α0501, α0601, α0701, α0302, α0403, α0503, α0504, α0604, α0704, α0605, α0705, α0706, α0908, α1008, α1108, α1208, α1308, α1009, α1109, α1209, α1309, α1110, α1210, α1310, α1211, α1311, β1, β6, β7, β8, β9, β10, β11, β12, β13, β1tilde, β6tilde, β7tilde, β8tilde, β9tilde, β10tilde, β11tilde, β12tilde, c2, c3, c4, c5, c6, c7, c8, c9, c10, c11, c12, c13 = cache + alg = unwrap_alg(integrator, true) + ν = alg.omega*dt + νsq = ν^2 + + c_ν = -0.19781108078634084 - (0.164050909125528499*νsq) + (0.042578310088756321*(νsq^2)) - (0.002300513610963998*(νsq^3)) + (0.000033467244551879287*(νsq^4)) - (7.8661142036921924*(1/100000000)*(νsq^5)) + d_ν = 1-(0.296457092123567400*νsq) + (0.0015793885907465726*(νsq^2)) - (0.00018913011771688527*(νsq^3)) + (0.000017089234650765179*(νsq^4)) - (1.2705211682518626*(1/10000000)*(νsq^5)) + + α0807 = c_ν/d_ν + α0801 = 0.026876256 + 0.0576576*α0807 + α0804 = 0.22464336 - 0.944944*α0807 + α0805 = 0.000369024 - 0.2061696*α0807 + α0806 = 0.21311136 + 0.093456*α0807 + α0901 = 0.07239997637512857 + 0.01913119863380767*α0807 + α0904 = -0.688400520601143 - 0.3135390887207368*α0807 + α0905 = -0.17301267570583073 - 0.06840852844816077*α0807 + α0906 = 0.1440060555560846 + 0.031009360422930017*α0807 + α0907 = 0.9982362892760762 + 0.33180705811215994*α0807 + α1001 = 0.16261514523236525 - 0.12125171966747463*α0807 + α1004 = -2.1255544052061124 + 1.9871809612169453*α0807 + α1005 = -0.216403903283323 + 0.43356675517460624*α0807 + α1006 = -0.060417230254934076 - 0.1965343807796979*α0807 + α1007 = 2.4846281621788395 - 2.102961615944379*α0807 + α1101 = -1.0320124180911034 + 1.061943768952537*α0807 + α1104 = 13.666683232895137 - 17.40407843561103*α0807 + α1105 = 0.25990355211486116 - 3.797253476860588*α0807 + α1106 = -5.759316475814002 + 1.7212824826428488*α0807 + α1107 = -12.822511612651839 + 18.41810566087623*α0807 + α1201 = 0.2478349764611783 - 0.06383934946543009*α0807 + α1204 = -4.593782880309185 + 1.046256005127882*α0807 + α1205 = -0.39566692537411896 + 0.22827403748244698*α0807 + α1206 = -3.0673550479691665 - 0.10347586863902129*α0807 + α1207 = 5.386688702227177 - 1.1072148245058775*α0807 + α1301 = 0.7332242174431163 - 0.5164807626867616*α0807 + α1304 = -10.196728938160977 + 8.464545832921925*α0807 + α1305 = -0.43865244706547707 + 1.846809999910238*α0807 + α1306 = 0.5693856884667226 - 0.8371528845746959*α0807 + α1307 = 10.52865228002416 - 8.957722185570706*α0807 + + k1 = integrator.fsalfirst + k2 = f(uprev+dt*α0201*k1, p, t + c2*dt) + k3 = f(uprev+dt*(α0301*k1+α0302*k2), p, t + c3*dt) + k4 = f(uprev+dt*(α0401*k1 +α0403*k3), p, t + c4*dt) + k5 = f(uprev+dt*(α0501*k1 +α0503*k3+α0504*k4), p, t + c5*dt) + k6 = f(uprev+dt*(α0601*k1 +α0604*k4+α0605*k5), p, t + c6*dt) + k7 = f(uprev+dt*(α0701*k1 +α0704*k4+α0705*k5+α0706*k6), p, t + c7*dt) + k8 = f(uprev+dt*(α0801*k1 +α0804*k4+α0805*k5+α0806*k6+α0807*k7), p, t + c8*dt) + k9 = f(uprev+dt*(α0901*k1 +α0904*k4+α0905*k5+α0906*k6+α0907*k7+α0908*k8), p, t + c9*dt) + k10 =f(uprev+dt*(α1001*k1 +α1004*k4+α1005*k5+α1006*k6+α1007*k7+α1008*k8+α1009*k9), p, t + c10*dt) + k11= f(uprev+dt*(α1101*k1 +α1104*k4+α1105*k5+α1106*k6+α1107*k7+α1108*k8+α1109*k9+α1110*k10), p, t + c11*dt) + k12= f(uprev+dt*(α1201*k1 +α1204*k4+α1205*k5+α1206*k6+α1207*k7+α1208*k8+α1209*k9+α1210*k10+α1211*k11), p, t+c12*dt) + k13= f(uprev+dt*(α1301*k1 +α1304*k4+α1305*k5+α1306*k6+α1307*k7+α1308*k8+α1309*k9+α1310*k10+α1311*k11), p, t+c13*dt) + integrator.destats.nf += 12 + u = uprev + dt*(β1*k1+β6*k6+β7*k7+β8*k8+β9*k9+β10*k10+β11*k11+β12*k12+β13*k13) + if integrator.opts.adaptive + utilde = dt*(β1tilde*k1 + β6tilde*k6 + β7tilde*k7 + β8tilde*k8 + β9tilde*k9 + β10tilde*k10 + β11tilde*k11 + β12tilde*k12) + atmp = calculate_residuals(utilde, uprev, u, integrator.opts.abstol, integrator.opts.reltol,integrator.opts.internalnorm,t) + integrator.EEst = integrator.opts.internalnorm(atmp,t) + end + integrator.fsallast = f(u, p, t+dt) + integrator.destats.nf += 1 + integrator.k[1] = integrator.fsalfirst + integrator.k[2] = integrator.fsallast + integrator.u = u +end + + +function initialize!(integrator, cache::PFRK87Cache) + integrator.fsalfirst = cache.fsalfirst + integrator.fsallast = cache.k + integrator.kshortsize = 2 + resize!(integrator.k, integrator.kshortsize) + integrator.k[1] = integrator.fsalfirst + integrator.k[2] = integrator.fsallast + integrator.f(integrator.fsalfirst, integrator.uprev, integrator.p, integrator.t) # Pre-start fsal + integrator.destats.nf += 1 +end + +@muladd function perform_step!(integrator, cache::PFRK87Cache, repeat_step=false) + @unpack t,dt,uprev,u,f,p = integrator + @unpack α0201, α0301, α0401, α0501, α0601, α0701, α0302, α0403, α0503, α0504, α0604, α0704, α0605, α0705, α0706, α0908, α1008, α1108, α1208, α1308, α1009, α1109, α1209, α1309, α1110, α1210, α1310, α1211, α1311, β1, β6, β7, β8, β9, β10, β11, β12, β13, β1tilde, β6tilde, β7tilde, β8tilde, β9tilde, β10tilde, β11tilde, β12tilde, c2, c3, c4, c5, c6, c7, c8, c9, c10, c11, c12, c13 = cache.tab + @unpack k2,k3,k4,k5,k6,k7,k8,k9,k10,k11,k12,k13,utilde,tmp,atmp,k = cache + + alg = unwrap_alg(integrator, true) + ν = alg.omega*dt + νsq = ν^2 + + c_ν = -0.19781108078634084 - (0.164050909125528499*νsq) + (0.042578310088756321*(νsq^2)) - (0.002300513610963998*(νsq^3)) + (0.000033467244551879287*(νsq^4)) - (7.8661142036921924*(10^(-8))*(νsq^5)) + d_ν = 1-(0.296457092123567400*νsq) + (0.0015793885907465726*(νsq^2)) - (0.00018913011771688527*(νsq^3)) + (0.000017089234650765179*(νsq^4)) - (1.2705211682518626*(10^(-7))*(νsq^5)) + + α0807 = c_ν/d_ν + α0801 = 0.026876256 + 0.0576576*α0807 + α0804 = 0.22464336 - 0.944944*α0807 + α0805 = 0.000369024 - 0.2061696*α0807 + α0806 = 0.21311136 + 0.093456*α0807 + α0901 = 0.07239997637512857 + 0.01913119863380767*α0807 + α0904 = -0.688400520601143 - 0.3135390887207368*α0807 + α0905 = -0.17301267570583073 - 0.06840852844816077*α0807 + α0906 = 0.1440060555560846 + 0.031009360422930017*α0807 + α0907 = 0.9982362892760762 + 0.33180705811215994*α0807 + α1001 = 0.16261514523236525 - 0.12125171966747463*α0807 + α1004 = -2.1255544052061124 + 1.9871809612169453*α0807 + α1005 = -0.216403903283323 + 0.43356675517460624*α0807 + α1006 = -0.060417230254934076 - 0.1965343807796979*α0807 + α1007 = 2.4846281621788395 - 2.102961615944379*α0807 + α1101 = -1.0320124180911034 + 1.061943768952537*α0807 + α1104 = 13.666683232895137 - 17.40407843561103*α0807 + α1105 = 0.25990355211486116 - 3.797253476860588*α0807 + α1106 = -5.759316475814002 + 1.7212824826428488*α0807 + α1107 = -12.822511612651839 + 18.41810566087623*α0807 + α1201 = 0.2478349764611783 - 0.06383934946543009*α0807 + α1204 = -4.593782880309185 + 1.046256005127882*α0807 + α1205 = -0.39566692537411896 + 0.22827403748244698*α0807 + α1206 = -3.0673550479691665 - 0.10347586863902129*α0807 + α1207 = 5.386688702227177 - 1.1072148245058775*α0807 + α1301 = 0.7332242174431163 - 0.5164807626867616*α0807 + α1304 = -10.196728938160977 + 8.464545832921925*α0807 + α1305 = -0.43865244706547707 + 1.846809999910238*α0807 + α1306 = 0.5693856884667226 - 0.8371528845746959*α0807 + α1307 = 10.52865228002416 - 8.957722185570706*α0807 + + k1 = cache.fsalfirst + f(k1, uprev, p, t) + @.. tmp = uprev+dt*α0201*k1 + f(k2, tmp, p, t + c2*dt) + @.. tmp = uprev+dt*(α0301*k1+α0302*k2) + f(k3, tmp, p, t + c3*dt) + @.. tmp = uprev+dt*(α0401*k1+α0403*k3) + f(k4, tmp, p, t + c4*dt) + @.. tmp = uprev+dt*(α0501*k1+α0503*k3+α0504*k4) + f(k5, tmp, p, t + c5*dt) + @.. tmp = uprev+dt*(α0601*k1+α0604*k4+α0605*k5) + f(k6, tmp, p, t + c6*dt) + @.. tmp = uprev+dt*(α0701*k1+α0704*k4+α0705*k5+α0706*k6) + f(k7, tmp, p, t + c7*dt) + @.. tmp = uprev+dt*(α0801*k1+α0804*k4+α0805*k5+α0806*k6+α0807*k7) + f(k8, tmp, p, t + c8*dt) + @.. tmp = uprev+dt*(α0901*k1+α0904*k4+α0905*k5+α0906*k6+α0907*k7+α0908*k8) + f(k9, tmp, p, t + c9*dt) + @.. tmp = uprev+dt*(α1001*k1+α1004*k4+α1005*k5+α1006*k6+α1007*k7+α1008*k8+α1009*k9) + f(k10, tmp, p, t + c10*dt) + @.. tmp = uprev+dt*(α1101*k1+α1104*k4+α1105*k5+α1106*k6+α1107*k7+α1108*k8+α1109*k9+α1110*k10) + f(k11, tmp, p, t + c11*dt) + @.. tmp = uprev+dt*(α1201*k1+α1204*k4+α1205*k5+α1206*k6+α1207*k7+α1208*k8+α1209*k9+α1210*k10+α1211*k11) + f(k12, tmp, p, t+c12*dt) + @.. tmp = uprev+dt*(α1301*k1+α1304*k4+α1305*k5+α1306*k6+α1307*k7+α1308*k8+α1309*k9+α1310*k10+α1311*k11) + f(k13, tmp, p, t+c13*dt) + @.. u = uprev + dt*(β1*k1+β6*k6+β7*k7+β8*k8+β9*k9+β10*k10+β11*k11+β12*k12+β13*k13) + integrator.destats.nf += 13 + if integrator.opts.adaptive + @.. utilde = dt*(β1tilde*k1 + β6tilde*k6 + β7tilde*k7 + β8tilde*k8 + β9tilde*k9 + β10tilde*k10 + β11tilde*k11 + β12tilde*k12) + calculate_residuals!(atmp, utilde, uprev, u, integrator.opts.abstol, integrator.opts.reltol,integrator.opts.internalnorm,t) + integrator.EEst = integrator.opts.internalnorm(atmp,t) + end + f(k, u, p, t+dt) + integrator.destats.nf += 1 + return nothing +end \ No newline at end of file diff --git a/src/tableaus/high_order_rk_tableaus.jl b/src/tableaus/high_order_rk_tableaus.jl index 07b54bf5a6..8705ba6f94 100644 --- a/src/tableaus/high_order_rk_tableaus.jl +++ b/src/tableaus/high_order_rk_tableaus.jl @@ -1020,3 +1020,155 @@ function DP8Interp_polyweights(T::Type) return d401,d406,d407,d408,d409,d410,d411,d412,d413,d414,d415,d416,d501,d506,d507,d508,d509,d510,d511,d512,d513,d514,d515,d516,d601,d606,d607,d608,d609,d610,d611,d612,d613,d614,d615,d616,d701,d706,d707,d708,d709,d710,d711,d712,d713,d714,d715,d716 end +struct PFRK87ConstantCache{T1,T2} <: OrdinaryDiffEqConstantCache + + α0201::T1 + α0301::T1 + α0401::T1 + α0501::T1 + α0601::T1 + α0701::T1 + + α0302::T1 + + α0403::T1 + α0503::T1 + + α0504::T1 + α0604::T1 + α0704::T1 + + α0605::T1 + α0705::T1 + + α0706::T1 + + α0908::T1 + α1008::T1 + α1108::T1 + α1208::T1 + α1308::T1 + + α1009::T1 + α1109::T1 + α1209::T1 + α1309::T1 + + α1110::T1 + α1210::T1 + α1310::T1 + + α1211::T1 + α1311::T1 + + β1::T1 + β6::T1 + β7::T1 + β8::T1 + β9::T1 + β10::T1 + β11::T1 + β12::T1 + β13::T1 + + β1tilde::T1 + β6tilde::T1 + β7tilde::T1 + β8tilde::T1 + β9tilde::T1 + β10tilde::T1 + β11tilde::T1 + β12tilde::T1 + + c2::T2 + c3::T2 + c4::T2 + c5::T2 + c6::T2 + c7::T2 + c8::T2 + c9::T2 + c10::T2 + c11::T2 + c12::T2 + c13::T2 +end + + +function PFRK87ConstantCache(T1::Type, T2::Type) + + # elements of Butcher Table + α0201 = convert(T1,1//18) + α0301 = convert(T1,1//48) + α0401 = convert(T1, 1//32) + α0501 = convert(T1, 5//16) + α0601 = convert(T1, 3//80) + α0701 = convert(T1, 29443841//614563906) + + + α0302 = convert(T1, 1//16) + + α0403 = convert(T1, 3//32) + α0503 = convert(T1, -75//64) + + α0504 = convert(T1, 75//64) + α0604 = convert(T1, 3//16) + α0704 = convert(T1, 77736538//692538347) + + α0605 = convert(T1, 3//20) + α0705 = convert(T1, -28693883//1125000000) + + α0706 = convert(T1, 23124283//1800000000) + + α0908 = convert(T1, 800635310//3783071287) + α1008 = convert(T1, 393006217//1396673457) + α1108 = convert(T1, 15336726248//1032824649) + α1208 = convert(T1, 5232866602//850066563) + α1308 = convert(T1, -13158990841//6184727034) + + α1009 = convert(T1, 123872331//1001029789) + α1109 = convert(T1, -45442868181//3398467696) + α1209 = convert(T1, -4093664535//808688257) + α1309 = convert(T1, 3936647629//1978049680) + + α1110 = convert(T1, 3065993473//597172653) + α1210 = convert(T1, 3962137247//1805957418) + α1310 = convert(T1, -160528059//685178525) + + α1211 = convert(T1, 65686358//487910083) + α1311 = convert(T1, 248638103//1413531060) + + β1 = convert(T1, 14005451//335480064) + β6 = convert(T1, -59238493//1068277825) + β7 = convert(T1, 181606767//758867731) + β8 = convert(T1, 561292985//797845732) + β9 = convert(T1, -1041891430//1371343529) + β10 = convert(T1, 760417239//1151165299) + β11 = convert(T1, 118820643//751138087) + β12 = convert(T1, -528747749//2220607170) + β13 = convert(T1, 1//4) + + β1tilde = convert(T1, 13451932//455176623) + β6tilde = convert(T1, -808719846//976000145) + β7tilde = convert(T1, 1757004468//5645159321) + β8tilde = convert(T1, 656045339//265891186) + β9tilde = convert(T1, -3867574721//1518517206) + β10tilde = convert(T1, 465885868//322736535) + β11tilde = convert(T1, 53011238//667516719) + β12tilde = convert(T1, 2//45) + + c2 = convert(T2, 1//18) + c3 = convert(T2, 1//12) + c4 = convert(T2, 1//8) + c5 = convert(T2, 5//16) + c6 = convert(T2, 3//8) + c7 = convert(T2, 59//400) + c8 = convert(T2, 93//200) + c9 = convert(T2, 5490023248//9719169821) + c10 = convert(T2, 13//20) + c11 = convert(T2, 1201146811//1299019798) + c12 = convert(T2, 1//1) + c13 = convert(T2, 1//1) + + PFRK87ConstantCache(α0201, α0301, α0401, α0501, α0601, α0701, α0302, α0403, α0503, α0504, α0604, α0704, α0605, α0705, α0706, α0908, α1008, α1108, α1208, α1308, α1009, α1109, α1209, α1309, α1110, α1210, α1310, α1211, α1311, β1, β6, β7, β8, β9, β10, β11, β12, β13, β1tilde, β6tilde, β7tilde, β8tilde, β9tilde, β10tilde, β11tilde, β12tilde, c2, c3, c4, c5, c6, c7, c8, c9, c10, c11, c12, c13) +end \ No newline at end of file diff --git a/test/algconvergence/ode_convergence_tests.jl b/test/algconvergence/ode_convergence_tests.jl index 2a18b50f73..855c44db41 100644 --- a/test/algconvergence/ode_convergence_tests.jl +++ b/test/algconvergence/ode_convergence_tests.jl @@ -11,6 +11,7 @@ dts1 = 1 .//2 .^(9:-1:5) dts2 = 1 .//2 .^(7:-1:3) dts3 = 1 .//2 .^(12:-1:7) dts4 = 1 .//2 .^(5:-1:3) +dts5 = 1 .//2 .^(3:-1:1) testTol = 0.2 @testset "Explicit Solver Convergence Tests ($(["out-of-place", "in-place"][i]))" for i in 1:2 @@ -40,6 +41,9 @@ testTol = 0.2 sim3 = test_convergence(dts4,prob,FRK65()) @test sim3.𝒪est[:l∞] ≈ 6 atol=0.6 + sim3 = test_convergence(dts5,prob,PFRK87()) + @test sim3.𝒪est[:l∞] ≈ 8.4 atol=0.2 + sim4 = test_convergence(dts,prob,BS3()) @test sim4.𝒪est[:l2] ≈ 3 atol=testTol sim5 = test_convergence(dts, prob, AB3())