Discontinuous Galerkin Library
#include "dg/algorithm.h"
dg::ShuOsherTableau< real_type > Struct Template Reference

Manage coefficients in Shu-Osher form. More...

Public Types

using value_type = real_type
 

Public Member Functions

 ShuOsherTableau ()=default
 No memory allocation. More...
 
 ShuOsherTableau (unsigned stages, unsigned order, const std::vector< real_type > &alpha_v, const std::vector< real_type > &beta_v)
 Construct a non-embedded explicit tableau. More...
 
 operator ButcherTableau< real_type > () const
 A Shu-Osher Tableau can be converted to a Butcher table. More...
 
real_type alpha (unsigned i, unsigned j)
 Read the alpha_ij coefficients. More...
 
real_type beta (unsigned i, unsigned j)
 Read the beta_ij coefficients. More...
 
unsigned num_stages () const
 The number of stages s. More...
 
unsigned order () const
 global order of accuracy for the method More...
 

Detailed Description

template<class real_type>
struct dg::ShuOsherTableau< real_type >

Manage coefficients in Shu-Osher form.

Currently only explicit tables that are marked with "Shu-Osher-Form" are available in this form

We follow the naming convention of the ARKode library https://sundials.readthedocs.io/en/latest/arkode/Butcher_link.html (They also provide nice stability plots for their methods) as NAME-S-P-Q or NAME-S-Q, where

  • NAME is the author or name of the method
  • S is the number of stages in the method
  • P is the global order of the embedding
  • Q is the global order of the method
Name Identifier Description
Euler dg::EXPLICIT_EULER_1_1 Explicit Euler
Midpoint-2-2 dg::MIDPOINT_2_2 Midpoint method
Kutta-3-3 dg::KUTTA_3_3 Kutta-3-3
Runge-Kutta-4-4 dg::CLASSIC_4_4 "The" Runge-Kutta method
SSPRK-2-2 dg::SSPRK_2_2 SSPRK (2,2) CFL_eff = 0.5 "Shu-Osher-Form"
SSPRK-3-2 dg::SSPRK_3_2 SSPRK (3,2) CFL_eff = 0.66 "Shu-Osher-Form"
SSPRK-3-3 dg::SSPRK_3_3 SSPRK (3,3) CFL_eff = 0.33 "Shu-Osher-Form"
SSPRK-5-3 dg::SSPRK_5_3 SSPRK (5,3) CFL_eff = 0.5 "Shu-Osher-Form"
SSPRK-5-4 dg::SSPRK_5_4 SSPRK (5,4) CFL_eff = 0.37 "Shu-Osher-Form"
Heun-Euler-2-1-2 dg::HEUN_EULER_2_1_2 Heun-Euler-2-1-2
Cavaglieri-3-1-2 (explicit) dg::CAVAGLIERI_3_1_2 Low-storage implicit/explicit Runge-Kutta schemes for the simulation of stiff high-dimensional ODE systems IMEXRKCB2 scheme
Fehlberg-3-2-3 dg::FEHLBERG_3_2_3 The original uses the embedding as the solution [Hairer, Noersett, Wanner, Solving ordinary differential Equations I, 1987]
Fehlberg-4-2-3 dg::FEHLBERG_4_2_3 The original uses the embedding as the solution [Hairer, Noersett, Wanner, Solving ordinary differential Equations I, 1987] (fsal)
Bogacki-Shampine-4-2-3 dg::BOGACKI_SHAMPINE_4_2_3 Bogacki-Shampine (fsal)
Cavaglieri-4-2-3 (explicit) dg::CAVAGLIERI_4_2_3 Low-storage implicit/explicit Runge-Kutta schemes for the simulation of stiff high-dimensional ODE systems The SSP scheme IMEXRKCB3c (explicit part)
ARK-4-2-3 (explicit) dg::ARK324L2SA_ERK_4_2_3 ARK-4-2-3 (explicit)
Zonneveld-5-3-4 dg::ZONNEVELD_5_3_4 Zonneveld-5-3-4
ARK-6-3-4 (explicit) dg::ARK436L2SA_ERK_6_3_4 ARK-6-3-4 (explicit)
Sayfy_Aburub-6-3-4 dg::SAYFY_ABURUB_6_3_4 Sayfy_Aburub_6_3_4
Cash_Karp-6-4-5 dg::CASH_KARP_6_4_5 Cash-Karp
Fehlberg-6-4-5 dg::FEHLBERG_6_4_5 Runge-Kutta-Fehlberg
Dormand-Prince-7-4-5 dg::DORMAND_PRINCE_7_4_5 Dormand-Prince method (fsal)
Tsitouras09-7-4-5 dg::TSITOURAS09_7_4_5 Tsitouras 5(4) method from 2009 (fsal), The default method in Julia
Tsitouras11-7-4-5 dg::TSITOURAS11_7_4_5 Tsitouras 5(4) method from 2011 (fsal) Further improves Tsitouras09 (Note that in the paper b-bt is printed instead of bt)
ARK-8-4-5 (explicit) dg::ARK548L2SA_ERK_8_4_5 Kennedy and Carpenter (2019) Optimum ARK_2 method (explicit part)
Verner-9-5-6 dg::VERNER_9_5_6 Verner-9-5-6 (fsal)
Verner-10-6-7 dg::VERNER_10_6_7 Verner-10-6-7
Fehlberg-13-7-8 dg::FEHLBERG_13_7_8 Fehlberg-13-7-8
Dormand-Prince-13-7-8 dg::DORMAND_PRINCE_13_7_8 [Hairer, Noersett, Wanner, Solving ordinary differential Equations I, 1987]
Feagin-17-8-10 dg::FEAGIN_17_8_10 Feagin (The RK10(8) method)
Note
A Shu-Osher Tableau can be uniquely converted to a ButcherTableau but the converse is not true
Template Parameters
real_typetype of the coefficients
See also
ShuOsher

Member Typedef Documentation

◆ value_type

template<class real_type >
using dg::ShuOsherTableau< real_type >::value_type = real_type

Constructor & Destructor Documentation

◆ ShuOsherTableau() [1/2]

template<class real_type >
dg::ShuOsherTableau< real_type >::ShuOsherTableau ( )
default

No memory allocation.

◆ ShuOsherTableau() [2/2]

template<class real_type >
dg::ShuOsherTableau< real_type >::ShuOsherTableau ( unsigned  stages,
unsigned  order,
const std::vector< real_type > &  alpha_v,
const std::vector< real_type > &  beta_v 
)
inline

Construct a non-embedded explicit tableau.

Parameters
stagesnumber of stages
order(global) order of the resulting method
alpha_vs*s/2 real numbers interpreted as a linearized triangle
beta_vs*s/2 real numbers interpreted as a linearized triangle

Member Function Documentation

◆ alpha()

template<class real_type >
real_type dg::ShuOsherTableau< real_type >::alpha ( unsigned  i,
unsigned  j 
)
inline

Read the alpha_ij coefficients.

Parameters
irow number 0<=i<s, i>=s results in undefined behaviour
jcol number 0<=j<s, j>=s results in undefined behaviour
Returns
alpha_ij

◆ beta()

template<class real_type >
real_type dg::ShuOsherTableau< real_type >::beta ( unsigned  i,
unsigned  j 
)
inline

Read the beta_ij coefficients.

Parameters
irow number 0<=i<s, i>=s results in undefined behaviour
jcol number 0<=j<s, j>=s results in undefined behaviour
Returns
beta_ij

◆ num_stages()

template<class real_type >
unsigned dg::ShuOsherTableau< real_type >::num_stages ( ) const
inline

The number of stages s.

◆ operator ButcherTableau< real_type >()

template<class real_type >
dg::ShuOsherTableau< real_type >::operator ButcherTableau< real_type > ( ) const
inline

A Shu-Osher Tableau can be converted to a Butcher table.

Returns
the corresponding Butcher Tableau

◆ order()

template<class real_type >
unsigned dg::ShuOsherTableau< real_type >::order ( ) const
inline

global order of accuracy for the method


The documentation for this struct was generated from the following file: