File imex-bdf2.hxx

class IMEXBDF2 : public Solver
#include <imex-bdf2.hxx>

IMEX-BDF2 time integration solver

Scheme taken from this paper: http://homepages.cwi.nl/~willem/DOCART/JCP07.pdf

W.Hundsdorfer, S.J.Ruuth “IMEX extensions of linear multistep methods with general

monotonicity and boundedness properties” JCP 225 (2007) 2016-2042

The method has been extended to variable order, variable timestep, and includes some adaptive capabilities

Public Functions

IMEXBDF2(Options *opt = nullptr)
~IMEXBDF2()
inline virtual BoutReal getCurrentTimestep() override

Returns the current internal timestep.

virtual int init(int nout, BoutReal tstep) override

Initialise solver. Must be called once and only once

Initialisation routine. Called once before solve.

Parameters
  • nout[in] Number of outputs

  • tstep[in] Time between outputs. NB: Not internal timestep

virtual int run() override

Run the simulation.

PetscErrorCode snes_function(Vec x, Vec f, bool linear)

Nonlinear function. This is called by PETSc SNES object via a static C-style function. For implicit time integration this function calculates:

f = (x - gamma*G(x)) - rhs

Parameters
  • x[in] The state vector

  • f[out] The vector for the result f(x)

  • linear[in] Specifies that the SNES solver is in a linear (KSP) inner loop, so the operator should be linearised if possible

PetscErrorCode precon(Vec x, Vec f)

Preconditioner. Called by PCapply via a C-style static function.

Parameters
  • x[in] The vector to be operated on

  • f[out] The result of the operation

Private Functions

void take_step(BoutReal curtime, BoutReal dt, int order = 2)

Take a full step at requested order

Take a full IMEX-BDF step of order “order”. Note that this assumes that enough time points are already available (in u and f).

Inputs: u* - Solution history f* - Non-stiff component history

Outputs: u - Latest Solution f1 - Non-stiff time derivative at current time

Parameters
  • curtime[in] The current simulation time

  • dt[in] The time step to take

  • order[in] The order of accuracy

void constructSNES(SNES *snesIn)

Setup a SNES object This includes creating, setting functions, options, and internal (KSP) solver, and Jacobian options including coloring.

void shuffleState()

Shuffle state along one step.

void calculateCoeffs(int order)

Populate the *Fac vectors and dtImp with appropriate coefficients for this order.

PetscErrorCode solve_implicit(BoutReal curtime, BoutReal gamma)
template<class Op>
void loopVars(BoutReal *u)

Loop over arrays, using template parameter to specify the operation to be performed at each point

void saveVars(BoutReal *u)

Save variables from BOUT++ fields into a pre-allocated array u

Copy data from fields into array

void loadVars(BoutReal *u)

Load variables from input vector u into BOUT++ fields.

Copy data from array into fields

void saveDerivs(BoutReal *u)

Save time derivatives from ddt() fields into a preallocated array u.

Copy time derivatives from fields into array

Private Members

int maxOrder

Specify the maximum order of the scheme to use (1/2/3)

BoutReal out_timestep

The output timestep.

int nsteps

Number of output steps.

BoutReal timestep

The internal timestep.

int ninternal

Number of internal steps per output.

int mxstep

Maximum number of internal steps between outputs.

bool adaptive

Use adaptive timestepping?

int nadapt

How often do we check the error.

int mxstepAdapt

Maximum no. consecutive times we try to reduce timestep.

BoutReal scaleCushUp

Don’t increase timestep if scale factor < 1.0+scaleCushUp.

BoutReal scaleCushDown

Don’t decrease timestep if scale factor > 1.0-scaleCushDown.

BoutReal adaptRtol

Target relative error for adaptivity.

BoutReal dtMin

Minimum timestep we want to use.

BoutReal dtMax

Maximum timestep we want to use.

BoutReal dtMinFatal

If timestep wants to drop below this we abort. Set -ve to deactivate

std::vector<BoutReal> uFac
std::vector<BoutReal> fFac
std::vector<BoutReal> gFac
BoutReal dtImp
int nlocal
int neq

Number of variables on local processor and in total.

Array<BoutReal> u

System state at current time.

std::vector<Array<BoutReal>> uV

The solution history.

std::vector<Array<BoutReal>> fV

The non-stiff solution history.

std::vector<BoutReal> timesteps

Timestep history.

Array<BoutReal> rhs
Array<BoutReal> err
BoutReal implicit_gamma
BoutReal implicit_curtime
int predictor

Predictor method.

PetscLib lib

Handles initialising, finalising PETSc.

Vec snes_f

Used by SNES to store function.

Vec snes_x

Result of SNES.

SNES snes

SNES context.

SNES snesAlt

Alternative SNES object for adaptive checks.

SNES snesUse

The snes object to use in solve stage. Allows easy switching.

Mat Jmf

Matrix-free Jacobian.

bool diagnose

Output diagnostics every timestep.

bool verbose

Gives a more verbose output for each timestep.

int linear_fails

Number of linear (KSP) convergence failures.

int nonlinear_fails

Numbef of nonlinear (SNES) convergence failures.

bool have_constraints

Are there any constraint variables?

Array<BoutReal> is_dae

If using constraints, 1 -> DAE, 0 -> AE.

MatFDColoring fdcoloring

Matrix coloring context, used for finite difference Jacobian evaluation

Private Static Attributes

static constexpr int MAX_SUPPORTED_ORDER = 4