diff --git a/.github/workflows/CI.yml b/.github/workflows/CI.yml index 62b0fe04..2adb3f9f 100644 --- a/.github/workflows/CI.yml +++ b/.github/workflows/CI.yml @@ -10,8 +10,8 @@ jobs: fail-fast: false matrix: os: [ubuntu-latest] - gcc_v: [10] # Version of GFortran we want to use. - python-version: [3.9] + gcc_v: [14] # Version of GFortran we want to use. + python-version: [3.12] env: FC: gfortran-${{ matrix.gcc_v }} GCC_V: ${{ matrix.gcc_v }} @@ -31,7 +31,7 @@ jobs: uses: ts-graphviz/setup-graphviz@v1 - name: Setup Fortran Package Manager - uses: fortran-lang/setup-fpm@v5 + uses: fortran-lang/setup-fpm@v7 with: github-token: ${{ secrets.GITHUB_TOKEN }} diff --git a/README.md b/README.md index 806e14e7..57eb4dd6 100644 --- a/README.md +++ b/README.md @@ -24,9 +24,9 @@ The focus of this library is single-step, explicit Runge-Kutta solvers for 1st o ### Available Runge-Kutta methods: - * Number of fixed-step methods: 27 + * Number of fixed-step methods: 28 * Number of variable-step methods: 48 - * Total number of methods: 75 + * Total number of methods: 76 ### Fixed-step methods: @@ -56,6 +56,7 @@ Name | Description| Properties | Order | Stages | Registers | CFL | Ref `rk8_10` | 10-stage, 8th order Runge-Kutta Shanks | | 8 | 10 | 10 | | [Shanks (1965)](http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19650022581.pdf) `rkcv8` | 11-stage, 8th order Runge-Kutta Cooper-Verner | | 8 | 11 | 11 | | [Cooper & Verner (1972)](https://epubs.siam.org/doi/abs/10.1137/0709037) `rk8_12` | 12-stage, 8th order Runge-Kutta Shanks | | 8 | 12 | 12 | | [Shanks (1965)](http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19650022581.pdf) +`rks10` | 10th order Runge-Kutta Stepanov | | 10 | 15 | 15 | | [Stepanov (2025)](https://arxiv.org/abs/2504.17329) `rkz10` | 10th order Runge-Kutta Zhang | | 10 | 16 | 16 | | [Zhang (2019)](https://arxiv.org/abs/1911.00318) `rko10` | 10th order Runge-Kutta Ono | | 10 | 17 | 17 | | [Ono (2003)](http://www.peterstone.name/Maplepgs/Maple/nmthds/RKcoeff/Runge_Kutta_schemes/RK10/RKcoeff10f_1.pdf) `rkh10` | 10th order Runge-Kutta Hairer | | 10 | 17 | 17 | | [Hairer (1978)](https://www.researchgate.net/publication/31221486_A_Runge-Kutta_Method_of_Order_10) @@ -182,7 +183,9 @@ xf = -0.1360372426E+01 0.1325538438E+01 ### Example performance comparison -Running the unit tests will generate some performance plots. The following is for the variable-step methods compiled with quadruple precision (e.g, `fpm test rk_test_variable_step --compiler ifort --flag "-DREAL128"`): [rk_test_variable_step_R16.pdf](media/rk_test_variable_step_R16.pdf) +Running the unit tests will generate some performance plots. The following is for the variable-step methods compiled with quadruple precision (e.g, `fpm test rk_test_variable_step --compiler ifort --flag "-DREAL128"`): + +![rk_test_variable_step_R16](media/rk_test_variable_step_R16.pdf) ### Compiling diff --git a/scripts/generate_files.py b/scripts/generate_files.py index db7dae12..084b4fe7 100644 --- a/scripts/generate_files.py +++ b/scripts/generate_files.py @@ -34,6 +34,7 @@ ('rk8_10' , '10-stage, 8th order Runge-Kutta Shanks' , ' ', 8 , 10 , 10 , None , '[Shanks (1965)](http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19650022581.pdf)'), ('rkcv8' , '11-stage, 8th order Runge-Kutta Cooper-Verner' , ' ', 8 , 11 , 11 , None , '[Cooper & Verner (1972)](https://epubs.siam.org/doi/abs/10.1137/0709037)'), ('rk8_12' , '12-stage, 8th order Runge-Kutta Shanks' , ' ', 8 , 12 , 12 , None , '[Shanks (1965)](http://ntrs.nasa.gov/archive/nasa/casi.ntrs.nasa.gov/19650022581.pdf)'), + ('rks10' , '10th order Runge-Kutta Stepanov' , ' ', 10 , 15 , 15 , None , '[Stepanov (2025)](https://arxiv.org/abs/2504.17329)'), ('rkz10' , '10th order Runge-Kutta Zhang' , ' ', 10 , 16 , 16 , None , '[Zhang (2019)](https://arxiv.org/abs/1911.00318)'), ('rko10' , '10th order Runge-Kutta Ono' , ' ', 10 , 17 , 17 , None , '[Ono (2003)](http://www.peterstone.name/Maplepgs/Maple/nmthds/RKcoeff/Runge_Kutta_schemes/RK10/RKcoeff10f_1.pdf)'), ('rkh10' , '10th order Runge-Kutta Hairer' , ' ', 10 , 17 , 17 , None , '[Hairer (1978)](https://www.researchgate.net/publication/31221486_A_Runge-Kutta_Method_of_Order_10)')] @@ -105,8 +106,6 @@ def readme_template(): ### Description -**This is a work in progress!** - The focus of this library is single-step, explicit Runge-Kutta solvers for 1st order differential equations. ### Novel features: @@ -153,7 +152,7 @@ def readme_template(): Running the unit tests will generate some performance plots. The following is for the variable-step methods compiled with quadruple precision (e.g, `fpm test rk_test_variable_step --compiler ifort --flag "-DREAL128"`): -![rk_test_variable_step_R16](media/rk_test_variable_step_R16.png) +![rk_test_variable_step_R16](media/rk_test_variable_step_R16.pdf) ### Compiling diff --git a/src/rklib_fixed_classes.inc b/src/rklib_fixed_classes.inc index 1cba8b3f..dbcc1066 100644 --- a/src/rklib_fixed_classes.inc +++ b/src/rklib_fixed_classes.inc @@ -168,6 +168,13 @@ procedure :: properties => rk8_12_properties end type rk8_12_class + type,extends(rk_fixed_step_class),public :: rks10_class + !! 10th order Runge-Kutta Stepanov + contains + procedure :: step => rks10 + procedure :: properties => rks10_properties + end type rks10_class + type,extends(rk_fixed_step_class),public :: rkz10_class !! 10th order Runge-Kutta Zhang contains diff --git a/src/rklib_fixed_properties.f90 b/src/rklib_fixed_properties.f90 index 3b12913a..ac20d58d 100644 --- a/src/rklib_fixed_properties.f90 +++ b/src/rklib_fixed_properties.f90 @@ -238,6 +238,15 @@ p%number_of_registers = 12 end procedure rk8_12_properties + module procedure rks10_properties + !! Returns the properties of the [[rks10]] method + p%short_name = 'rks10' + p%long_name = '10th order Runge-Kutta Stepanov' + p%order = 10 + p%number_of_stages = 15 + p%number_of_registers = 15 + end procedure rks10_properties + module procedure rkz10_properties !! Returns the properties of the [[rkz10]] method p%short_name = 'rkz10' diff --git a/src/rklib_fixed_property_interfaces.inc b/src/rklib_fixed_property_interfaces.inc index 1699bd4e..7832eaab 100644 --- a/src/rklib_fixed_property_interfaces.inc +++ b/src/rklib_fixed_property_interfaces.inc @@ -142,6 +142,12 @@ pure module function rk8_12_properties(me) result(p) type(rklib_properties) :: p !! properties of the method end function rk8_12_properties +pure module function rks10_properties(me) result(p) + implicit none + class(rks10_class),intent(in) :: me + type(rklib_properties) :: p !! properties of the method +end function rks10_properties + pure module function rkz10_properties(me) result(p) implicit none class(rkz10_class),intent(in) :: me diff --git a/src/rklib_fixed_step_interfaces.inc b/src/rklib_fixed_step_interfaces.inc index fbc955db..7979c1c7 100644 --- a/src/rklib_fixed_step_interfaces.inc +++ b/src/rklib_fixed_step_interfaces.inc @@ -216,6 +216,15 @@ real(wp),dimension(me%n),intent(out) :: xf !! state at time `t+h` end subroutine rk8_12 + module subroutine rks10(me,t,x,h,xf) + implicit none + class(rks10_class),intent(inout) :: me + real(wp),intent(in) :: t !! initial time + real(wp),dimension(me%n),intent(in) :: x !! initial state + real(wp),intent(in) :: h !! time step + real(wp),dimension(me%n),intent(out) :: xf !! state at time `t+h` + end subroutine rks10 + module subroutine rkz10(me,t,x,h,xf) implicit none class(rkz10_class),intent(inout) :: me diff --git a/src/rklib_fixed_steps.f90 b/src/rklib_fixed_steps.f90 index 26fc1f05..6bad11b3 100644 --- a/src/rklib_fixed_steps.f90 +++ b/src/rklib_fixed_steps.f90 @@ -1842,6 +1842,170 @@ end procedure rkh10 !***************************************************************************************** +!***************************************************************************************** +!> +! Stepanov 10th order (15 stage) method. +! +!### References +! * Misha Stepanov, "On Runge-Kutta methods of order 10" +! arXiv:2504.17329, 24 Apr 2025 + + module procedure rks10 + + ! rounded decimal form (90 digits) of the coefficients provided by Stepanov. + + real(wp),parameter :: a2 = +0.133333333333333333333333333333333333333333333333333333333333333333333333333333333333333333_wp + real(wp),parameter :: a3 = +0.266666666666666666666666666666666666666666666666666666666666666666666666666666666666666667_wp + real(wp),parameter :: a4 = +0.400000000000000000000000000000000000000000000000000000000000000000000000000000000000000000_wp + real(wp),parameter :: a5 = +0.571428571428571428571428571428571428571428571428571428571428571428571428571428571428571429_wp + real(wp),parameter :: a6 = +0.778740761536291800442363524550407562438954342313508170658506769545868148577329534757889703_wp + real(wp),parameter :: a7 = +0.642615758240322548157075497020439535959501736363212695909875208263848965457099799090837864_wp + real(wp),parameter :: a8 = +0.882527661964732346425501486979669075182867844268052119663791177918527658519413257061748635_wp + real(wp),parameter :: a9 = +0.117472338035267653574498513020330924817132155731947880336208822081472341480586742938251365_wp + real(wp),parameter :: a10 = +0.117472338035267653574498513020330924817132155731947880336208822081472341480586742938251365_wp + real(wp),parameter :: a11 = +0.357384241759677451842924502979560464040498263636787304090124791736151034542900200909162136_wp + real(wp),parameter :: a12 = +0.357384241759677451842924502979560464040498263636787304090124791736151034542900200909162136_wp + real(wp),parameter :: a13 = +0.642615758240322548157075497020439535959501736363212695909875208263848965457099799090837864_wp + real(wp),parameter :: a14 = +0.882527661964732346425501486979669075182867844268052119663791177918527658519413257061748635_wp + + real(wp),parameter :: b21 = +0.133333333333333333333333333333333333333333333333333333333333333333333333333333333333333333_wp + real(wp),parameter :: b32 = +0.266666666666666666666666666666666666666666666666666666666666666666666666666666666666666667_wp + real(wp),parameter :: b41 = +0.100000000000000000000000000000000000000000000000000000000000000000000000000000000000000000_wp + real(wp),parameter :: b43 = +0.300000000000000000000000000000000000000000000000000000000000000000000000000000000000000000_wp + real(wp),parameter :: b51 = +0.134110787172011661807580174927113702623906705539358600583090379008746355685131195335276968_wp + real(wp),parameter :: b53 = +0.087463556851311953352769679300291545189504373177842565597667638483965014577259475218658892_wp + real(wp),parameter :: b54 = +0.349854227405247813411078717201166180758017492711370262390670553935860058309037900874635569_wp + real(wp),parameter :: b61 = +0.328658938997791240721282668696245180100949935428058492348520659513222104304668996831025694_wp + real(wp),parameter :: b63 = -0.843123044711636120617820434267319755889422538822419939805710279613810068291349217990455336_wp + real(wp),parameter :: b64 = +1.230383064721160870667126440845906263177038685065274297872024392536130940503747853752249173_wp + real(wp),parameter :: b65 = +0.062821802528975809671774849275575875050388260642595320243671997110325172060261902165070172_wp + real(wp),parameter :: b71 = +0.130186990844533481397415277676876710451288921648802895948985973227427327921848463186124942_wp + real(wp),parameter :: b74 = +0.580738972917374107341609289970199905398640941037865280565410031214758908228344711848158528_wp + real(wp),parameter :: b75 = -0.132061017042928219570339015811796150446818327131853650470115559422384511796983327996330108_wp + real(wp),parameter :: b76 = +0.063750811521343178988389945185159070556390200808398169865594763244047241103889952052884502_wp + real(wp),parameter :: b81 = +0.116969144664032355608536137983019767284447494556651903864723553337776482921362147717802148_wp + real(wp),parameter :: b84 = +0.804791321585223079797049560718389631085997501079130322088670870267840258541270110417387198_wp + real(wp),parameter :: b85 = -1.032442407303393452179074692383412947850773751656102435103865475343877801453723942152851429_wp + real(wp),parameter :: b86 = +0.141237507577425159424271186625123531674337461239364982255644624276392633209837338226267645_wp + real(wp),parameter :: b87 = +0.851972095441445203774719294036549092988859139049007346558617605380396085300667602853143073_wp + real(wp),parameter :: b91 = +0.072301376815318053808169355285248577196274719038173036516820302048727944927561329493507003_wp + real(wp),parameter :: b94 = +0.278449353793857950736187552433349927872809528174204139729412569505988600192254101364671140_wp + real(wp),parameter :: b95 = -0.720355742529399485982635010432963793377590743017129914556141009831958711185177678815313348_wp + real(wp),parameter :: b96 = +0.071478309723007569077180467110101963232268708471224068127334024952865428728115593534703972_wp + real(wp),parameter :: b97 = +0.480546127100390254710165280839106235150089447220646221776588301884109281650681775624438024_wp + real(wp),parameter :: b98 = -0.064947086867906688774569132214511985256719504155169671257805366478260202832848378263755427_wp + real(wp),parameter :: b101 = -0.033604580075425596431640044865448312697754125457571033966446479388085421715398642048231498_wp + real(wp),parameter :: b104 = -0.857669388668674924032812732229566068998144989807984562742043274869246367307599404275438590_wp + real(wp),parameter :: b105 = +2.054105395119576054328336908587945831582620175191578897288368999063154807831596816873316215_wp + real(wp),parameter :: b106 = -0.209909042661849185066865444237374444876978972036734389443755323786107747791309251353961678_wp + real(wp),parameter :: b107 = -1.333972235951611409051125524672704502068159831145103714187256702629992309144152509633176930_wp + real(wp),parameter :: b108 = +0.180806348384581339693160596923480075192005804935788561945038863986674416934017278251008686_wp + real(wp),parameter :: b109 = +0.317715841888671374135444753513998346683544094051974121442302739705074962673432455124735159_wp + real(wp),parameter :: b111 = +0.037810053861933727583877526864046322952305373764286754688878816133302931065788321936975530_wp + real(wp),parameter :: b114 = +0.045591558641579644208161372472470642588279099631100286855624026273895766568226149784507397_wp + real(wp),parameter :: b115 = +0.374858247294747357801811789348440678077393709386181864992340497763246020363342448831864027_wp + real(wp),parameter :: b116 = -0.018981904165238088665273073099056486944960468054380557714800392623159772511027805727133276_wp + real(wp),parameter :: b117 = -0.336038862637559379967733799098120208915305739878288300584069940705912757234040695695257751_wp + real(wp),parameter :: b118 = +0.042207716375559153022580951807545725160812786625465879746013555646386854001549779891730475_wp + real(wp),parameter :: b119 = +0.051900360983690832037790128319440699859433125187462257771821672390157423238366045436562647_wp + real(wp),parameter :: b1110 = +0.160037071404964205821709606364793091262540376974959118334316556858234569050695956449913087_wp + real(wp),parameter :: b121 = +0.025378474265434522073341785109062596527860323276820540705890698859323818498838966304096456_wp + real(wp),parameter :: b124 = +0.362306327543474080870356990535692676714482207942118981416939928949791728367123473328004120_wp + real(wp),parameter :: b125 = -2.130022066259949711216588139786172085385079053283141571779776466074951544604054464037893691_wp + real(wp),parameter :: b126 = +0.167271594859437715452103295041620121357577378494618057193516457140141705837480748503078071_wp + real(wp),parameter :: b127 = +1.582750941864663213591839501273807297171503649177669091900056915258923612981224332245177665_wp + real(wp),parameter :: b128 = -0.204720185139081567587174987631558173178320360472369067149216918644283831827088288528863705_wp + real(wp),parameter :: b129 = +0.319265328144806137168538521838050367024648119039965645533219052127751799436983299054606914_wp + real(wp),parameter :: b1210 = -0.199485256254696927298131093258193563570518501610697726802791703448058055735145490139989292_wp + real(wp),parameter :: b1211 = +0.434639082735589988788638629857251227378344501071803353072286827567511801587537624180945598_wp + real(wp),parameter :: b131 = +0.021029555123492649231004219196592342935742304160151528452676708946550307466281435746860028_wp + real(wp),parameter :: b134 = -0.208531388344693562066394580214099278549298495126493308746134843353719483236004257692460587_wp + real(wp),parameter :: b135 = +0.326863498654767652665730419629359718148991455519984894344217473947412982044929978522731887_wp + real(wp),parameter :: b136 = -0.040291595063466341977981502499810508092167265608412827254081246436545949894568768120086225_wp + real(wp),parameter :: b137 = -0.062168570110085104851127006245165266954000899589268210982639466796612080783267242335332813_wp + real(wp),parameter :: b138 = +0.019401511635634607470015991594882193791752892731789828413457858797773885042894717414314220_wp + real(wp),parameter :: b139 = +0.033618934322149698701521634179531768990520248587976979368806692500394465464703080284812928_wp + real(wp),parameter :: b1310 = +0.184057539066668989188409585677311173226149887887488899665117691435303709700110659619387592_wp + real(wp),parameter :: b1311 = +0.126719042450063065291362679244566355420313628999369516069616211043815410052049067488162884_wp + real(wp),parameter :: b1312 = +0.241917230505790894504534056457271037041497978800625396578838128179475719599971128162447951_wp + real(wp),parameter :: b141 = +0.059015439601144226763876106824299822651056803159364026683828389696455698859560833618210307_wp + real(wp),parameter :: b144 = +0.346753529400439458709721774328857076918185484739163574624551262985723353031416083991491991_wp + real(wp),parameter :: b145 = -0.543520439250941361637466822737163641756559230304227658933756185028034640683436881957640271_wp + real(wp),parameter :: b146 = +0.066998320513439386022349163118524383402374304344205128122905426297240312187142777357799209_wp + real(wp),parameter :: b147 = -0.209105926273890705326433037391968985907198847287185344061438195214198838012365442723920080_wp + real(wp),parameter :: b148 = +0.034196819215377393468216514309803928830340701888167708646927556886609080880789258681682710_wp + real(wp),parameter :: b149 = +0.656808691640117064523673744875271638362693883629496564722278370329924383275858838375005043_wp + real(wp),parameter :: b1410 = -0.530023086491271629685194176217211641242245756479495833589514996293186125976018346545708951_wp + real(wp),parameter :: b1411 = +0.682370713621575520444713157488350028626121649514411954149882119768342984225842629755650051_wp + real(wp),parameter :: b1412 = -0.524316794649154201841271672103049178047360868487598331909107445957588634113426671093395740_wp + real(wp),parameter :: b1413 = +0.843350394637897194983316734483955643345459719551750331207234874447240084844050177602574367_wp + real(wp),parameter :: b151 = -0.040440866328662678612786513864723972854343311900228087268921187067545862661799086441064095_wp + real(wp),parameter :: b154 = -0.936174412160392166960880644803504126455580435000519336040696957832662353800600021951728948_wp + real(wp),parameter :: b155 = +1.467410954959016212947954342116057777884286077126485056662059170859341898803083590224129382_wp + real(wp),parameter :: b156 = -0.180883849778987051699697942242124022922323870654198891654454984812850334688440411791417790_wp + real(wp),parameter :: b157 = +1.167589102274410575304328368479226265476383761610677844233764974596479007878520383510866872_wp + real(wp),parameter :: b158 = -0.280113604810353377498686142206136702751291429624401159362004219851994095769258196924953998_wp + real(wp),parameter :: b159 = -1.871668376348516420302371567652811494220817967866998659532896736116592408041514651376908974_wp + real(wp),parameter :: b1510 = +2.241588116125051681252015179294157675223848655884912766574651463604345999471163122136357746_wp + real(wp),parameter :: b1511 = -1.572572082691418559566908951416310106551328354276059431821310834756397117889598376197133068_wp + real(wp),parameter :: b1512 = +2.151700399601058696148477230250852973515302452974155466886336399845840759500171462284652967_wp + real(wp),parameter :: b1513 = -1.813340450902764402851969633637645578401662879405552525500029323314042291416407578145425752_wp + real(wp),parameter :: b1514 = +0.666905070061557491840526275682961312057527301131726956823502234846076798614679764672625658_wp + + real(wp),parameter :: c1 = +0.033333333333333333333333333333333333333333333333333333333333333333333333333333333333333333_wp + real(wp),parameter :: c9 = +0.135169627249231064398790288647151661598687390677589878805138875485701561628210909193328656_wp + real(wp),parameter :: c10 = +0.054067850899692425759516115458860664639474956271035951522055550194280624651284363677331462_wp + real(wp),parameter :: c11 = +0.215778257736022470617613537547175598111058915336253983819589520767421262523355528508005093_wp + real(wp),parameter :: c12 = +0.061650930781720705890746725013478742317445404381786852519882720219263217863815865288001455_wp + real(wp),parameter :: c13 = +0.277429188517743176508360262560654340428504319718040836339472240986684480387171393796006548_wp + real(wp),parameter :: c14 = +0.189237478148923490158306404106012326238162346948625830327194425679982186279495272870660119_wp + real(wp),parameter :: c15 = +0.033333333333333333333333333333333333333333333333333333333333333333333333333333333333333333_wp + + associate (f1 => me%funcs(:,1), & + f2 => me%funcs(:,2), & + f3 => me%funcs(:,3), & + f4 => me%funcs(:,4), & + f5 => me%funcs(:,5), & + f6 => me%funcs(:,6), & + f7 => me%funcs(:,7), & + f8 => me%funcs(:,8), & + f9 => me%funcs(:,9), & + f10 => me%funcs(:,10), & + f11 => me%funcs(:,11), & + f12 => me%funcs(:,12), & + f13 => me%funcs(:,13), & + f14 => me%funcs(:,14), & + f15 => me%funcs(:,15)) + + call me%f(t,x,f1) + call me%f(t+a2*h,x+h*(b21*f1),f2) + call me%f(t+a3*h,x+h*(b32*f2),f3) + call me%f(t+a4*h,x+h*(b41*f1+b43*f3),f4) + call me%f(t+a5*h,x+h*(b51*f1+b53*f3+b54*f4),f5) + call me%f(t+a6*h,x+h*(b61*f1+b63*f3+b64*f4+b65*f5),f6) + call me%f(t+a7*h,x+h*(b71*f1+b74*f4+b75*f5+b76*f6),f7) + call me%f(t+a8*h,x+h*(b81*f1+b84*f4+b85*f5+b86*f6+b87*f7),f8) + call me%f(t+a9*h,x+h*(b91*f1+b94*f4+b95*f5+b96*f6+b97*f7+b98*f8),f9) + call me%f(t+a10*h,x+h*(b101*f1+b104*f4+b105*f5+b106*f6+b107*f7+& + b108*f8+b109*f9),f10) + call me%f(t+a11*h,x+h*(b111*f1+b114*f4+b115*f5+b116*f6+b117*f7+& + b118*f8+b119*f9+b1110*f10),f11) + call me%f(t+a12*h,x+h*(b121*f1+b124*f4+b125*f5+b126*f6+b127*f7+& + b128*f8+b129*f9+b1210*f10+b1211*f11),f12) + call me%f(t+a13*h,x+h*(b131*f1+b134*f4+b135*f5+b136*f6+b137*f7+& + b138*f8+b139*f9+b1310*f10+b1311*f11+b1312*f12),f13) + call me%f(t+a14*h,x+h*(b141*f1+b144*f4+b145*f5+b146*f6+b147*f7+& + b148*f8+b149*f9+b1410*f10+b1411*f11+b1412*f12+b1413*f13),f14) + call me%f(t,x+h*(b151*f1+b154*f4+b155*f5+b156*f6+b157*f7+b158*f8+& + b159*f9+b1510*f10+b1511*f11+b1512*f12+b1513*f13+b1514*f14),f15) + + xf = x+h*(c1*f1+c9*f9+c10*f10+c11*f11+c12*f12+c13*f13+c14*f14+c15*f15) + + end associate + + end procedure rks10 +!***************************************************************************************** + !***************************************************************************************** end submodule rklib_fixed_steps !***************************************************************************************** \ No newline at end of file diff --git a/test/rk_test.f90 b/test/rk_test.f90 index 5d5607ab..6dc25720 100644 --- a/test/rk_test.f90 +++ b/test/rk_test.f90 @@ -60,9 +60,11 @@ program rk_test allocate(rk8_10_class :: s); call run_all_tests([122, 52, 235]) allocate(rk8_12_class :: s); call run_all_tests([229, 52, 235]) allocate(rkcv8_class :: s); call run_all_tests([217, 163, 163]) + allocate(rkz10_class :: s); call run_all_tests([222, 115, 73]) allocate(rko10_class :: s); call run_all_tests([200, 200, 200]) allocate(rkh10_class :: s); call run_all_tests([180, 180, 180]) + allocate(rks10_class :: s); call run_all_tests([50, 50, 50]) ! save plot: write(rstr,'(I3)') wp diff --git a/test/rk_test_variable_step.f90 b/test/rk_test_variable_step.f90 index 20bf06a3..4c1c0c40 100644 --- a/test/rk_test_variable_step.f90 +++ b/test/rk_test_variable_step.f90 @@ -117,6 +117,9 @@ program rk_test_variable_step call plt%savefig(figfile='rk_test_variable_step_R'//trim(adjustl(rstr))//'.pdf',istat=istat) call plt2%savefig(figfile='rk_test_variable_step_FIXED_MODE_R'//trim(adjustl(rstr))//'.pdf',istat=istat) + call plt%savefig(figfile='rk_test_variable_step_R'//trim(adjustl(rstr))//'.png',istat=istat) + call plt2%savefig(figfile='rk_test_variable_step_FIXED_MODE_R'//trim(adjustl(rstr))//'.png',istat=istat) + contains !***************************************************************************************** diff --git a/test/rklib_allocate_and_test.inc b/test/rklib_allocate_and_test.inc index 1dc2d609..f8f68265 100644 --- a/test/rklib_allocate_and_test.inc +++ b/test/rklib_allocate_and_test.inc @@ -22,6 +22,7 @@ allocate(rk8_10_class :: s); call run_test() allocate(rkcv8_class :: s); call run_test() allocate(rk8_12_class :: s); call run_test() + allocate(rks10_class :: s); call run_test() allocate(rkz10_class :: s); call run_test() allocate(rko10_class :: s); call run_test() allocate(rkh10_class :: s); call run_test()