Time Steppers

Pydro offers a few different time steppers, including many strong-stability-preserving (SSP) time steppers. The linear Runge-Kutta (RK) time steppers, while simple, should be avoided since they can only achieve the promised order of accuracy for linear systems. The SSP RK3 and SSP RK4 are both robust choices for nonlinear hydrodynamics simulations. The SSP RK4 has a maximum step size roughly 50% larger than an Euler step and the SSP RK3, making it a more efficient method than the SSP RK3 There are also (mostly complete) Adams-Bashforth time steppers up to fourth order, which still require a self-starting procedure in order to truly reach the high-order accuracy.

class TimeStepper.LinearRk4Ssp(time_deriv, initial_state, initial_time)

A fourth-order SSP linear Runge-Kutta method from [6]

Warning

This time stepper is a linear method and should not be used on nonlinear systems.

Warning

This time stepper passes None as the time.

Parameters
  • time_deriv (funcion) – The time derivative function which must be invokable with two arguments, the evolved variables and the time.

  • initial_state (list) – A list of the values of the variables at the initial time.

  • time (double) – Initial time.

get_cfl_coefficient()

Returns the CFL factor required for stability.

get_evolved_vars()

Returns an list of the evolved variables.

get_time()

Returns the current time.

take_step(dt)

Take a time step.

Parameters

dt (double) – The time step size to use already including the CFL coefficient for the time stepper.

class TimeStepper.LinearRk6Ssp(time_deriv, initial_state, initial_time)

A sixth-order SSP linear Runge-Kutta method from [6]

Warning

This time stepper is a linear method and should not be used on nonlinear systems.

Warning

This time stepper passes None as the time.

Parameters
  • time_deriv (funcion) – The time derivative function which must be invokable with two arguments, the evolved variables and the time.

  • initial_state (list) – A list of the values of the variables at the initial time.

  • time (double) – Initial time.

get_cfl_coefficient()

Returns the CFL factor required for stability.

get_evolved_vars()

Returns an list of the evolved variables.

get_time()

Returns the current time.

take_step(dt)

Take a time step.

Parameters

dt (double) – The time step size to use already including the CFL coefficient for the time stepper.

class TimeStepper.LinearRk8Ssp(time_deriv, initial_state, initial_time)

A eighth-order SSP linear Runge-Kutta method from [6]

Warning

This time stepper is a linear method and should not be used on nonlinear systems.

Warning

This time stepper passes None as the time.

Parameters
  • time_deriv (funcion) – The time derivative function which must be invokable with two arguments, the evolved variables and the time.

  • initial_state (list) – A list of the values of the variables at the initial time.

  • time (double) – Initial time.

get_cfl_coefficient()

Returns the CFL factor required for stability.

get_evolved_vars()

Returns an list of the evolved variables.

get_time()

Returns the current time.

take_step(dt)

Take a time step.

Parameters

dt (double) – The time step size to use already including the CFL coefficient for the time stepper.

class TimeStepper.Rk3Ssp(time_deriv, initial_state, initial_time)

A third-order strong-stability-preserving nonlinear Runge-Kutta time stepper [7]

Denoting the time derivative operator by \(\mathcal{L}\), the stepper is given by

\[\begin{split}\begin{align} v^{(1)} &=q^n+\Delta t \mathcal{L}(u^n, t^n) \\ v^{(2)} &= \frac{1}{4}\left[3q^n + v^{(1)} + \Delta t \mathcal{L}\left(v^{(1)}, t^n + \Delta t\right)\right] \\ q^{n+1} &=\frac{1}{3}\left[q^n + 2 v^{(2)} + 2 \Delta t \mathcal{L}\left(v^{(2)},t^n + \frac{1}{2}\Delta t \right)\right] \end{align}\end{split}\]
Parameters
  • time_deriv (funcion) – The time derivative function which must be invokable with two arguments, the evolved variables and the time.

  • initial_state (list) – A list of the values of the variables at the initial time.

  • time (double) – Initial time.

get_cfl_coefficient()

Returns the CFL factor required for stability.

get_evolved_vars()

Returns an list of the evolved variables.

get_time()

Returns the current time.

take_step(dt)

Take a time step.

Parameters

dt (double) – The time step size to use already including the CFL coefficient for the time stepper.

class TimeStepper.Rk4Ssp(time_deriv, initial_state, initial_time)

A fourth-order strong-stability-preserving nonlinear Runge-Kutta time stepper [7]

Denoting the time derivative operator by \(\mathcal{L}\), the stepper is given by

\[\begin{split}\begin{align} v^{(1)} &=q^n + 0.39175222700392 \Delta t \mathcal{L}(u^n,t^n) \\ v^{(2)} &= 0.44437049406734 q^n + 0.55562950593266 v^{(1)} \\ &+ 0.36841059262959 \Delta t \mathcal{L}\left( v^{(1)}, t^n + 0.39175222700392 \Delta t \right) \\ v^{(3)} &= 0.6201018513854 q^n + 0.3798981486146 v^{(2)} \\ &+ 0.25189177424738 \Delta t \mathcal{L}\left( v^{(2)}, t^n + 0.5860796889678 \Delta t \right) \\ v^{(4)} &= 0.17807995410773 q^n + 0.82192004589227 v^{(3)} \\ &+ 0.54497475021237 \Delta t \mathcal{L}\left( v^{(3)}, t^n + 0.47454236302687 \Delta t \right) \\ q^{n+1} &= 0.00683325884039 q^n + 0.51723167208978 v^{(2)} \\ &+ 0.12759831133288 v^{(3)} + 0.34833675773694 v^{(4)} \\ &+ 0.08460416338212 \Delta t \mathcal{L}\left( v^{(3)}, t^n + 0.47454236302687 \Delta t \right) \\ &+ 0.22600748319395 \Delta t \mathcal{L}\left( v^{(4)}, t^n + 0.93501063100924 \Delta t \right) \\ \end{align}\end{split}\]
Parameters
  • time_deriv (funcion) – The time derivative function which must be invokable with two arguments, the evolved variables and the time.

  • initial_state (list) – A list of the values of the variables at the initial time.

  • time (double) – Initial time.

get_cfl_coefficient()

Returns the CFL factor required for stability.

get_evolved_vars()

Returns an list of the evolved variables.

get_time()

Returns the current time.

take_step(dt)

Take a time step.

Parameters

dt (double) – The time step size to use already including the CFL coefficient for the time stepper.