diff --git a/CHANGELOG.md b/CHANGELOG.md index 45619e61e3..c27937e5a4 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -14,6 +14,24 @@ Added "Resize" capability, as well as missing `SetRootDirection` and Deprecated `ARKStepSetOptimalParams` function; added instructions to user guide for users who wish to retain the current functionality. +Added the following Runge-Kutta Butcher tables +* `ARKODE_FORWARD_EULER_1_1` +* `ARKODE_RALSTON_EULER_2_1_2` +* `ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2` +* `ARKODE_BACKWARD_EULER_1_1` +* `ARKODE_IMPLICIT_MIDPOINT_1_2` +* `ARKODE_IMPLICIT_TRAPEZOIDAL_2_2` + +Added the following MRI coupling tables +* `ARKODE_MRI_GARK_FORWARD_EULER` +* `ARKODE_MRI_GARK_RALSTON2` +* `ARKODE_MRI_GARK_RALSTON3` +* `ARKODE_MRI_GARK_BACKWARD_EULER` +* `ARKODE_MRI_GARK_IMPLICIT_MIDPOINT` +* `ARKODE_IMEX_MRI_GARK_EULER` +* `ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL` +* `ARKODE_IMEX_MRI_GARK_MIDPOINT` + Updated the CMake variable `HIP_PLATFORM` default to `amd` as the previous default, `hcc`, is no longer recognized in ROCm 5.7.0 or newer. The new default is also valid in older version of ROCm (at least back to version 4.3.1). diff --git a/doc/arkode/guide/source/Butcher.rst b/doc/arkode/guide/source/Butcher.rst index 015b500a04..3d44dbecdc 100644 --- a/doc/arkode/guide/source/Butcher.rst +++ b/doc/arkode/guide/source/Butcher.rst @@ -150,6 +150,37 @@ specified via a unique ID and name: with values specified for each method below (e.g., ``ARKODE_HEUN_EULER_2_1_2``). +.. _Butcher.Forward_Euler: + +Forward-Euler-1-1 +^^^^^^^^^^^^^^^^^ + +.. index:: Forward-Euler-1-1 ERK method + +Accessible via the constant ``ARKODE_FORWARD_EULER_1_1`` to +:c:func:`ARKStepSetTableNum`, :c:func:`ERKStepSetTableNum` or +:c:func:`ARKodeButcherTable_LoadERK`. +Accessible via the string ``"ARKODE_FORWARD_EULER_1_1"`` to +:c:func:`ARKStepSetTableName`, :c:func:`ERKStepSetTableName` or +:c:func:`ARKodeButcherTable_LoadERKByName`. +This is the default 1st order explicit method (from :cite:p:`Euler:68`). + +.. math:: + + \renewcommand{\arraystretch}{1.5} + \begin{array}{r|c} + 0 & 0 \\ + \hline + 1 & 1 + \end{array} + +.. figure:: /figs/arkode/forward_euler_erk_stab_region.png + :scale: 50 % + :align: center + + Linear stability region for the forward Euler method. + + .. _Butcher.Heun_Euler: Heun-Euler-2-1-2 @@ -184,6 +215,74 @@ This is the default 2nd order explicit method. region is outlined in blue; the embedding's region is in red. +.. _Butcher.Ralston_Euler: + +Ralston-Euler-2-1-2 +^^^^^^^^^^^^^^^^^^^^ + +.. index:: Ralston-Euler-2-1-2 ERK method + +Accessible via the constant ``ARKODE_RALSTON_EULER_2_1_2`` to +:c:func:`ARKStepSetTableNum`, :c:func:`ERKStepSetTableNum` or +:c:func:`ARKodeButcherTable_LoadERK`. +Accessible via the string ``"ARKODE_RALSTON_EULER_2_1_2"`` to +:c:func:`ARKStepSetTableName`, :c:func:`ERKStepSetTableName` or +:c:func:`ARKodeButcherTable_LoadERKByName` +(primary method from :cite:p:`Ralston:62`). + +.. math:: + + \renewcommand{\arraystretch}{1.5} + \begin{array}{r|cc} + 0 & 0 & 0 \\ + \frac{2}{3} & \frac{2}{3} & 0 \\ + \hline + 2 & \frac{1}{4} & \frac{3}{4} \\ + 1 & 1 & 0 + \end{array} + +.. figure:: /figs/arkode/ralston_euler_erk_stab_region.png + :scale: 50 % + :align: center + + Linear stability region for the Ralston-Euler method. The method's + region is outlined in blue; the embedding's region is in red. + + +.. _Butcher.Explicit_Midpoint_Euler: + +Explicit-Midpoint-Euler-2-1-2 +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ + +.. index:: Explicit-Midpoint-Euler-2-1-2 ERK method + +Accessible via the constant ``ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2`` to +:c:func:`ARKStepSetTableNum`, :c:func:`ERKStepSetTableNum` or +:c:func:`ARKodeButcherTable_LoadERK`. +Accessible via the string ``"ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2"`` to +:c:func:`ARKStepSetTableName`, :c:func:`ERKStepSetTableName` or +:c:func:`ARKodeButcherTable_LoadERKByName` +(primary method from :cite:p:`Runge:95`). + +.. math:: + + \renewcommand{\arraystretch}{1.5} + \begin{array}{r|cc} + 0 & 0 & 0 \\ + \frac{1}{2} & \frac{1}{2} & 0 \\ + \hline + 2 & 0 & 1 \\ + 1 & 1 & 0 + \end{array} + +.. figure:: /figs/arkode/explicit_midpoint_euler_erk_stab_region.png + :scale: 50 % + :align: center + + Linear stability region for the Explicit-Midpoint-Euler method. The method's + region is outlined in blue; the embedding's region is in red. + + .. _Butcher.ARK2_ERK: ARK2-ERK-3-1-2 @@ -1123,6 +1222,37 @@ specified via a unique ID and name: with values specified for each method below (e.g., ``ARKODE_SDIRK_2_1_2``). +.. _Butcher.Backward-Euler: + +Backward-Euler-1-1 +^^^^^^^^^^^^^^^^^^ + +.. index:: Backward-Euler-1-1 method + +Accessible via the constant ``ARKODE_BACKWARD_EULER_1_1`` to +:c:func:`ARKStepSetTableNum` or +:c:func:`ARKodeButcherTable_LoadDIRK`. +Accessible via the string ``"ARKODE_BACKWARD_EULER_1_1"`` to +:c:func:`ARKStepSetTableName` or +:c:func:`ARKodeButcherTable_LoadDIRKByName`. +This is the default 1st order implicit method. The method is A-, L-, and B-stable. + +.. math:: + + \renewcommand{\arraystretch}{1.5} + \begin{array}{r|c} + 1 & 1 \\ + \hline + 1 & 1 + \end{array} + +.. figure:: /figs/arkode/backward_euler_dirk_stab_region.png + :scale: 50 % + :align: center + + Linear stability region for the backward Euler method. + + .. _Butcher.SDIRK-2-1: SDIRK-2-1-2 @@ -1194,6 +1324,69 @@ implicit portion of the ARK2 method from :cite:p:`giraldo2013implicit`). region is outlined in blue; the embedding's region is in red. +.. _Butcher.Implicit_Midpoint: + +Implicit-Midpoint-1-2 +^^^^^^^^^^^^^^^^^^^^^ + +.. index:: Implicit-Midpoint-1-2 method + +Accessible via the constant ``ARKODE_IMPLICIT_MIDPOINT_1_2`` to +:c:func:`ARKStepSetTableNum` or +:c:func:`ARKodeButcherTable_LoadDIRK`. +Accessible via the string ``"ARKODE_IMPLICIT_MIDPOINT_1_2"`` to +:c:func:`ARKStepSetTableName` or +:c:func:`ARKodeButcherTable_LoadDIRKByName`. +The method is A- and B-stable. + +.. math:: + + \renewcommand{\arraystretch}{1.5} + \begin{array}{r|c} + \frac{1}{2} & \frac{1}{2} \\ + \hline + 2 & 1 + \end{array} + +.. figure:: /figs/arkode/implicit_midpoint_dirk_stab_region.png + :scale: 50 % + :align: center + + Linear stability region for the implicit midpoint method. + + +.. _Butcher.Implicit_Trapezoidal: + +Implicit-Trapezoidal-2-2 +^^^^^^^^^^^^^^^^^^^^^^^^ + +.. index:: Implicit-Trapezoidal-2-2 method + +Accessible via the constant ``ARKODE_IMPLICIT_TRAPEZOIDAL_2_2`` to +:c:func:`ARKStepSetTableNum` or +:c:func:`ARKodeButcherTable_LoadDIRK`. +Accessible via the string ``"ARKODE_IMPLICIT_TRAPEZOIDAL_2_2"`` to +:c:func:`ARKStepSetTableName` or +:c:func:`ARKodeButcherTable_LoadDIRKByName`. +The method is A-stable. + +.. math:: + + \renewcommand{\arraystretch}{1.5} + \begin{array}{r|cc} + 0 & 0 & 0 \\ + 1 & \frac{1}{2} & \frac{1}{2} \\ + \hline + 2 & \frac{1}{2} & \frac{1}{2} + \end{array} + +.. figure:: /figs/arkode/implicit_trapezoidal_dirk_stab_region.png + :scale: 50 % + :align: center + + Linear stability region for the implicit trapezoidal method. + + .. _Butcher.Billington: Billington-3-3-2 diff --git a/doc/arkode/guide/source/Constants.rst b/doc/arkode/guide/source/Constants.rst index 69affcfaa9..088dbbdda8 100644 --- a/doc/arkode/guide/source/Constants.rst +++ b/doc/arkode/guide/source/Constants.rst @@ -74,8 +74,16 @@ contains the ARKODE output constants. +-----------------------------------------------+------------------------------------------------------------+ | **Explicit Butcher table specification** | | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_FORWARD_EULER_1_1` | Use the Forward-Euler-1-1 ERK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_HEUN_EULER_2_1_2` | Use the Heun-Euler-2-1-2 ERK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_RALSTON_EULER_2_1_2` | Use the Ralston-Euler-2-1-2 ERK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2` | Use the Explicit-Midpoint-Euler-2-1-2 ERK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_ARK2_ERK_3_1_2` | Use the ARK2-ERK-3-1-2 ERK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_BOGACKI_SHAMPINE_4_2_3` | Use the Bogacki-Shampine-4-2-3 ERK method. | +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_ARK324L2SA_ERK_4_2_3` | Use the ARK-4-2-3 ERK method. | @@ -116,6 +124,9 @@ contains the ARKODE output constants. +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_VERNER_16_8_9` | Use the Verner-16-8-9 ERK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKSTEP_DEFAULT_ERK_1` | Use ARKStep's default first-order ERK method | + | | (ARKODE_FORWARD_EULER_1_1). | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKSTEP_DEFAULT_ERK_2` | Use ARKStep's default second-order ERK method | | | (ARKODE_HEUN_EULER_2_1_2). | +-----------------------------------------------+------------------------------------------------------------+ @@ -140,6 +151,9 @@ contains the ARKODE output constants. | :index:`ARKSTEP_DEFAULT_ERK_9` | Use ARKStep's default ninth-order ERK method | | | (ARKODE_VERNER_16_8_9). | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ERKSTEP_DEFAULT_1` | Use ERKStep's default first-order ERK method | + | | (ARKODE_FORWARD_EULER_1_1). | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ERKSTEP_DEFAULT_2` | Use ERKStep's default second-order ERK method | | | (ARKODE_HEUN_EULER_2_1_2). | +-----------------------------------------------+------------------------------------------------------------+ @@ -168,8 +182,16 @@ contains the ARKODE output constants. +-----------------------------------------------+------------------------------------------------------------+ | **Implicit Butcher table specification** | | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_BACKWARD_EULER_1_1` | Use the Backward-Euler-1-1 SDIRK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_SDIRK_2_1_2` | Use the SDIRK-2-1-2 SDIRK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_ARK2_DIRK_3_1_2` | Use the ARK2-DIRK-3-1-2 SDIRK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_IMPLICIT_MIDPOINT_1_2` | Use the Implicit-Midpoint-1-2 SDIRK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_IMPLICIT_TRAPEZOIDAL_2_2` | Use the Implicit-Trapezoidal-2-2 ESDIRK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_BILLINGTON_3_3_2` | Use the Billington-3-3-2 SDIRK method. | +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_ESDIRK324L2SA_4_2_3` | Use the ESDIRK324L2SA-4-2-3 ESDIRK method. | @@ -214,6 +236,9 @@ contains the ARKODE output constants. +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_ESDIRK547L2SA2_7_4_5` | Use the ESDIRK547L2SA2-7-4-5 ESDIRK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKSTEP_DEFAULT_DIRK_1` | Use ARKStep's default first-order DIRK method | + | | (ARKODE_BACKWARD_EULER_1_1). | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKSTEP_DEFAULT_DIRK_2` | Use ARKStep's default second-order DIRK method | | | (ARKODE_SDIRK_2_1_2). | +-----------------------------------------------+------------------------------------------------------------+ @@ -230,6 +255,9 @@ contains the ARKODE output constants. +-----------------------------------------------+------------------------------------------------------------+ | **ImEx Butcher table specification** | | +-----------------------------------------------+------------------------------------------------------------+ + | ARKODE_ARK2_ERK_3_1_2 & | Use the :index:`ARK-3-1-2 ARK method`. | + | ARKODE_ARK2_DIRK_3_1_2 | | + +-----------------------------------------------+------------------------------------------------------------+ | ARKODE_ARK324L2SA_ERK_4_2_3 & | Use the :index:`ARK-4-2-3 ARK method`. | | ARKODE_ARK324L2SA_DIRK_4_2_3 | | +-----------------------------------------------+------------------------------------------------------------+ @@ -299,30 +327,59 @@ contains the ARKODE output constants. +-----------------------------------------------+------------------------------------------------------------+ | **MRI coupling table specification** | | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_FORWARD_EULER` | Use the forward Euler MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_ERK22b` | Use the ERK22b MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_ERK22a` | Use the ERK22a MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_RALSTON2` | Use the second order Ralston MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_MIS_MW3` | Use the Knoth-Wolke-3 MIS method. | +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_MRI_GARK_ERK33a` | Use the ERK33a MRI-GARK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_RALSTON3` | Use the third order Ralston MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_MRI_GARK_ERK45a` | Use the ERK45a MRI-GARK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_BACKWARD_EULER` | Use the backward Euler MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_MRI_GARK_IRK21a` | Use the IRK21a MRI-GARK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_MRI_GARK_IMPLICIT_MIDPOINT` | Use the implicit midpoint MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_MRI_GARK_ESDIRK34a` | Use the ESDIRK34a MRI-GARK method. | +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_MRI_GARK_ESDIRK46a` | Use the ESDIRK46a MRI-GARK method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_IMEX_MRI_GARK_EULER` | Use the Euler IMEX-MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL` | Use the trapezoidal rule IMEX-MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`ARKODE_IMEX_MRI_GARK_MIDPOINT` | Use the midpoint rule IMEX-MRI-GARK method. | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_IMEX_MRI_GARK3a` | Use the IMEX-MRI-GARK3a method. | +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_IMEX_MRI_GARK3b` | Use the IMEX-MRI-GARK3b method. | +-----------------------------------------------+------------------------------------------------------------+ | :index:`ARKODE_IMEX_MRI_GARK4` | Use the IMEX-MRI-GARK4 method. | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`MRISTEP_DEFAULT_EXPL_TABLE_1` | Use MRIStep's default 1st-order explicit method | + | | (MRI_GARK_FORWARD_EULER). | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`MRISTEP_DEFAULT_EXPL_TABLE_2` | Use MRIStep's default 2nd-order explicit method | + | | (MRI_GARK_ERK22b). | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`MRISTEP_DEFAULT_EXPL_TABLE_3` | Use MRIStep's default 3rd-order explicit method | | | (MIS_MW3). | +-----------------------------------------------+------------------------------------------------------------+ | :index:`MRISTEP_DEFAULT_EXPL_TABLE_4` | Use MRIStep's default 4th-order explicit method | | | (MRI_GARK_ERK45a). | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`MRISTEP_DEFAULT_IMPL_SD_TABLE_1` | Use MRIStep's default 1st-order solve-decoupled implicit | + | | method (MRI_GARK_BACKWARD_EULER). | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`MRISTEP_DEFAULT_IMPL_SD_TABLE_2` | Use MRIStep's default 2nd-order solve-decoupled implicit | | | method (MRI_GARK_IRK21a). | +-----------------------------------------------+------------------------------------------------------------+ @@ -332,6 +389,12 @@ contains the ARKODE output constants. | :index:`MRISTEP_DEFAULT_IMPL_SD_TABLE_4` | Use MRIStep's default 4th-order solve-decoupled implicit | | | method (MRI_GARK_ESDIRK46a). | +-----------------------------------------------+------------------------------------------------------------+ + | :index:`MRISTEP_DEFAULT_IMEX_SD_TABLE_1` | Use MRIStep's default 1st-order solve-decoupled ImEx | + | | method (IMEX_MRI_GARK_EULER). | + +-----------------------------------------------+------------------------------------------------------------+ + | :index:`MRISTEP_DEFAULT_IMEX_SD_TABLE_2` | Use MRIStep's default 2nd-order solve-decoupled ImEx | + | | method (ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL). | + +-----------------------------------------------+------------------------------------------------------------+ | :index:`MRISTEP_DEFAULT_IMEX_SD_TABLE_3` | Use MRIStep's default 3rd-order solve-decoupled ImEx | | | method (IMEX_MRI_GARK3b). | +-----------------------------------------------+------------------------------------------------------------+ diff --git a/doc/arkode/guide/source/Usage/MRIStep_c_interface/MRIStepCoupling.rst b/doc/arkode/guide/source/Usage/MRIStep_c_interface/MRIStepCoupling.rst index cb8a725b4d..a948bd6dbb 100644 --- a/doc/arkode/guide/source/Usage/MRIStep_c_interface/MRIStepCoupling.rst +++ b/doc/arkode/guide/source/Usage/MRIStep_c_interface/MRIStepCoupling.rst @@ -304,36 +304,46 @@ with values specified for each method below (e.g., ``ARKODE_MIS_KW3``). .. table:: Explicit MRI-GARK coupling tables. The default method for each order is marked with an asterisk (:math:`^*`). - ========================== =========== ===================== - Table name Order Reference - ========================== =========== ===================== - ``ARKODE_MIS_KW3`` :math:`3^*` :cite:p:`Schlegel:09` - ``ARKODE_MRI_GARK_ERK33a`` 3 :cite:p:`Sandu:19` - ``ARKODE_MRI_GARK_ERK45a`` :math:`4^*` :cite:p:`Sandu:19` - ========================== =========== ===================== + ================================= =========== ===================== + Table name Order Reference + ================================= =========== ===================== + ``ARKODE_MRI_GARK_FORWARD_EULER`` :math:`1^*` + ``ARKODE_MRI_GARK_ERK22b`` :math:`2^*` :cite:p:`Sandu:19` + ``ARKODE_MRI_GARK_ERK22a`` 2 :cite:p:`Sandu:19` + ``ARKODE_MRI_GARK_RALSTON2`` 2 :cite:p:`Roberts:22` + ``ARKODE_MIS_KW3`` :math:`3^*` :cite:p:`Schlegel:09` + ``ARKODE_MRI_GARK_ERK33a`` 3 :cite:p:`Sandu:19` + ``ARKODE_MRI_GARK_RALSTON3`` 3 :cite:p:`Roberts:22` + ``ARKODE_MRI_GARK_ERK45a`` :math:`4^*` :cite:p:`Sandu:19` + ================================= =========== ===================== .. table:: Diagonally-implicit, solve-decoupled MRI-GARK coupling tables. The default method for each order is marked with an asterisk (:math:`^*`). - ============================= =========== =============== ================== - Table name Order Implicit Solves Reference - ============================= =========== =============== ================== - ``ARKODE_MRI_GARK_IRK21a`` :math:`2^*` 1 :cite:p:`Sandu:19` - ``ARKODE_MRI_GARK_ESDIRK34a`` :math:`3^*` 3 :cite:p:`Sandu:19` - ``ARKODE_MRI_GARK_ESDIRK46a`` :math:`4^*` 5 :cite:p:`Sandu:19` - ============================= =========== =============== ================== + ===================================== =========== =============== ================== + Table name Order Implicit Solves Reference + ===================================== =========== =============== ================== + ``ARKODE_MRI_GARK_BACKWARD_EULER`` :math:`1^*` 1 + ``ARKODE_MRI_GARK_IRK21a`` :math:`2^*` 1 :cite:p:`Sandu:19` + ``ARKODE_MRI_GARK_IMPLICIT_MIDPOINT`` 2 2 + ``ARKODE_MRI_GARK_ESDIRK34a`` :math:`3^*` 3 :cite:p:`Sandu:19` + ``ARKODE_MRI_GARK_ESDIRK46a`` :math:`4^*` 5 :cite:p:`Sandu:19` + ===================================== =========== =============== ================== .. table:: Diagonally-implicit, solve-decoupled IMEX-MRI-GARK coupling tables. The default method for each order is marked with an asterisk (:math:`^*`). - =========================== =========== =============== =================== - Table name Order Implicit Solves Reference - =========================== =========== =============== =================== - ``ARKODE_IMEX_MRI_GARK3a`` :math:`3^*` 2 :cite:p:`ChiRen:21` - ``ARKODE_IMEX_MRI_GARK3b`` 3 2 :cite:p:`ChiRen:21` - ``ARKODE_IMEX_MRI_GARK4`` :math:`4^*` 5 :cite:p:`ChiRen:21` - =========================== =========== =============== =================== + ==================================== =========== =============== =================== + Table name Order Implicit Solves Reference + ==================================== =========== =============== =================== + ``ARKODE_IMEX_MRI_GARK_EULER`` :math:`1^*` 1 + ``ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL`` :math:`2^*` 1 + ``ARKODE_IMEX_MRI_GARK_MIDPOINT`` 2 2 + ``ARKODE_IMEX_MRI_GARK3a`` :math:`3^*` 2 :cite:p:`ChiRen:21` + ``ARKODE_IMEX_MRI_GARK3b`` 3 2 :cite:p:`ChiRen:21` + ``ARKODE_IMEX_MRI_GARK4`` :math:`4^*` 5 :cite:p:`ChiRen:21` + ==================================== =========== =============== =================== diff --git a/doc/shared/RecentChanges.rst b/doc/shared/RecentChanges.rst index 9212d1ff37..ef560561bf 100644 --- a/doc/shared/RecentChanges.rst +++ b/doc/shared/RecentChanges.rst @@ -15,6 +15,26 @@ for users who wish to retain the current functionality. Added CMake infrastructure that enables externally maintained addons/plugins to be *optionally* built with SUNDIALS. See :ref:`Contributing` for details. +Added the following Runge-Kutta Butcher tables + +* ``ARKODE_FORWARD_EULER_1_1`` +* ``ARKODE_RALSTON_EULER_2_1_2`` +* ``ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2`` +* ``ARKODE_BACKWARD_EULER_1_1`` +* ``ARKODE_IMPLICIT_MIDPOINT_1_2`` +* ``ARKODE_IMPLICIT_TRAPEZOIDAL_2_2`` + +Added the following MRI coupling tables + +* ``ARKODE_MRI_GARK_FORWARD_EULER`` +* ``ARKODE_MRI_GARK_RALSTON2`` +* ``ARKODE_MRI_GARK_RALSTON3`` +* ``ARKODE_MRI_GARK_BACKWARD_EULER`` +* ``ARKODE_MRI_GARK_IMPLICIT_MIDPOINT`` +* ``ARKODE_IMEX_MRI_GARK_EULER`` +* ``ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL`` +* ``ARKODE_IMEX_MRI_GARK_MIDPOINT`` + **Bug Fixes** Updated the CMake variable ``HIP_PLATFORM`` default to ``amd`` as the previous diff --git a/doc/shared/figs/arkode/backward_euler_dirk_stab_region.png b/doc/shared/figs/arkode/backward_euler_dirk_stab_region.png new file mode 100644 index 0000000000..388c85f0e0 Binary files /dev/null and b/doc/shared/figs/arkode/backward_euler_dirk_stab_region.png differ diff --git a/doc/shared/figs/arkode/explicit_midpoint_euler_erk_stab_region.png b/doc/shared/figs/arkode/explicit_midpoint_euler_erk_stab_region.png new file mode 100644 index 0000000000..c3f7e7b937 Binary files /dev/null and b/doc/shared/figs/arkode/explicit_midpoint_euler_erk_stab_region.png differ diff --git a/doc/shared/figs/arkode/forward_euler_erk_stab_region.png b/doc/shared/figs/arkode/forward_euler_erk_stab_region.png new file mode 100644 index 0000000000..e5bedb74c7 Binary files /dev/null and b/doc/shared/figs/arkode/forward_euler_erk_stab_region.png differ diff --git a/doc/shared/figs/arkode/implicit_midpoint_dirk_stab_region.png b/doc/shared/figs/arkode/implicit_midpoint_dirk_stab_region.png new file mode 100644 index 0000000000..94695afff9 Binary files /dev/null and b/doc/shared/figs/arkode/implicit_midpoint_dirk_stab_region.png differ diff --git a/doc/shared/figs/arkode/implicit_trapezoidal_dirk_stab_region.png b/doc/shared/figs/arkode/implicit_trapezoidal_dirk_stab_region.png new file mode 100644 index 0000000000..da40013312 Binary files /dev/null and b/doc/shared/figs/arkode/implicit_trapezoidal_dirk_stab_region.png differ diff --git a/doc/shared/figs/arkode/ralston_euler_erk_stab_region.png b/doc/shared/figs/arkode/ralston_euler_erk_stab_region.png new file mode 100644 index 0000000000..e05a95f3c5 Binary files /dev/null and b/doc/shared/figs/arkode/ralston_euler_erk_stab_region.png differ diff --git a/doc/shared/sundials.bib b/doc/shared/sundials.bib index 3516b83111..b4b5cf0343 100644 --- a/doc/shared/sundials.bib +++ b/doc/shared/sundials.bib @@ -1795,6 +1795,15 @@ @article{DorPri:80 doi = {10.1016/0771-050X(80)90013-3} } +@book{Euler:68 + author = {Leonhard Euler}, + title = {Institutiones calculi integralis}, + volume = {Volumen Primum}, + year = {1768}, + publisher = {B. G. Teubner Verlag}, + note = {reprinted in Opera Omnia Series 1, Volume 11} +} + @article{HaSo:05, title={Explicit, time reversible, adaptive step size control}, author={Hairer, Ernst and S{\"o}derlind, Gustaf}, @@ -1943,6 +1952,42 @@ @article{Mclachlan:92 publisher = {IOP Publishing} } +@article{Ralston:62, + author = {Ralston, Anthony}, + title = {{Runge--Kutta} methods with minimum error bounds}, + journal = {Mathematics of Computation}, + volume = {16}, + number = {80}, + pages = {431–437}, + year = {1962}, + publisher = {American Mathematical Society}, + doi = {10.1090/s0025-5718-1962-0150954-0} +} + +@article{Roberts:22, + author = {Roberts, Steven and Popov, Andrey A and Sarshar, Arash and Sandu, Adrian}, + title = {A Fast Time-Stepping Strategy for Dynamical Systems Equipped with a Surrogate Model}, + journal = {SIAM Journal on Scientific Computing}, + volume = {44}, + number = {3}, + pages = {A1405--A1427}, + year = {2022}, + publisher = {SIAM}, + doi = {10.1137/20M1386281} +} + +@Article{Runge:95, + author = {Runge, C.}, + title = {Ueber die numerische Aufl{\"o}sung von Differentialgleichungen}, + journal = {Mathematische Annalen}, + volume = {46}, + number = {2}, + pages = {167-178}, + year = {1895}, + publisher = {Springer}, + doi = {10.1007/BF01446807} +} + @article{Sandu:19, author = {Sandu, A.}, title = {A Class of Multirate Infinitesimal GARK Methods}, diff --git a/include/arkode/arkode_arkstep.h b/include/arkode/arkode_arkstep.h index 8c43ee198b..f7a6f112f5 100644 --- a/include/arkode/arkode_arkstep.h +++ b/include/arkode/arkode_arkstep.h @@ -35,6 +35,7 @@ extern "C" { /* Default Butcher tables for each method/order */ /* explicit */ +static const int ARKSTEP_DEFAULT_ERK_1 = ARKODE_FORWARD_EULER_1_1; static const int ARKSTEP_DEFAULT_ERK_2 = ARKODE_HEUN_EULER_2_1_2; static const int ARKSTEP_DEFAULT_ERK_3 = ARKODE_BOGACKI_SHAMPINE_4_2_3; static const int ARKSTEP_DEFAULT_ERK_4 = ARKODE_ZONNEVELD_5_3_4; @@ -45,6 +46,7 @@ static const int ARKSTEP_DEFAULT_ERK_8 = ARKODE_FEHLBERG_13_7_8; static const int ARKSTEP_DEFAULT_ERK_9 = ARKODE_VERNER_16_8_9; /* implicit */ +static const int ARKSTEP_DEFAULT_DIRK_1 = ARKODE_BACKWARD_EULER_1_1; static const int ARKSTEP_DEFAULT_DIRK_2 = ARKODE_SDIRK_2_1_2; static const int ARKSTEP_DEFAULT_DIRK_3 = ARKODE_ARK324L2SA_DIRK_4_2_3; static const int ARKSTEP_DEFAULT_DIRK_4 = ARKODE_SDIRK_5_3_4; diff --git a/include/arkode/arkode_butcher_dirk.h b/include/arkode/arkode_butcher_dirk.h index fa60103e89..ffa4908c91 100644 --- a/include/arkode/arkode_butcher_dirk.h +++ b/include/arkode/arkode_butcher_dirk.h @@ -51,7 +51,10 @@ typedef enum ARKODE_ESDIRK547L2SA_7_4_5, ARKODE_ESDIRK547L2SA2_7_4_5, ARKODE_ARK2_DIRK_3_1_2, - ARKODE_MAX_DIRK_NUM = ARKODE_ARK2_DIRK_3_1_2 + ARKODE_BACKWARD_EULER_1_1, + ARKODE_IMPLICIT_MIDPOINT_1_2, + ARKODE_IMPLICIT_TRAPEZOIDAL_2_2, + ARKODE_MAX_DIRK_NUM = ARKODE_IMPLICIT_TRAPEZOIDAL_2_2 } ARKODE_DIRKTableID; /* Accessor routine to load built-in DIRK table */ diff --git a/include/arkode/arkode_butcher_erk.h b/include/arkode/arkode_butcher_erk.h index bb99368b41..764b67bbfe 100644 --- a/include/arkode/arkode_butcher_erk.h +++ b/include/arkode/arkode_butcher_erk.h @@ -49,7 +49,10 @@ typedef enum ARKODE_VERNER_10_6_7, ARKODE_VERNER_13_7_8, ARKODE_VERNER_16_8_9, - ARKODE_MAX_ERK_NUM = ARKODE_VERNER_16_8_9 + ARKODE_FORWARD_EULER_1_1, + ARKODE_RALSTON_EULER_2_1_2, + ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2, + ARKODE_MAX_ERK_NUM = ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2 } ARKODE_ERKTableID; /* Accessor routine to load built-in ERK table */ diff --git a/include/arkode/arkode_erkstep.h b/include/arkode/arkode_erkstep.h index 73f2a474c6..ed93a240b9 100644 --- a/include/arkode/arkode_erkstep.h +++ b/include/arkode/arkode_erkstep.h @@ -32,6 +32,7 @@ extern "C" { /* Default Butcher tables for each order */ +static const int ERKSTEP_DEFAULT_1 = ARKODE_FORWARD_EULER_1_1; static const int ERKSTEP_DEFAULT_2 = ARKODE_HEUN_EULER_2_1_2; static const int ERKSTEP_DEFAULT_3 = ARKODE_BOGACKI_SHAMPINE_4_2_3; static const int ERKSTEP_DEFAULT_4 = ARKODE_ZONNEVELD_5_3_4; diff --git a/include/arkode/arkode_mristep.h b/include/arkode/arkode_mristep.h index 1191e0a39d..000bae135a 100644 --- a/include/arkode/arkode_mristep.h +++ b/include/arkode/arkode_mristep.h @@ -52,17 +52,32 @@ typedef enum ARKODE_IMEX_MRI_GARK3a, ARKODE_IMEX_MRI_GARK3b, ARKODE_IMEX_MRI_GARK4, - ARKODE_MAX_MRI_NUM = ARKODE_IMEX_MRI_GARK4 + ARKODE_MRI_GARK_FORWARD_EULER, + ARKODE_MRI_GARK_RALSTON2, + ARKODE_MRI_GARK_ERK22a, + ARKODE_MRI_GARK_ERK22b, + ARKODE_MRI_GARK_RALSTON3, + ARKODE_MRI_GARK_BACKWARD_EULER, + ARKODE_MRI_GARK_IMPLICIT_MIDPOINT, + ARKODE_IMEX_MRI_GARK_EULER, + ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL, + ARKODE_IMEX_MRI_GARK_MIDPOINT, + ARKODE_MAX_MRI_NUM = ARKODE_IMEX_MRI_GARK_MIDPOINT, } ARKODE_MRITableID; /* Default MRI coupling tables for each order */ +static const int MRISTEP_DEFAULT_EXPL_1 = ARKODE_MRI_GARK_FORWARD_EULER; +static const int MRISTEP_DEFAULT_EXPL_2 = ARKODE_MRI_GARK_ERK22b; +static const int MRISTEP_DEFAULT_EXPL_3 = ARKODE_MIS_KW3; +static const int MRISTEP_DEFAULT_EXPL_4 = ARKODE_MRI_GARK_ERK45a; -static const int MRISTEP_DEFAULT_3 = ARKODE_MIS_KW3; -static const int MRISTEP_DEFAULT_EXPL_3 = ARKODE_MIS_KW3; -static const int MRISTEP_DEFAULT_EXPL_4 = ARKODE_MRI_GARK_ERK45a; +static const int MRISTEP_DEFAULT_IMPL_SD_1 = ARKODE_MRI_GARK_BACKWARD_EULER; static const int MRISTEP_DEFAULT_IMPL_SD_2 = ARKODE_MRI_GARK_IRK21a; static const int MRISTEP_DEFAULT_IMPL_SD_3 = ARKODE_MRI_GARK_ESDIRK34a; static const int MRISTEP_DEFAULT_IMPL_SD_4 = ARKODE_MRI_GARK_ESDIRK46a; + +static const int MRISTEP_DEFAULT_IMEX_SD_1 = ARKODE_IMEX_MRI_GARK_EULER; +static const int MRISTEP_DEFAULT_IMEX_SD_2 = ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL; static const int MRISTEP_DEFAULT_IMEX_SD_3 = ARKODE_IMEX_MRI_GARK3b; static const int MRISTEP_DEFAULT_IMEX_SD_4 = ARKODE_IMEX_MRI_GARK4; diff --git a/src/arkode/arkode_arkstep.c b/src/arkode/arkode_arkstep.c index 5e181982c1..70275484ec 100644 --- a/src/arkode/arkode_arkstep.c +++ b/src/arkode/arkode_arkstep.c @@ -2159,6 +2159,7 @@ int arkStep_SetButcherTables(ARKodeMem ark_mem) { switch (step_mem->q) { + case (1): itable = ARKSTEP_DEFAULT_DIRK_1; break; case (2): itable = ARKSTEP_DEFAULT_DIRK_2; break; case (3): itable = ARKSTEP_DEFAULT_DIRK_3; break; case (4): itable = ARKSTEP_DEFAULT_DIRK_4; break; @@ -2176,6 +2177,7 @@ int arkStep_SetButcherTables(ARKodeMem ark_mem) { switch (step_mem->q) { + case (1): itable = ARKSTEP_DEFAULT_ERK_1; break; case (2): etable = ARKSTEP_DEFAULT_ERK_2; break; case (3): etable = ARKSTEP_DEFAULT_ERK_3; break; case (4): etable = ARKSTEP_DEFAULT_ERK_4; break; diff --git a/src/arkode/arkode_butcher_dirk.def b/src/arkode/arkode_butcher_dirk.def index 2ce5444f8c..40d13e60bd 100644 --- a/src/arkode/arkode_butcher_dirk.def +++ b/src/arkode/arkode_butcher_dirk.def @@ -34,7 +34,11 @@ imeth type A-stable L-stable QP ----------------------------------------------------------------- + ARKODE_BACKWARD_EULER_1_1 SDIRK Y Y Y ARKODE_SDIRK_2_1_2 SDIRK Y N Y + ARKODE_ARK2_DIRK_3_1_2 ESDIRK Y Y Y + ARKODE_IMPLICIT_MIDPOINT_1_2 SDIRK Y N Y + ARKODE_IMPLICIT_TRAPEZOIDAL_2_2 ESDIRK Y N Y ARKODE_BILLINGTON_3_3_2 SDIRK N N N ARKODE_TRBDF2_3_3_2 ESDIRK N N Y ARKODE_KVAERNO_4_2_3 ESDIRK Y Y N @@ -57,7 +61,6 @@ ARKODE_ARK548L2SAb_DIRK_8_4_5* ESDIRK Y Y N ARKODE_ESDIRK547L2SA_7_4_5 ESDIRK Y Y N ARKODE_ESDIRK547L2SA2_7_4_5 ESDIRK Y Y N - ARKODE_ARK2_DIRK_3_1_2 ESDIRK Y Y Y ----------------------------------------------------------------- */ @@ -65,6 +68,19 @@ ARK_BUTCHER_TABLE(ARKODE_DIRK_NONE, { return NULL; }) +ARK_BUTCHER_TABLE(ARKODE_BACKWARD_EULER_1_1, { /* Backward Euler (L,B stable) */ + ARKodeButcherTable B = ARKodeButcherTable_Alloc(1, SUNFALSE); + B->q = 1; + B->p = 0; + + B->A[0][0] = SUN_RCONST(1.0); + + B->b[0] = SUN_RCONST(1.0); + + B->c[0] = SUN_RCONST(1.0); + return B; + }) + ARK_BUTCHER_TABLE(ARKODE_SDIRK_2_1_2, { /* SDIRK-2-1 (A,B stable) */ ARKodeButcherTable B = ARKodeButcherTable_Alloc(2, SUNTRUE); B->q = 2; @@ -119,6 +135,34 @@ ARK_BUTCHER_TABLE(ARKODE_ARK2_DIRK_3_1_2, { /* ARK2 Implicit Table (A,L stable) return B; }) +ARK_BUTCHER_TABLE(ARKODE_IMPLICIT_MIDPOINT_1_2, { /* Implicit Midpoint Rule (A,B stable) */ + ARKodeButcherTable B = ARKodeButcherTable_Alloc(1, SUNFALSE); + B->q = 2; + B->p = 0; + + B->A[0][0] = SUN_RCONST(0.5); + + B->b[0] = SUN_RCONST(1.0); + + B->c[0] = SUN_RCONST(0.5); + return B; + }) + +ARK_BUTCHER_TABLE(ARKODE_IMPLICIT_TRAPEZOIDAL_2_2, { /* Implicit Trapezoidal Rule (A stable) */ + ARKodeButcherTable B = ARKodeButcherTable_Alloc(2, SUNFALSE); + B->q = 2; + B->p = 0; + + B->A[1][0] = SUN_RCONST(0.5); + B->A[1][1] = SUN_RCONST(0.5); + + B->b[0] = SUN_RCONST(0.5); + B->b[1] = SUN_RCONST(0.5); + + B->c[1] = SUN_RCONST(1.0); + return B; + }) + ARK_BUTCHER_TABLE(ARKODE_BILLINGTON_3_3_2, { /* Billington-SDIRK */ ARKodeButcherTable B = ARKodeButcherTable_Alloc(3, SUNTRUE); diff --git a/src/arkode/arkode_butcher_erk.def b/src/arkode/arkode_butcher_erk.def index 5fe439bbab..a29292e27d 100644 --- a/src/arkode/arkode_butcher_erk.def +++ b/src/arkode/arkode_butcher_erk.def @@ -36,38 +36,50 @@ are known precisely enough for use in quad precision (128-bit) calculations. - imeth QP - -------------------------------------- - ARKODE_HEUN_EULER_2_1_2 Y - ARKODE_BOGACKI_SHAMPINE_4_2_3 Y - ARKODE_ARK324L2SA_ERK_4_2_3* N - ARKODE_SHU_OSHER_3_2_3 Y - ARKODE_SOFRONIOU_SPALETTA_5_3_4 Y - ARKODE_ZONNEVELD_5_3_4 Y - ARKODE_ARK436L2SA_ERK_6_3_4* N - ARKODE_ARK437L2SA_ERK_7_3_4* N - ARKODE_SAYFY_ABURUB_6_3_4 N - ARKODE_CASH_KARP_6_4_5 Y - ARKODE_FEHLBERG_6_4_5 Y - ARKODE_DORMAND_PRINCE_7_4_5 Y - ARKODE_ARK548L2SA_ERK_8_4_5* N - ARKODE_ARK548L2SAb_ERK_8_4_5* N - ARKODE_VERNER_8_5_6 Y - ARKODE_FEHLBERG_13_7_8 Y - ARKODE_ARK2_ERK_3_1_2 Y - ARKODE_VERNER_9_5_6 Y - ARKODE_VERNER_10_6_7 Y - ARKODE_VERNER_13_7_8 Y - ARKODE_VERNER_16_8_9 Y - -------------------------------------- - ARKODE_KNOTH_WOLKE_3_3^ Y - -------------------------------------- + imeth QP + --------------------------------------- + ARKODE_FORWARD_EULER_1_1 Y + ARKODE_HEUN_EULER_2_1_2 Y + ARKODE_RALSTON_EULER_2_1_2 Y + ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2 Y + ARKODE_BOGACKI_SHAMPINE_4_2_3 Y + ARKODE_ARK324L2SA_ERK_4_2_3* N + ARKODE_SHU_OSHER_3_2_3 Y + ARKODE_SOFRONIOU_SPALETTA_5_3_4 Y + ARKODE_ZONNEVELD_5_3_4 Y + ARKODE_ARK436L2SA_ERK_6_3_4* N + ARKODE_ARK437L2SA_ERK_7_3_4* N + ARKODE_SAYFY_ABURUB_6_3_4 N + ARKODE_CASH_KARP_6_4_5 Y + ARKODE_FEHLBERG_6_4_5 Y + ARKODE_DORMAND_PRINCE_7_4_5 Y + ARKODE_ARK548L2SA_ERK_8_4_5* N + ARKODE_ARK548L2SAb_ERK_8_4_5* N + ARKODE_VERNER_8_5_6 Y + ARKODE_FEHLBERG_13_7_8 Y + ARKODE_ARK2_ERK_3_1_2 Y + ARKODE_VERNER_9_5_6 Y + ARKODE_VERNER_10_6_7 Y + ARKODE_VERNER_13_7_8 Y + ARKODE_VERNER_16_8_9 Y + --------------------------------------- + ARKODE_KNOTH_WOLKE_3_3^ Y + --------------------------------------- */ ARK_BUTCHER_TABLE(ARKODE_ERK_NONE, { return NULL; }) +ARK_BUTCHER_TABLE(ARKODE_FORWARD_EULER_1_1, { /* Euler-ERK */ + ARKodeButcherTable B = ARKodeButcherTable_Alloc(1, SUNFALSE); + B->q = 1; + B->p = 0; + + B->b[0] = SUN_RCONST(1.0); + return B; + }) + ARK_BUTCHER_TABLE(ARKODE_HEUN_EULER_2_1_2, { /* Heun-Euler-ERK */ ARKodeButcherTable B = ARKodeButcherTable_Alloc(2, SUNTRUE); B->q = 2; @@ -84,6 +96,37 @@ ARK_BUTCHER_TABLE(ARKODE_HEUN_EULER_2_1_2, { /* Heun-Euler-ERK */ return B; }) +ARK_BUTCHER_TABLE(ARKODE_RALSTON_EULER_2_1_2, { /* Ralston-Euler-ERK */ + ARKodeButcherTable B = ARKodeButcherTable_Alloc(2, SUNTRUE); + B->q = 2; + B->p = 1; + + B->A[1][0] = SUN_RCONST(2.0) / SUN_RCONST(3.0); + + B->b[0] = SUN_RCONST(1.0) / SUN_RCONST(4.0); + B->b[1] = SUN_RCONST(3.0) / SUN_RCONST(4.0); + + B->d[0] = SUN_RCONST(1.0); + + B->c[1] = SUN_RCONST(2.0) / SUN_RCONST(3.0); + return B; + }) + +ARK_BUTCHER_TABLE(ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2, { /* Explicit-Midpoint-Euler-ERK */ + ARKodeButcherTable B = ARKodeButcherTable_Alloc(2, SUNTRUE); + B->q = 2; + B->p = 1; + + B->A[1][0] = SUN_RCONST(1.0) / SUN_RCONST(2.0); + + B->b[1] = SUN_RCONST(1.0); + + B->d[0] = SUN_RCONST(1.0); + + B->c[1] = SUN_RCONST(1.0) / SUN_RCONST(2.0); + return B; + }) + ARK_BUTCHER_TABLE(ARKODE_ARK2_ERK_3_1_2, { /* ARK2 Explicit Table */ ARKodeButcherTable B = ARKodeButcherTable_Alloc(3, SUNTRUE); diff --git a/src/arkode/arkode_erkstep.c b/src/arkode/arkode_erkstep.c index adf22daadb..2432784596 100644 --- a/src/arkode/arkode_erkstep.c +++ b/src/arkode/arkode_erkstep.c @@ -838,6 +838,7 @@ int erkStep_SetButcherTable(ARKodeMem ark_mem) /* select method based on order */ switch (step_mem->q) { + case (1): etable = ERKSTEP_DEFAULT_1; break; case (2): etable = ERKSTEP_DEFAULT_2; break; case (3): etable = ERKSTEP_DEFAULT_3; break; case (4): etable = ERKSTEP_DEFAULT_4; break; diff --git a/src/arkode/arkode_mri_tables.def b/src/arkode/arkode_mri_tables.def index 13c5635b43..999ef71f95 100644 --- a/src/arkode/arkode_mri_tables.def +++ b/src/arkode/arkode_mri_tables.def @@ -31,24 +31,49 @@ are known precisely enough for use in quad precision (128-bit) calculations. - imeth order type QP - ------------------------------------------------ - ARKODE_MIS_KW3 3 E Y - ARKODE_MRI_GARK_ERK33a 3 E Y - ARKODE_MRI_GARK_ERK45a 4 E Y - ARKODE_MRI_GARK_IRK21a 2 ID Y - ARKODE_MRI_GARK_ESDIRK34a 3 ID Y - ARKODE_MRI_GARK_ESDIRK46a 4 ID Y - ARKODE_IMEX_MRI_GARK3a 3 ID Y - ARKODE_IMEX_MRI_GARK3b 3 ID Y - ARKODE_IMEX_MRI_GARK4 4 ID Y - ------------------------------------------------ + imeth order type QP + ----------------------------------------------------- + ARKODE_MRI_GARK_FORWARD_EULER 1 E Y + ARKODE_MRI_GARK_RALSTON2 2 E Y + ARKODE_MIS_KW3 3 E Y + ARKODE_MRI_GARK_ERK22a 2 E Y + ARKODE_MRI_GARK_ERK22b 2 E Y + ARKODE_MRI_GARK_ERK33a 3 E Y + ARKODE_MRI_GARK_RALSTON3 3 E Y + ARKODE_MRI_GARK_ERK45a 4 E Y + ARKODE_MRI_GARK_BACKWARD_EULER 1 ID Y + ARKODE_MRI_GARK_IRK21a 2 ID Y + ARKODE_MRI_GARK_IMPLICIT_MIDPOINT 2 ID Y + ARKODE_MRI_GARK_ESDIRK34a 3 ID Y + ARKODE_MRI_GARK_ESDIRK46a 4 ID Y + ARKODE_IMEX_MRI_GARK_EULER 1 ID Y + ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL 2 ID Y + ARKODE_IMEX_MRI_GARK_MIDPOINT 2 ID Y + ARKODE_IMEX_MRI_GARK3a 3 ID Y + ARKODE_IMEX_MRI_GARK3b 3 ID Y + ARKODE_IMEX_MRI_GARK4 4 ID Y + ----------------------------------------------------- */ + ARK_MRI_TABLE(ARKODE_MRI_NONE, { return NULL; }) +ARK_MRI_TABLE(ARKODE_MRI_GARK_FORWARD_EULER, { + ARKodeButcherTable B = ARKodeButcherTable_LoadERK(ARKODE_FORWARD_EULER_1_1); + MRIStepCoupling C = MRIStepCoupling_MIStoMRI(B, 1, 0); + ARKodeButcherTable_Free(B); + return C; + }) + +ARK_MRI_TABLE(ARKODE_MRI_GARK_RALSTON2, { /* Roberts et al., SISC 44:A1405 - A1427, 2022 */ + ARKodeButcherTable B = ARKodeButcherTable_LoadERK(ARKODE_RALSTON_EULER_2_1_2); + MRIStepCoupling C = MRIStepCoupling_MIStoMRI(B, 2, 0); + ARKodeButcherTable_Free(B); + return C; + }) + ARK_MRI_TABLE(ARKODE_MIS_KW3, { /* Schlegel et al., JCAM 226:345-357, 2009 */ ARKodeButcherTable B = ARKodeButcherTable_LoadERK(ARKODE_KNOTH_WOLKE_3_3); MRIStepCoupling C = MRIStepCoupling_MIStoMRI(B, 3, 0); @@ -56,6 +81,20 @@ ARK_MRI_TABLE(ARKODE_MIS_KW3, { /* Schlegel et al., JCAM 226:345-357, 2009 */ return C; }) +ARK_MRI_TABLE(ARKODE_MRI_GARK_ERK22a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ + ARKodeButcherTable B = ARKodeButcherTable_LoadERK(ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2); + MRIStepCoupling C = MRIStepCoupling_MIStoMRI(B, 2, 0); + ARKodeButcherTable_Free(B); + return C; + }) + +ARK_MRI_TABLE(ARKODE_MRI_GARK_ERK22b, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ + ARKodeButcherTable B = ARKodeButcherTable_LoadERK(ARKODE_HEUN_EULER_2_1_2); + MRIStepCoupling C = MRIStepCoupling_MIStoMRI(B, 2, 0); + ARKodeButcherTable_Free(B); + return C; + }) + ARK_MRI_TABLE(ARKODE_MRI_GARK_ERK33a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ MRIStepCoupling C = MRIStepCoupling_Alloc(2, 4, MRISTEP_EXPLICIT); @@ -77,6 +116,31 @@ ARK_MRI_TABLE(ARKODE_MRI_GARK_ERK33a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ return C; }) +ARK_MRI_TABLE(ARKODE_MRI_GARK_RALSTON3, { /* Roberts et al., SISC 44:A1405 - A1427, 2022 */ + MRIStepCoupling C = MRIStepCoupling_Alloc(2, 4, MRISTEP_EXPLICIT); + + C->q = 3; + C->p = 0; + + C->c[1] = ONE/TWO; + C->c[2] = SUN_RCONST(3.0)/SUN_RCONST(4.0); + C->c[3] = ONE; + + C->W[0][1][0] = ONE/TWO; + C->W[0][2][0] = -SUN_RCONST(11.0)/SUN_RCONST(4.0); + C->W[0][2][1] = SUN_RCONST(3.0); + C->W[0][3][0] = SUN_RCONST(47.0)/SUN_RCONST(36.0); + C->W[0][3][1] = -ONE/SUN_RCONST(6.0); + C->W[0][3][2] = -SUN_RCONST(8.0)/SUN_RCONST(9.0); + + C->W[1][2][0] = SUN_RCONST(9.0)/TWO; + C->W[1][2][1] = -SUN_RCONST(9.0)/TWO; + C->W[1][3][0] = -SUN_RCONST(13.0)/SUN_RCONST(6.0); + C->W[1][3][1] = -ONE/TWO; + C->W[1][3][2] = SUN_RCONST(8.0)/SUN_RCONST(3.0); + return C; + }) + ARK_MRI_TABLE(ARKODE_MRI_GARK_ERK45a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ MRIStepCoupling C = MRIStepCoupling_Alloc(2, 6, MRISTEP_EXPLICIT); @@ -122,6 +186,21 @@ ARK_MRI_TABLE(ARKODE_MRI_GARK_ERK45a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ return C; }) +ARK_MRI_TABLE(ARKODE_MRI_GARK_BACKWARD_EULER, { + MRIStepCoupling C = MRIStepCoupling_Alloc(1, 3, MRISTEP_IMPLICIT); + + C->q = 1; + C->p = 0; + + C->c[1] = ONE; + C->c[2] = ONE; + + C->G[0][1][0] = ONE; + C->G[0][2][0] = -ONE; + C->G[0][2][2] = ONE; + return C; + }) + ARK_MRI_TABLE(ARKODE_MRI_GARK_IRK21a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ MRIStepCoupling C; ARKodeButcherTable B = ARKodeButcherTable_Alloc(3, SUNFALSE); @@ -143,6 +222,23 @@ ARK_MRI_TABLE(ARKODE_MRI_GARK_IRK21a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ return C; }) +ARK_MRI_TABLE(ARKODE_MRI_GARK_IMPLICIT_MIDPOINT, { + MRIStepCoupling C = MRIStepCoupling_Alloc(1, 4, MRISTEP_IMPLICIT); + + C->q = 2; + C->p = 0; + + C->c[1] = ONE/TWO; + C->c[2] = ONE/TWO; + C->c[3] = ONE; + + C->G[0][1][0] = ONE/TWO; + C->G[0][2][0] = -ONE/TWO; + C->G[0][2][2] = ONE/TWO; + C->G[0][3][2] = ONE/TWO; + return C; + }) + ARK_MRI_TABLE(ARKODE_MRI_GARK_ESDIRK34a, { /* A. Sandu, SINUM 57:2300-2327, 2019 */ MRIStepCoupling C = MRIStepCoupling_Alloc(1, 7, MRISTEP_IMPLICIT); sunrealtype beta = SUN_RCONST(0.4358665215084589994160194511935568425); @@ -257,6 +353,64 @@ ARK_MRI_TABLE(ARKODE_MRI_GARK_ESDIRK46a, { /* A. Sandu, SINUM 57:2300-2327, 2019 return C; }) +ARK_MRI_TABLE(ARKODE_IMEX_MRI_GARK_EULER, { + MRIStepCoupling C = MRIStepCoupling_Alloc(1, 3, MRISTEP_IMEX); + + C->q = 1; + C->p = 0; + + C->c[1] = ONE; + C->c[2] = ONE; + + C->W[0][1][0] = ONE; + + C->G[0][1][0] = ONE; + C->G[0][2][0] = -ONE; + C->G[0][2][2] = ONE; + return C; + }) + +ARK_MRI_TABLE(ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL, { + MRIStepCoupling C = MRIStepCoupling_Alloc(1, 4, MRISTEP_IMEX); + + C->q = 2; + C->p = 0; + + C->c[1] = ONE; + C->c[2] = ONE; + C->c[3] = ONE; + + C->W[0][1][0] = ONE; + C->W[0][3][0] = -ONE/TWO; + C->W[0][3][2] = ONE/TWO; + + C->G[0][1][0] = ONE; + C->G[0][2][0] = -ONE/TWO; + C->G[0][2][2] = ONE/TWO; + return C; + }) + +ARK_MRI_TABLE(ARKODE_IMEX_MRI_GARK_MIDPOINT, { + MRIStepCoupling C = MRIStepCoupling_Alloc(1, 4, MRISTEP_IMEX); + + C->q = 2; + C->p = 0; + + C->c[1] = ONE/TWO; + C->c[2] = ONE/TWO; + C->c[3] = ONE; + + C->W[0][1][0] = ONE/TWO; + C->W[0][3][0] = -ONE/TWO; + C->W[0][3][2] = ONE; + + C->G[0][1][0] = ONE/TWO; + C->G[0][2][0] = -ONE/TWO; + C->G[0][2][2] = ONE/TWO; + C->G[0][3][2] = ONE/TWO; + return C; + }) + ARK_MRI_TABLE(ARKODE_IMEX_MRI_GARK3a, { /* R. Chinomona & D. Reynolds SINUM 43(5):A3082-A3113, 2021 */ MRIStepCoupling C = MRIStepCoupling_Alloc(1, 8, MRISTEP_IMEX); sunrealtype beta = SUN_RCONST(0.4358665215084589994160194511935568425); diff --git a/src/arkode/arkode_mristep.c b/src/arkode/arkode_mristep.c index cea390373e..78f34da788 100644 --- a/src/arkode/arkode_mristep.c +++ b/src/arkode/arkode_mristep.c @@ -90,8 +90,7 @@ void* MRIStepCreate(ARKRhsFn fse, ARKRhsFn fsi, sunrealtype t0, N_Vector y0, } /* Allocate ARKodeMRIStepMem structure, and initialize to zero */ - step_mem = NULL; - step_mem = (ARKodeMRIStepMem)malloc(sizeof(struct ARKodeMRIStepMemRec)); + step_mem = (ARKodeMRIStepMem)calloc(1, sizeof(*step_mem)); if (step_mem == NULL) { arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, @@ -99,7 +98,6 @@ void* MRIStepCreate(ARKRhsFn fse, ARKRhsFn fsi, sunrealtype t0, N_Vector y0, ARKodeFree((void**)&ark_mem); return (NULL); } - memset(step_mem, 0, sizeof(struct ARKodeMRIStepMemRec)); /* Attach step_mem structure and function pointers to ark_mem */ ark_mem->step_attachlinsol = mriStep_AttachLinsol; @@ -917,23 +915,23 @@ int mriStep_Init(ARKodeMem ark_mem, int init_type) return (ARK_ILL_INPUT); } - /* Retrieve/store method and embedding orders now that tables are finalized */ - step_mem->stages = step_mem->MRIC->stages; - step_mem->q = step_mem->MRIC->q; - step_mem->p = step_mem->MRIC->p; - /* allocate/fill derived quantities from MRIC structure */ /* stage map */ if (step_mem->stage_map) { free(step_mem->stage_map); - step_mem->stage_map = NULL; ark_mem->liw -= step_mem->stages; } - step_mem->stage_map = (int*)calloc(step_mem->stages, sizeof(int)); - ark_mem->liw += step_mem->stages; - + step_mem->stage_map = (int*)calloc(step_mem->MRIC->stages, + sizeof(*step_mem->stage_map)); + if (step_mem->stage_map == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + MSG_ARK_MEM_FAIL); + return (ARK_MEM_FAIL); + } + ark_mem->liw += step_mem->MRIC->stages; retval = mriStepCoupling_GetStageMap(step_mem->MRIC, step_mem->stage_map, &(step_mem->nstages_active)); if (retval != ARK_SUCCESS) @@ -947,12 +945,18 @@ int mriStep_Init(ARKodeMem ark_mem, int init_type) if (step_mem->stagetypes) { free(step_mem->stagetypes); - step_mem->stagetypes = NULL; ark_mem->liw -= step_mem->stages; } - step_mem->stagetypes = (int*)calloc(step_mem->stages, sizeof(int)); - ark_mem->liw += step_mem->stages; - for (j = 0; j < step_mem->stages; j++) + step_mem->stagetypes = (int*)calloc(step_mem->MRIC->stages, + sizeof(*step_mem->stagetypes)); + if (step_mem->stagetypes == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + MSG_ARK_MEM_FAIL); + return (ARK_MEM_FAIL); + } + ark_mem->liw += step_mem->MRIC->stages; + for (j = 0; j < step_mem->MRIC->stages; j++) { step_mem->stagetypes[j] = mriStepCoupling_GetStageType(step_mem->MRIC, j); } @@ -961,23 +965,70 @@ int mriStep_Init(ARKodeMem ark_mem, int init_type) if (step_mem->Ae_row) { free(step_mem->Ae_row); - step_mem->Ae_row = NULL; ark_mem->lrw -= step_mem->stages; } - step_mem->Ae_row = (sunrealtype*)calloc(step_mem->stages, - sizeof(sunrealtype)); - ark_mem->lrw += step_mem->stages; + step_mem->Ae_row = (sunrealtype*)calloc(step_mem->MRIC->stages, + sizeof(*step_mem->Ae_row)); + if (step_mem->Ae_row == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + MSG_ARK_MEM_FAIL); + return (ARK_MEM_FAIL); + } + ark_mem->lrw += step_mem->MRIC->stages; /* implicit RK coefficient row */ if (step_mem->Ai_row) { free(step_mem->Ai_row); - step_mem->Ai_row = NULL; ark_mem->lrw -= step_mem->stages; } - step_mem->Ai_row = (sunrealtype*)calloc(step_mem->stages, - sizeof(sunrealtype)); - ark_mem->lrw += step_mem->stages; + step_mem->Ai_row = (sunrealtype*)calloc(step_mem->MRIC->stages, + sizeof(*step_mem->Ai_row)); + if (step_mem->Ai_row == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + MSG_ARK_MEM_FAIL); + return (ARK_MEM_FAIL); + } + ark_mem->lrw += step_mem->MRIC->stages; + + /* Allocate reusable arrays for fused vector interface */ + if (step_mem->cvals) + { + free(step_mem->cvals); + ark_mem->lrw -= step_mem->nfusedopvecs; + } + if (step_mem->Xvecs) + { + free(step_mem->Xvecs); + ark_mem->liw -= step_mem->nfusedopvecs; + } + step_mem->nfusedopvecs = 2 * step_mem->MRIC->stages + 2; + step_mem->cvals = (sunrealtype*)calloc(step_mem->nfusedopvecs, + sizeof(*step_mem->cvals)); + if (step_mem->cvals == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + MSG_ARK_MEM_FAIL); + return (ARK_MEM_FAIL); + } + ark_mem->lrw += step_mem->nfusedopvecs; + + step_mem->Xvecs = (N_Vector*)calloc(step_mem->nfusedopvecs, + sizeof(*step_mem->Xvecs)); + if (step_mem->Xvecs == NULL) + { + arkProcessError(ark_mem, ARK_MEM_FAIL, __LINE__, __func__, __FILE__, + MSG_ARK_MEM_FAIL); + return (ARK_MEM_FAIL); + } + ark_mem->liw += step_mem->nfusedopvecs; + + /* Retrieve/store method and embedding orders now that tables are finalized */ + step_mem->stages = step_mem->MRIC->stages; + step_mem->q = step_mem->MRIC->q; + step_mem->p = step_mem->MRIC->p; /* Allocate MRI RHS vector memory, update storage requirements */ /* Allocate Fse[0] ... Fse[nstages_active - 1] and */ @@ -1054,23 +1105,6 @@ int mriStep_Init(ARKodeMem ark_mem, int init_type) step_mem->lmem = NULL; } - /* Allocate reusable arrays for fused vector interface */ - step_mem->nfusedopvecs = 2 * step_mem->stages + 2; - if (step_mem->cvals == NULL) - { - step_mem->cvals = (sunrealtype*)calloc(step_mem->nfusedopvecs, - sizeof(sunrealtype)); - if (step_mem->cvals == NULL) { return (ARK_MEM_FAIL); } - ark_mem->lrw += (step_mem->nfusedopvecs); - } - if (step_mem->Xvecs == NULL) - { - step_mem->Xvecs = (N_Vector*)calloc(step_mem->nfusedopvecs, - sizeof(N_Vector)); - if (step_mem->Xvecs == NULL) { return (ARK_MEM_FAIL); } - ark_mem->liw += (step_mem->nfusedopvecs); /* pointers */ - } - /* Allocate inner stepper data */ retval = mriStepInnerStepper_AllocVecs(step_mem->stepper, step_mem->MRIC->nmat, ark_mem->ewt); @@ -1693,6 +1727,8 @@ int mriStep_SetCoupling(ARKodeMem ark_mem) { ARKodeMRIStepMem step_mem; sunindextype Cliw, Clrw; + int q_actual; + ARKODE_MRITableID table_id; /* access ARKodeMRIStepMem structure */ if (ark_mem->step_mem == NULL) @@ -1702,72 +1738,57 @@ int mriStep_SetCoupling(ARKodeMem ark_mem) return (ARK_MEM_NULL); } step_mem = (ARKodeMRIStepMem)ark_mem->step_mem; + q_actual = step_mem->q; /* if coupling has already been specified, just return */ if (step_mem->MRIC != NULL) { return (ARK_SUCCESS); } + if (q_actual < 1 || q_actual > 4) + { + arkProcessError(ark_mem, ARK_ILL_INPUT, __LINE__, __func__, __FILE__, + "No MRI method at requested order, using q=3."); + q_actual = 3; + } + /* select method based on order and type */ /**** ImEx methods ****/ if (step_mem->implicit_rhs && step_mem->explicit_rhs) { - switch (step_mem->q) + switch (q_actual) { - case 3: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMEX_SD_3); - break; - case 4: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMEX_SD_4); - break; - default: - arkProcessError(ark_mem, ARK_ILL_INPUT, __LINE__, __func__, __FILE__, - "No MRI method at requested order, using q=3."); - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMEX_SD_3); - break; + case 1: table_id = MRISTEP_DEFAULT_IMEX_SD_1; break; + case 2: table_id = MRISTEP_DEFAULT_IMEX_SD_2; break; + case 3: table_id = MRISTEP_DEFAULT_IMEX_SD_2; break; + case 4: table_id = MRISTEP_DEFAULT_IMEX_SD_2; break; } /**** implicit methods ****/ } else if (step_mem->implicit_rhs) { - switch (step_mem->q) + switch (q_actual) { - case 2: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMPL_SD_3); - break; - case 3: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMPL_SD_3); - break; - case 4: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMPL_SD_4); - break; - default: - arkProcessError(ark_mem, ARK_ILL_INPUT, __LINE__, __func__, __FILE__, - "No MRI method at requested order, using q=3."); - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_IMPL_SD_3); - break; + case 1: table_id = MRISTEP_DEFAULT_IMPL_SD_1; break; + case 2: table_id = MRISTEP_DEFAULT_IMPL_SD_2; break; + case 3: table_id = MRISTEP_DEFAULT_IMPL_SD_3; break; + case 4: table_id = MRISTEP_DEFAULT_IMPL_SD_4; break; } /**** explicit methods ****/ } else { - switch (step_mem->q) + switch (q_actual) { - case 3: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_EXPL_3); - break; - case 4: - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_EXPL_4); - break; - default: - arkProcessError(ark_mem, ARK_ILL_INPUT, __LINE__, __func__, __FILE__, - "No MRI method at requested order, using q=3."); - step_mem->MRIC = MRIStepCoupling_LoadTable(MRISTEP_DEFAULT_EXPL_3); - break; + case 1: table_id = MRISTEP_DEFAULT_EXPL_1; break; + case 2: table_id = MRISTEP_DEFAULT_EXPL_2; break; + case 3: table_id = MRISTEP_DEFAULT_EXPL_3; break; + case 4: table_id = MRISTEP_DEFAULT_EXPL_4; break; } } + step_mem->MRIC = MRIStepCoupling_LoadTable(table_id); if (step_mem->MRIC == NULL) { arkProcessError(ark_mem, ARK_INVALID_TABLE, __LINE__, __func__, __FILE__, @@ -2882,7 +2903,7 @@ int mriStepInnerStepper_AllocVecs(MRIStepInnerStepper stepper, int count, /* Allocate fused operation workspace arrays */ if (stepper->vecs == NULL) { - stepper->vecs = (N_Vector*)calloc(count + 1, sizeof(N_Vector)); + stepper->vecs = (N_Vector*)calloc(count + 1, sizeof(*stepper->vecs)); if (stepper->vecs == NULL) { mriStepInnerStepper_FreeVecs(stepper); @@ -2892,7 +2913,7 @@ int mriStepInnerStepper_AllocVecs(MRIStepInnerStepper stepper, int count, if (stepper->vals == NULL) { - stepper->vals = (sunrealtype*)calloc(count + 1, sizeof(sunrealtype)); + stepper->vals = (sunrealtype*)calloc(count + 1, sizeof(*stepper->vals)); if (stepper->vals == NULL) { mriStepInnerStepper_FreeVecs(stepper); diff --git a/src/arkode/arkode_mristep_io.c b/src/arkode/arkode_mristep_io.c index e4ca3abaa2..d803141044 100644 --- a/src/arkode/arkode_mristep_io.c +++ b/src/arkode/arkode_mristep_io.c @@ -333,7 +333,7 @@ int mriStep_SetOrder(ARKodeMem ark_mem, int ord) if (retval) { return (retval); } /* check for illegal inputs */ - if (ord < 3 || ord > 4) { step_mem->q = 3; } + if (ord <= 0) { step_mem->q = 3; } else { step_mem->q = ord; } /* Clear tables, the user is requesting a change in method or a reset to diff --git a/src/arkode/fmod/farkode_mod.f90 b/src/arkode/fmod/farkode_mod.f90 index beac179931..f6252c31f4 100644 --- a/src/arkode/fmod/farkode_mod.f90 +++ b/src/arkode/fmod/farkode_mod.f90 @@ -315,7 +315,10 @@ module farkode_mod enumerator :: ARKODE_ESDIRK547L2SA_7_4_5 enumerator :: ARKODE_ESDIRK547L2SA2_7_4_5 enumerator :: ARKODE_ARK2_DIRK_3_1_2 - enumerator :: ARKODE_MAX_DIRK_NUM = ARKODE_ARK2_DIRK_3_1_2 + enumerator :: ARKODE_BACKWARD_EULER_1_1 + enumerator :: ARKODE_IMPLICIT_MIDPOINT_1_2 + enumerator :: ARKODE_IMPLICIT_TRAPEZOIDAL_2_2 + enumerator :: ARKODE_MAX_DIRK_NUM = ARKODE_IMPLICIT_TRAPEZOIDAL_2_2 end enum integer, parameter, public :: ARKODE_DIRKTableID = kind(ARKODE_DIRK_NONE) public :: ARKODE_DIRK_NONE, ARKODE_MIN_DIRK_NUM, ARKODE_SDIRK_2_1_2, ARKODE_BILLINGTON_3_3_2, ARKODE_TRBDF2_3_3_2, & @@ -324,7 +327,7 @@ module farkode_mod ARKODE_ARK437L2SA_DIRK_7_3_4, ARKODE_ARK548L2SAb_DIRK_8_4_5, ARKODE_ESDIRK324L2SA_4_2_3, ARKODE_ESDIRK325L2SA_5_2_3, & ARKODE_ESDIRK32I5L2SA_5_2_3, ARKODE_ESDIRK436L2SA_6_3_4, ARKODE_ESDIRK43I6L2SA_6_3_4, ARKODE_QESDIRK436L2SA_6_3_4, & ARKODE_ESDIRK437L2SA_7_3_4, ARKODE_ESDIRK547L2SA_7_4_5, ARKODE_ESDIRK547L2SA2_7_4_5, ARKODE_ARK2_DIRK_3_1_2, & - ARKODE_MAX_DIRK_NUM + ARKODE_BACKWARD_EULER_1_1, ARKODE_IMPLICIT_MIDPOINT_1_2, ARKODE_IMPLICIT_TRAPEZOIDAL_2_2, ARKODE_MAX_DIRK_NUM public :: FARKodeButcherTable_LoadDIRK public :: FARKodeButcherTable_LoadDIRKByName ! typedef enum ARKODE_ERKTableID @@ -353,7 +356,10 @@ module farkode_mod enumerator :: ARKODE_VERNER_10_6_7 enumerator :: ARKODE_VERNER_13_7_8 enumerator :: ARKODE_VERNER_16_8_9 - enumerator :: ARKODE_MAX_ERK_NUM = ARKODE_VERNER_16_8_9 + enumerator :: ARKODE_FORWARD_EULER_1_1 + enumerator :: ARKODE_RALSTON_EULER_2_1_2 + enumerator :: ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2 + enumerator :: ARKODE_MAX_ERK_NUM = ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2 end enum integer, parameter, public :: ARKODE_ERKTableID = kind(ARKODE_ERK_NONE) public :: ARKODE_ERK_NONE, ARKODE_MIN_ERK_NUM, ARKODE_HEUN_EULER_2_1_2, ARKODE_BOGACKI_SHAMPINE_4_2_3, & @@ -361,7 +367,8 @@ module farkode_mod ARKODE_CASH_KARP_6_4_5, ARKODE_FEHLBERG_6_4_5, ARKODE_DORMAND_PRINCE_7_4_5, ARKODE_ARK548L2SA_ERK_8_4_5, & ARKODE_VERNER_8_5_6, ARKODE_FEHLBERG_13_7_8, ARKODE_KNOTH_WOLKE_3_3, ARKODE_ARK437L2SA_ERK_7_3_4, & ARKODE_ARK548L2SAb_ERK_8_4_5, ARKODE_ARK2_ERK_3_1_2, ARKODE_SOFRONIOU_SPALETTA_5_3_4, ARKODE_SHU_OSHER_3_2_3, & - ARKODE_VERNER_9_5_6, ARKODE_VERNER_10_6_7, ARKODE_VERNER_13_7_8, ARKODE_VERNER_16_8_9, ARKODE_MAX_ERK_NUM + ARKODE_VERNER_9_5_6, ARKODE_VERNER_10_6_7, ARKODE_VERNER_13_7_8, ARKODE_VERNER_16_8_9, ARKODE_FORWARD_EULER_1_1, & + ARKODE_RALSTON_EULER_2_1_2, ARKODE_EXPLICIT_MIDPOINT_EULER_2_1_2, ARKODE_MAX_ERK_NUM public :: FARKodeButcherTable_LoadERK public :: FARKodeButcherTable_LoadERKByName ! typedef enum ARKODE_SPRKMethodID diff --git a/src/arkode/fmod/farkode_mristep_mod.f90 b/src/arkode/fmod/farkode_mristep_mod.f90 index 1b7b8aef81..782c3c1df9 100644 --- a/src/arkode/fmod/farkode_mristep_mod.f90 +++ b/src/arkode/fmod/farkode_mristep_mod.f90 @@ -47,18 +47,35 @@ module farkode_mristep_mod enumerator :: ARKODE_IMEX_MRI_GARK3a enumerator :: ARKODE_IMEX_MRI_GARK3b enumerator :: ARKODE_IMEX_MRI_GARK4 - enumerator :: ARKODE_MAX_MRI_NUM = ARKODE_IMEX_MRI_GARK4 + enumerator :: ARKODE_MRI_GARK_FORWARD_EULER + enumerator :: ARKODE_MRI_GARK_RALSTON2 + enumerator :: ARKODE_MRI_GARK_ERK22a + enumerator :: ARKODE_MRI_GARK_ERK22b + enumerator :: ARKODE_MRI_GARK_RALSTON3 + enumerator :: ARKODE_MRI_GARK_BACKWARD_EULER + enumerator :: ARKODE_MRI_GARK_IMPLICIT_MIDPOINT + enumerator :: ARKODE_IMEX_MRI_GARK_EULER + enumerator :: ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL + enumerator :: ARKODE_IMEX_MRI_GARK_MIDPOINT + enumerator :: ARKODE_MAX_MRI_NUM = ARKODE_IMEX_MRI_GARK_MIDPOINT end enum integer, parameter, public :: ARKODE_MRITableID = kind(ARKODE_MRI_NONE) public :: ARKODE_MRI_NONE, ARKODE_MIN_MRI_NUM, ARKODE_MIS_KW3, ARKODE_MRI_GARK_ERK33a, ARKODE_MRI_GARK_ERK45a, & ARKODE_MRI_GARK_IRK21a, ARKODE_MRI_GARK_ESDIRK34a, ARKODE_MRI_GARK_ESDIRK46a, ARKODE_IMEX_MRI_GARK3a, & - ARKODE_IMEX_MRI_GARK3b, ARKODE_IMEX_MRI_GARK4, ARKODE_MAX_MRI_NUM - integer(C_INT), parameter, public :: MRISTEP_DEFAULT_3 = ARKODE_MIS_KW3 + ARKODE_IMEX_MRI_GARK3b, ARKODE_IMEX_MRI_GARK4, ARKODE_MRI_GARK_FORWARD_EULER, ARKODE_MRI_GARK_RALSTON2, & + ARKODE_MRI_GARK_ERK22a, ARKODE_MRI_GARK_ERK22b, ARKODE_MRI_GARK_RALSTON3, ARKODE_MRI_GARK_BACKWARD_EULER, & + ARKODE_MRI_GARK_IMPLICIT_MIDPOINT, ARKODE_IMEX_MRI_GARK_EULER, ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL, & + ARKODE_IMEX_MRI_GARK_MIDPOINT, ARKODE_MAX_MRI_NUM + integer(C_INT), parameter, public :: MRISTEP_DEFAULT_EXPL_1 = ARKODE_MRI_GARK_FORWARD_EULER + integer(C_INT), parameter, public :: MRISTEP_DEFAULT_EXPL_2 = ARKODE_MRI_GARK_ERK22b integer(C_INT), parameter, public :: MRISTEP_DEFAULT_EXPL_3 = ARKODE_MIS_KW3 integer(C_INT), parameter, public :: MRISTEP_DEFAULT_EXPL_4 = ARKODE_MRI_GARK_ERK45a + integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMPL_SD_1 = ARKODE_MRI_GARK_BACKWARD_EULER integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMPL_SD_2 = ARKODE_MRI_GARK_IRK21a integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMPL_SD_3 = ARKODE_MRI_GARK_ESDIRK34a integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMPL_SD_4 = ARKODE_MRI_GARK_ESDIRK46a + integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMEX_SD_1 = ARKODE_IMEX_MRI_GARK_EULER + integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMEX_SD_2 = ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMEX_SD_3 = ARKODE_IMEX_MRI_GARK3b integer(C_INT), parameter, public :: MRISTEP_DEFAULT_IMEX_SD_4 = ARKODE_IMEX_MRI_GARK4 diff --git a/test/answers b/test/answers index ea6ac15fdc..073b119355 160000 --- a/test/answers +++ b/test/answers @@ -1 +1 @@ -Subproject commit ea6ac15fdcd8615e25e52692bcd453e2f003f46b +Subproject commit 073b119355058ac88a2e4abc87acfb007bc211f6 diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_0.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_0.out index 0b5e75274c..0632b1dd9d 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_0.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_0.out @@ -818,6 +818,111 @@ Fe RHS evals: expected: 68 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + ======================== Test implicit RK methods ======================== @@ -1797,6 +1902,123 @@ Fi RHS evals: expected: 17 -------------------- +======================== +DIRK Table ID 124 + stages: 1 + order: 1 + explicit 1st stage: 0 + stiffly accurate: 1 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + +======================== +DIRK Table ID 125 + stages: 1 + order: 2 + explicit 1st stage: 0 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- +Dense Output +Fi RHS evals: + actual: 10 + expected: 10 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 12 + expected: 12 +-------------------- + +======================== +DIRK Table ID 126 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 1 + first same as last: 1 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + ===================== Test IMEX ARK methods ===================== diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_1.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_1.out index 63acf12fd8..2594d1cd48 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_1.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_0_1.out @@ -818,6 +818,111 @@ Fe RHS evals: expected: 64 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + ======================== Test implicit RK methods ======================== @@ -1797,6 +1902,123 @@ Fi RHS evals: expected: 17 -------------------- +======================== +DIRK Table ID 124 + stages: 1 + order: 1 + explicit 1st stage: 0 + stiffly accurate: 1 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +DIRK Table ID 125 + stages: 1 + order: 2 + explicit 1st stage: 0 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +DIRK Table ID 126 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 1 + first same as last: 1 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + ===================== Test IMEX ARK methods ===================== diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_0.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_0.out index 02af9aea19..c4acdc8074 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_0.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_0.out @@ -818,6 +818,111 @@ Fe RHS evals: expected: 68 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + ======================== Test implicit RK methods ======================== @@ -1797,6 +1902,123 @@ Fi RHS evals: expected: 17 -------------------- +======================== +DIRK Table ID 124 + stages: 1 + order: 1 + explicit 1st stage: 0 + stiffly accurate: 1 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + +======================== +DIRK Table ID 125 + stages: 1 + order: 2 + explicit 1st stage: 0 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- +Dense Output +Fi RHS evals: + actual: 10 + expected: 10 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 12 + expected: 12 +-------------------- + +======================== +DIRK Table ID 126 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 1 + first same as last: 1 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + ===================== Test IMEX ARK methods ===================== diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_1.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_1.out index 8e54df6754..6d307f1e20 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_1.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_1_1.out @@ -818,6 +818,111 @@ Fe RHS evals: expected: 64 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + ======================== Test implicit RK methods ======================== @@ -1797,6 +1902,123 @@ Fi RHS evals: expected: 17 -------------------- +======================== +DIRK Table ID 124 + stages: 1 + order: 1 + explicit 1st stage: 0 + stiffly accurate: 1 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +DIRK Table ID 125 + stages: 1 + order: 2 + explicit 1st stage: 0 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +DIRK Table ID 126 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 1 + first same as last: 1 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + ===================== Test IMEX ARK methods ===================== diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_0.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_0.out index 363b66b2ac..bf6e07f6c5 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_0.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_0.out @@ -818,6 +818,111 @@ Fe RHS evals: expected: 68 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + ======================== Test implicit RK methods ======================== @@ -1797,6 +1902,123 @@ Fi RHS evals: expected: 17 -------------------- +======================== +DIRK Table ID 124 + stages: 1 + order: 1 + explicit 1st stage: 0 + stiffly accurate: 1 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + +======================== +DIRK Table ID 125 + stages: 1 + order: 2 + explicit 1st stage: 0 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- +Dense Output +Fi RHS evals: + actual: 10 + expected: 10 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 12 + expected: 12 +-------------------- + +======================== +DIRK Table ID 126 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 1 + first same as last: 1 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + ===================== Test IMEX ARK methods ===================== diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_1.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_1.out index bf8edcf74b..0269da8178 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_1.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_ark_2_1.out @@ -818,6 +818,111 @@ Fe RHS evals: expected: 64 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + ======================== Test implicit RK methods ======================== @@ -1797,6 +1902,123 @@ Fi RHS evals: expected: 17 -------------------- +======================== +DIRK Table ID 124 + stages: 1 + order: 1 + explicit 1st stage: 0 + stiffly accurate: 1 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +DIRK Table ID 125 + stages: 1 + order: 2 + explicit 1st stage: 0 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fi RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +DIRK Table ID 126 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 1 + first same as last: 1 +======================== +-------------------- +Steps: 1 +NLS iters: 1 +Fi RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 2 +NLS iters: 2 +Fi RHS evals: + actual: 5 + expected: 5 +-------------------- +Steps: 3 +NLS iters: 3 +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Dense Output +Fi RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +NLS iters: 4 +Fi RHS evals: + actual: 9 + expected: 9 +-------------------- + ===================== Test IMEX ARK methods ===================== diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_0.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_0.out index 7e2256037c..0a8aaf3158 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_0.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_0.out @@ -815,5 +815,110 @@ Fe RHS evals: expected: 68 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 7 + expected: 7 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + All tests passed! diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_1.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_1.out index 74bc1b6a8d..3851797978 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_1.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_erk_1.out @@ -815,5 +815,110 @@ Fe RHS evals: expected: 64 -------------------- +======================== +ERK Table ID 22 + stages: 1 + order: 1 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 1 + expected: 1 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Dense Output +Fe RHS evals: + actual: 3 + expected: 3 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- + +======================== +ERK Table ID 23 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + +======================== +ERK Table ID 24 + stages: 2 + order: 2 + explicit 1st stage: 1 + stiffly accurate: 0 + first same as last: 0 +======================== +-------------------- +Steps: 1 +Fe RHS evals: + actual: 2 + expected: 2 +-------------------- +Steps: 2 +Fe RHS evals: + actual: 4 + expected: 4 +-------------------- +Steps: 3 +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Dense Output +Fe RHS evals: + actual: 6 + expected: 6 +-------------------- +Steps: 4 +Fe RHS evals: + actual: 8 + expected: 8 +-------------------- + All tests passed! diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.cpp b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.cpp index 579e76ddef..7e3d20ca84 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.cpp +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.cpp @@ -24,6 +24,7 @@ #include #include #include +#include #include #include #include @@ -242,10 +243,8 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, // Evolve with various IMEX MRI methods // ------------------------------------ - // Methods to test (order most stages to least since reinit does not realloc) - int num_methods; - ARKODE_MRITableID* methods = nullptr; - bool* stiffly_accurate = nullptr; + // Methods to test paired with whether they are stiffly accurate + std::map methods; if (type == MRISTEP_EXPLICIT) { @@ -253,12 +252,14 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, cout << "Test explicit MRI methods\n"; cout << "=========================\n"; - num_methods = 3; - methods = new ARKODE_MRITableID[num_methods]; - - methods[0] = ARKODE_MIS_KW3; - methods[1] = ARKODE_MRI_GARK_ERK33a; - methods[2] = ARKODE_MRI_GARK_ERK45a; + methods.insert({{"ARKODE_MRI_GARK_FORWARD_EULER", false}, + {"ARKODE_MRI_GARK_ERK22a", false}, + {"ARKODE_MRI_GARK_ERK22b", false}, + {"ARKODE_MRI_GARK_RALSTON2", false}, + {"ARKODE_MIS_KW3", false}, + {"ARKODE_MRI_GARK_ERK33a", false}, + {"ARKODE_MRI_GARK_RALSTON3", false}, + {"ARKODE_MRI_GARK_ERK45a", false}}); } else if (type == MRISTEP_IMPLICIT) { @@ -266,18 +267,11 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, cout << "Test implicit MRI methods\n"; cout << "=========================\n"; - num_methods = 3; - methods = new ARKODE_MRITableID[num_methods]; - stiffly_accurate = new bool[num_methods]; - - methods[0] = ARKODE_MRI_GARK_IRK21a; - stiffly_accurate[0] = true; - - methods[1] = ARKODE_MRI_GARK_ESDIRK34a; - stiffly_accurate[1] = true; - - methods[2] = ARKODE_MRI_GARK_ESDIRK46a; - stiffly_accurate[2] = true; + methods.insert({{"ARKODE_MRI_GARK_BACKWARD_EULER", true}, + {"ARKODE_MRI_GARK_IRK21a", true}, + {"ARKODE_MRI_GARK_IMPLICIT_MIDPOINT", false}, + {"ARKODE_MRI_GARK_ESDIRK34a", true}, + {"ARKODE_MRI_GARK_ESDIRK46a", true}}); } else if (type == MRISTEP_IMEX) { @@ -285,31 +279,27 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, cout << "Test IMEX MRI methods\n"; cout << "=====================\n"; - num_methods = 3; - methods = new ARKODE_MRITableID[num_methods]; - stiffly_accurate = new bool[num_methods]; - - methods[0] = ARKODE_IMEX_MRI_GARK3a; - stiffly_accurate[0] = false; - - methods[1] = ARKODE_IMEX_MRI_GARK3b; - stiffly_accurate[1] = false; - - methods[2] = ARKODE_IMEX_MRI_GARK4; - stiffly_accurate[2] = false; + methods.insert({{"ARKODE_IMEX_MRI_GARK_EULER", true}, + {"ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL", false}, + {"ARKODE_IMEX_MRI_GARK_MIDPOINT", false}, + {"ARKODE_IMEX_MRI_GARK3a", false}, + {"ARKODE_IMEX_MRI_GARK3b", false}, + {"ARKODE_IMEX_MRI_GARK4", false}}); } else { return 1; } - for (int i = 0; i < num_methods; i++) + for (const auto& pair : methods) { - cout << "\nTesting method " << i << "\n"; + std::string id = pair.first; + bool stiffly_accurate = pair.second; + cout << "\nTesting method " << id << "\n"; // ------------- // Select method // ------------- // Load method table - MRIStepCoupling C = MRIStepCoupling_LoadTable(methods[i]); + MRIStepCoupling C = MRIStepCoupling_LoadTableByName(id.c_str()); if (check_flag((void*)C, "MRIStepCoupling_LoadTable", 0)) { return 1; } MRIStepCoupling_Write(C, stdout); @@ -415,10 +405,12 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, cout << "\nComparing Solver Statistics:\n"; + int nstages_evaluated = nstages_stored; + if (stiffly_accurate) nstages_evaluated--; long int fe_evals = 0; if (type == MRISTEP_EXPLICIT || type == MRISTEP_IMEX) { - fe_evals = mri_nst * nstages_stored; + fe_evals = mri_nst * nstages_evaluated; } if (mri_nfse != fe_evals) @@ -430,17 +422,7 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, long int fi_evals = 0; if (type == MRISTEP_IMPLICIT || type == MRISTEP_IMEX) { - if (stiffly_accurate[i]) - { - // The last stage is implicit so it does not correspond to a column of - // zeros in the coupling matrix and is counted in "nstages_stored" - // however we do not evaluate the RHS functions after the solve since - // the methods is "FSAL" (the index map value and allocated space is - // used in the nonlinear for this stage). The RHS functions will be - // evaluated and stored at the start of the next step. - fi_evals = mri_nst * (nstages_stored - 1) + mri_nni; - } - else { fi_evals = mri_nst * nstages_stored + mri_nni; } + fi_evals = mri_nst * nstages_evaluated + mri_nni; } if (mri_nfsi != fi_evals) @@ -493,8 +475,6 @@ int run_tests(MRISTEP_METHOD_TYPE type, sunrealtype t0, int nsteps, SUNMatDestroy(A); } N_VDestroy(y); - delete[] methods; - delete[] stiffly_accurate; return numfails; } diff --git a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.out b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.out index 5e7243e3e8..e5f2a9b602 100644 --- a/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.out +++ b/test/unit_tests/arkode/CXX_serial/ark_test_dahlquist_mri.out @@ -12,7 +12,7 @@ Dahlquist ODE test problem: Test explicit MRI methods ========================= -Testing method 0 +Testing method ARKODE_MIS_KW3 nmat = 1 stages = 4 method order (q) = 3 @@ -38,7 +38,57 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed -Testing method 1 +Testing method ARKODE_MRI_GARK_ERK22a + nmat = 1 + stages = 3 + method order (q) = 2 + embedding order (p) = 0 + c = 0 0.5 1 + W[0] = + 0 0 0 + 0.5 0 0 + -0.5 1 0 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.980199 + Error = -4.95647e-07 + Steps = 1 + Fe evals = 2 + Fi evals = 0 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_MRI_GARK_ERK22b + nmat = 1 + stages = 3 + method order (q) = 2 + embedding order (p) = 0 + c = 0 1 1 + W[0] = + 0 0 0 + 1 0 0 + -0.5 0.5 0 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.980199 + Error = -4.95856e-07 + Steps = 1 + Fe evals = 2 + Fi evals = 0 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_MRI_GARK_ERK33a nmat = 2 stages = 4 method order (q) = 3 @@ -70,7 +120,7 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed -Testing method 2 +Testing method ARKODE_MRI_GARK_ERK45a nmat = 2 stages = 6 method order (q) = 4 @@ -106,28 +156,109 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed +Testing method ARKODE_MRI_GARK_FORWARD_EULER + nmat = 1 + stages = 2 + method order (q) = 1 + embedding order (p) = 0 + c = 0 1 + W[0] = + 0 0 + 1 0 + + Stored stages = 1 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.9801 + Error = 9.90058e-05 + Steps = 1 + Fe evals = 1 + Fi evals = 0 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_MRI_GARK_RALSTON2 + nmat = 1 + stages = 3 + method order (q) = 2 + embedding order (p) = 0 + c = 0 0.6666666666666666 1 + W[0] = + 0 0 0 + 0.6666666666666666 0 0 + -0.4166666666666666 0.75 0 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.980199 + Error = -4.9567e-07 + Steps = 1 + Fe evals = 2 + Fi evals = 0 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_MRI_GARK_RALSTON3 + nmat = 2 + stages = 4 + method order (q) = 3 + embedding order (p) = 0 + c = 0 0.5 0.75 1 + W[0] = + 0 0 0 0 + 0.5 0 0 0 + -2.75 3 0 0 + 1.305555555555556 -0.1666666666666667 -0.8888888888888888 0 + + W[1] = + 0 0 0 0 + 0 0 0 0 + 4.5 -4.5 0 0 + -2.166666666666667 -0.5 2.666666666666667 0 + + Stored stages = 3 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.980199 + Error = 1.72143e-09 + Steps = 1 + Fe evals = 3 + Fi evals = 0 + +Comparing Solver Statistics: +All checks passed + ========================= Test implicit MRI methods ========================= -Testing method 0 +Testing method ARKODE_MRI_GARK_BACKWARD_EULER nmat = 1 stages = 3 - method order (q) = 2 + method order (q) = 1 embedding order (p) = 0 c = 0 1 1 G[0] = 0 0 0 1 0 0 - -0.5 0 0.5 + -1 0 1 Stored stages = 2 MRIStep Statistics: Time = 0.01 y(t) = 0.980199 - y_n = 0.980199 - Error = -8.22598e-10 + y_n = 0.980297 + Error = -9.80272e-05 Steps = 1 Fe evals = 0 Fi evals = 2 @@ -140,7 +271,7 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed -Testing method 1 +Testing method ARKODE_MRI_GARK_ESDIRK34a nmat = 1 stages = 7 method order (q) = 3 @@ -174,7 +305,7 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed -Testing method 2 +Testing method ARKODE_MRI_GARK_ESDIRK46a nmat = 2 stages = 11 method order (q) = 4 @@ -225,11 +356,72 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed +Testing method ARKODE_MRI_GARK_IMPLICIT_MIDPOINT + nmat = 1 + stages = 4 + method order (q) = 2 + embedding order (p) = 0 + c = 0 0.5 0.5 1 + G[0] = + 0 0 0 0 + 0.5 0 0 0 + -0.5 0 0.5 0 + 0 0 0.5 0 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.980199 + Error = 1.2304e-07 + Steps = 1 + Fe evals = 0 + Fi evals = 3 + NLS iters = 1 + NLS fails = 0 + LS setups = 1 + LS Fi evals = 0 + Ji evals = 1 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_MRI_GARK_IRK21a + nmat = 1 + stages = 3 + method order (q) = 2 + embedding order (p) = 0 + c = 0 1 1 + G[0] = + 0 0 0 + 1 0 0 + -0.5 0 0.5 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.980199 + y_n = 0.980199 + Error = -8.22598e-10 + Steps = 1 + Fe evals = 0 + Fi evals = 2 + NLS iters = 1 + NLS fails = 0 + LS setups = 1 + LS Fi evals = 0 + Ji evals = 1 + +Comparing Solver Statistics: +All checks passed + ===================== Test IMEX MRI methods ===================== -Testing method 0 +Testing method ARKODE_IMEX_MRI_GARK3a nmat = 1 stages = 8 method order (q) = 3 @@ -274,7 +466,7 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed -Testing method 1 +Testing method ARKODE_IMEX_MRI_GARK3b nmat = 1 stages = 8 method order (q) = 3 @@ -319,7 +511,7 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed -Testing method 2 +Testing method ARKODE_IMEX_MRI_GARK4 nmat = 2 stages = 12 method order (q) = 4 @@ -400,5 +592,114 @@ MRIStep Statistics: Comparing Solver Statistics: All checks passed +Testing method ARKODE_IMEX_MRI_GARK_EULER + nmat = 1 + stages = 3 + method order (q) = 1 + embedding order (p) = 0 + c = 0 1 1 + W[0] = + 0 0 0 + 1 0 0 + 0 0 0 + + G[0] = + 0 0 0 + 1 0 0 + -1 0 1 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.970446 + y_n = 0.970445 + Error = 4.82806e-07 + Steps = 1 + Fe evals = 1 + Fi evals = 2 + NLS iters = 1 + NLS fails = 0 + LS setups = 1 + LS Fi evals = 0 + Ji evals = 1 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_IMEX_MRI_GARK_MIDPOINT + nmat = 1 + stages = 4 + method order (q) = 2 + embedding order (p) = 0 + c = 0 0.5 0.5 1 + W[0] = + 0 0 0 0 + 0.5 0 0 0 + 0 0 0 0 + -0.5 0 1 0 + + G[0] = + 0 0 0 0 + 0.5 0 0 0 + -0.5 0 0.5 0 + 0 0 0.5 0 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.970446 + y_n = 0.970446 + Error = -8.01486e-07 + Steps = 1 + Fe evals = 2 + Fi evals = 3 + NLS iters = 1 + NLS fails = 0 + LS setups = 1 + LS Fi evals = 0 + Ji evals = 1 + +Comparing Solver Statistics: +All checks passed + +Testing method ARKODE_IMEX_MRI_GARK_TRAPEZOIDAL + nmat = 1 + stages = 4 + method order (q) = 2 + embedding order (p) = 0 + c = 0 1 1 1 + W[0] = + 0 0 0 0 + 1 0 0 0 + 0 0 0 0 + -0.5 0 0.5 0 + + G[0] = + 0 0 0 0 + 1 0 0 0 + -0.5 0 0.5 0 + 0 0 0 0 + + Stored stages = 2 + +MRIStep Statistics: + Time = 0.01 + y(t) = 0.970446 + y_n = 0.970447 + Error = -9.8759e-07 + Steps = 1 + Fe evals = 2 + Fi evals = 3 + NLS iters = 1 + NLS fails = 0 + LS setups = 1 + LS Fi evals = 0 + Ji evals = 1 + +Comparing Solver Statistics: +All checks passed + All tests passed!