File solver.hxx#

Defines

BOUT_NO_USING_NAMESPACE_BOUTGLOBALS

Typedefs

using rhsfunc = int (*)(BoutReal): RHS function pointer.

using PhysicsPrecon = int (*)(BoutReal t, BoutReal gamma, BoutReal delta)#: User-supplied preconditioner function.

using Jacobian = int (*)(BoutReal t)#: User-supplied Jacobian function.

using TimestepMonitorFunc = int (*)(Solver *solver, BoutReal simtime, BoutReal lastdt)#: Solution monitor, called each timestep.

using SolverType = std::string#

using RegisterSolver = SolverFactory::RegisterInFactory<DerivedType>#

Simpler name for Factory registration helper class

Usage:

#include <bout/solverfactory.hxx>
namespace {
RegisterSolver<MySolver> registersolvermine("mysolver");
}

using RegisterUnavailableSolver = SolverFactory::RegisterUnavailableInFactory#

Enums

enum SOLVER_VAR_OP#

Values:

enumerator LOAD_VARS#

enumerator LOAD_DERIVS#

enumerator SET_ID#

enumerator SAVE_VARS#

enumerator SAVE_DERIVS#

enum MonitorPosition#

A type to set where in the list monitors are added.

Values:

enumerator BACK#

enumerator FRONT#

Variables

constexpr auto SOLVERCVODE = "cvode"#

constexpr auto SOLVERPVODE = "pvode"#

constexpr auto SOLVERIDA = "ida"#

constexpr auto SOLVERPETSC = "petsc"#

constexpr auto SOLVERSLEPC = "slepc"#

constexpr auto SOLVERRK4 = "rk4"#

constexpr auto SOLVEREULER = "euler"#

constexpr auto SOLVERRK3SSP = "rk3ssp"#

constexpr auto SOLVERPOWER = "power"#

constexpr auto SOLVERARKODE = "arkode"#

constexpr auto SOLVERIMEXBDF2 = "imexbdf2"#

constexpr auto SOLVERSNES = "snes"#

constexpr auto SOLVERRKGENERIC = "rkgeneric"#

class SolverFactory : public Factory<Solver, SolverFactory, Options*>#

Public Static Attributes

static constexpr auto type_name = "Solver"#

static constexpr auto section_name = "solver"#

static constexpr auto option_name = "type"#

static constexpr auto default_type = SOLVERCVODE #

class Solver#

#include <solver.hxx>

Interface to integrators, mainly for time integration

Creation

Solver is a base class and can’t be created directly:

Solver *solver = Solver(); // Error

Instead, use the create() static function:

auto solver = Solver::create(); // ok

By default this will use the options in the “solver” section of the options, equivalent to:

Options *opts = Options::getRoot()->getSection("solver");
auto solver = Solver::create(opts);

To use a different set of options, for example if there are multiple solvers, use a different option section:

Options *opts = Options::getRoot()->getSection("anothersolver");
auto anothersolver = Solver::create(opts);

Problem specification

The equations to be solved are specified in a PhysicsModel object

class MyProblem : public PhysicsModel {
  protected:
    // This function called once at beginning
    int init(bool restarting) {
      SOLVE_FOR(f);   // Specify variables to solve
      // Set f to initial value if needed
      // otherwise read from input options
      return 0;
    }

    // This function called to evaluate time derivatives
    int rhs(BoutReal t) {
       ddt(f) = 1.0; // Calculate time derivatives
       return 0;
    }
  private:
    Field3D f; // A variable to evolve
}

The init() and rhs() functions must be defined, but there are other functions which can be defined. See PhysicsModel documentation for details.

Create an object, then add to the solver:

MyProblem *prob = MyProblem();
solver->setModel(prob);

Running simulation

To run a calculation

solver->solve();

This will use NOUT and TIMESTEP in the solver options (specified during creation). If these are not present then the global options will be used.

To specify NOUT and TIMESTEP, pass the values to solve:

solver->solve(NOUT, TIMESTEP);

Subclassed by AdamsBashforthSolver, ArkodeSolver, CvodeSolver, EulerSolver, IdaSolver, IMEXBDF2, PetscSolver, PowerSolver, PvodeSolver, RK3SSP, RK4Solver, RKGenericSolver, SlepcSolver, SNESSolver, SplitRK

Public Functions

Solver(Options *opts = nullptr)#

virtual ~Solver() = default#

virtual void setModel(PhysicsModel *model)#: Specify physics model to solve. Currently only one model can be evolved by a Solver.

void addMonitor(Monitor *monitor, MonitorPosition pos = MonitorPosition::FRONT)#

Add a monitor to be called regularly

The frequency at which the monitor is called is set by Monitor::timestep. By default this is every output timestep. When a new Monitor with a smaller timestep is added, the solver attemps to adjust its internal timestep to match, as well as adjusting the timesteps of the current set of Monitors to be multiples of the new timestep. If this is not possible, addMonitor will throw an exception.

Adding new Monitors after the Solver has been initialised is only possible if their timestep is a multiple of the Solver’s timestep. Smaller timesteps will throw an exception.

void removeMonitor(Monitor *monitor)#: Remove a monitor function previously added.

void addTimestepMonitor(TimestepMonitorFunc monitor)#: Add a monitor function to be called every timestep.

void removeTimestepMonitor(TimestepMonitorFunc monitor)#: Remove a previously added timestep monitor.

virtual void add(Field2D &v, const std::string &name, const std::string &description = "")#: Add a variable to be solved. This must be done in the initialisation stage, before the simulation starts.

virtual void add(Field3D &v, const std::string &name, const std::string &description = "")#

virtual void add(Vector2D &v, const std::string &name, const std::string &description = "")#

virtual void add(Vector3D &v, const std::string &name, const std::string &description = "")#

inline virtual bool constraints()#: Returns true if constraints available.

virtual void constraint(Field2D &v, Field2D &C_v, std::string name)#: Add constraint functions (optional). These link a variable v to a control parameter C_v such that v is adjusted to keep C_v = 0.

virtual void constraint(Field3D &v, Field3D &C_v, std::string name)#

virtual void constraint(Vector2D &v, Vector2D &C_v, std::string name)#

virtual void constraint(Vector3D &v, Vector3D &C_v, std::string name)#

inline virtual void setMaxTimestep(MAYBE_UNUSED(BoutReal dt))#: Set a maximum internal timestep (only for explicit schemes)

inline virtual BoutReal getCurrentTimestep()#: Return the current internal timestep.

int solve(int nout = -1, BoutReal timestep = 0.0)#

Start the solver. By default solve() uses options to determine the number of steps and the output timestep. If nout and dt are specified here then the options are not used

Parameters:

nout – [in] Number of output timesteps to run for
timestep – [in] The time between outputs

virtual int init()#: Initialise the solver.

virtual int run() = 0#

Run the solver, calling monitors nout times, at intervals of tstep. This function is called by solve(), and is specific to each solver type

This should probably be protected, since it shouldn’t be called by users.

inline virtual void resetInternalFields()#: Should wipe out internal field vector and reset from current field object data.

inline virtual int n2Dvars() const#: Number of 2D variables. Vectors count as 3.

inline virtual int n3Dvars() const#: Number of 3D variables. Vectors count as 3.

int resetRHSCounter()#: Get and reset the number of calls to the RHS function.

int resetRHSCounter_e()#: Same but for explicit timestep counter - for IMEX.

int resetRHSCounter_i()#: Same but fur implicit timestep counter - for IMEX.

bool splitOperator()#: Test if this solver supports split operators (e.g. implicit/explicit)

virtual void outputVars(Options &output_options, bool save_repeat = true)#

Add evolving variables to output (dump) file or restart file

Parameters:

outputfile – [inout] The Options to add variable to
save_repeat – [in] If true, add variables with time dimension

virtual void readEvolvingVariablesFromOptions(Options &options)#: Copy evolving variables out of options.

std::string getRunID() const#: A unique identifier for this run. Throws if the identifier hasn’t been set yet.

std::string getRunRestartFrom() const#: The run from which this was restarted. Throws if the identifier hasn’t been set yet.

inline int getIterationCounter() const#: Get the number of completed output steps.

inline int incrementIterationCounter()#: Add one to the iteration count, used by BoutMonitor, but could be called by a.

void writeToModelOutputFile(const Options &options)#: Write options to the model’s output file.

Public Members

bool canReset = {false}#

Public Static Functions

static std::unique_ptr<Solver> create(Options *opts = nullptr)#: Create a Solver object. This uses the “type” option in the given Option section to determine which solver type to create.

static std::unique_ptr<Solver> create(const SolverType &type, Options *opts = nullptr)#: Create a Solver object, specifying the type.

static inline void setArgs(int &c, char **&v)#: Pass the command-line arguments. This static function is called by BoutInitialise, and puts references into protected variables. These may then be used by Solvers to control behavior

Public Static Attributes

static constexpr MonitorPosition BACK = MonitorPosition::BACK #

static constexpr MonitorPosition FRONT = MonitorPosition::FRONT #

Protected Functions

int getLocalN()#: Calculate the number of evolving variables on this processor.

template<class T> inline bool contains(const std::vector<VarStr<T>> &vars, const std::string &name)#: Does vars contain a field with name?

inline void add_int_diagnostic(int &i, const std::string &name, const std::string &description = "")#

inline void add_BoutReal_diagnostic(BoutReal &r, const std::string &name, const std::string &description = "")#

int run_rhs(BoutReal t, bool linear = false)#: Run the user’s RHS function.

int run_convective(BoutReal t, bool linear = false)#

Calculate only the convective parts.

NOTE: This calls add_mms_sources.

int run_diffusive(BoutReal t, bool linear = false)#: Calculate only the diffusive parts.

inline void resetIterationCounter(int value = 0)#: Reset the iteration counter.

int call_monitors(BoutReal simtime, int iter, int NOUT)#

Calls all monitor functions

There are two important things to note about how iter is passed along to each monitor:

The solvers all start their iteration numbering from zero, so the initial state is calculated at iter = -1
Secondly, iter is passed along to each monitor relative to that monitor’s period

In practice, this means that each monitor is called like:

monitor->call(solver, simulation_time,
              ((iter + 1) / monitor->period) - 1,
              NOUT / monitor->period);

e.g. for a monitor with period 10, passing iter = 9 will result in it being called with a value of (9 + 1)/10 - 1 == 0

int call_timestep_monitors(BoutReal simtime, BoutReal lastdt)#

bool hasPreconditioner()#: Do we have a user preconditioner?

int runPreconditioner(BoutReal time, BoutReal gamma, BoutReal delta)#: Run the user preconditioner.

bool hasJacobian()#: Do we have a user Jacobian?

int runJacobian(BoutReal time)#: Run the user Jacobian.

void load_vars(BoutReal *udata)#

void load_derivs(BoutReal *udata)#

void save_vars(BoutReal *udata)#

void save_derivs(BoutReal *dudata)#

void set_id(BoutReal *udata)#

Field3D globalIndex(int localStart)#: Returns a Field3D containing the global indices.

inline auto getMonitors() const -> const std::list<MonitorInfo>&#: Get the list of monitors.

inline int getNumberOutputSteps() const#: Get the currently set number of output steps requested.

inline void setNumberOutputSteps(int nout)#: Change the number of requested output steps.

inline BoutReal getOutputTimestep() const#: Get the currently set output timestep.

Protected Attributes

Options *options = {nullptr}#: Settings to use during initialisation (set by constructor)

int NPES = {1}#: Number of processors.

int MYPE = {0}#: This processor’s index.

std::vector<VarStr<Field2D>> f2d#: Vectors of variables to evolve.

std::vector<VarStr<Field3D>> f3d#

std::vector<VarStr<Vector2D>> v2d#

std::vector<VarStr<Vector3D>> v3d#

std::vector<VarStr<int>> diagnostic_int#: Vectors of diagnostic variables to save.

std::vector<VarStr<BoutReal>> diagnostic_BoutReal#

bool has_constraints = {false}#: Can this solver handle constraints? Set to true if so.

bool initialised = {false}#: Has init been called yet?

BoutReal simtime = {0.0}#: Current simulation time.

bool monitor_timestep = {false}#: Should timesteps be monitored?

BoutReal max_dt = {-1.0}#: Maximum internal timestep.

Protected Static Attributes

static int *pargc = nullptr#: Number of command-line arguments.

static char ***pargv = nullptr#: Command-line arguments.

Private Functions

std::string createRunID() const#: Generate a random UUID (version 4) and broadcast it to all processors.

void add_mms_sources(BoutReal t)#

void calculate_mms_error(BoutReal t)#

void pre_rhs(BoutReal t)#: Should be run before user RHS is called.

void post_rhs(BoutReal t)#: Should be run after user RHS is called.

void loop_vars_op(Ind2D i2d, BoutReal *udata, int &p, SOLVER_VAR_OP op, bool bndry)#

Loading data from BOUT++ to/from solver.

Perform an operation at a given Ind2D (jx,jy) location, moving data between BOUT++ and CVODE.

void loop_vars(BoutReal *udata, SOLVER_VAR_OP op)#: Loop over variables and domain. Used for all data operations for consistency.

bool varAdded(const std::string &name)#: Check if a variable has already been added.

BoutReal adjustMonitorPeriods(Monitor *monitor)#: (Possibly) adjust the periods of monitor, and the monitors timesteps, returning the new Solver timestep

void finaliseMonitorPeriods(int &NOUT, BoutReal &output_timestep)#: Fix all the monitor periods based on output_timestep, as well as adjusting NOUT and output_timestep to be consistent

Private Members

std::string run_id = default_run_id #: Randomly generated run ID Initialise with 36 characters so the allocated array is the right size

std::string run_restart_from = "yyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyyy"#: The run from which this was restarted.

bool save_repeat_run_id = {false}#: Save run_id and run_restart_from every output.

int iteration = {0}#: Current iteration (output time-step) number.

int rhs_ncalls = {0}#: Number of calls to the RHS function.

int rhs_ncalls_e = {0}#: Number of calls to the explicit (convective) RHS function.

int rhs_ncalls_i = {0}#: Number of calls to the implicit (diffusive) RHS function.

int default_monitor_period = {1}#: Default sampling rate at which to call monitors - same as output to screen.

BoutReal internal_timestep = {-1}#: timestep - shouldn’t be changed after init is called.

PhysicsModel *model = {nullptr}#: Physics model being evolved.

bool is_nonsplit_model_diffusive = {true}#: Should non-split physics models be treated as diffusive?

bool mms = {false}#: Enable sources and solutions for Method of Manufactured Solutions.

bool mms_initialise = {false}#: Initialise variables to the manufactured solution.

std::list<MonitorInfo> monitors#: List of monitor functions.

std::list<TimestepMonitorFunc> timestep_monitors#: List of timestep monitor functions.

int number_output_steps#: Number of requested output steps.

BoutReal output_timestep#: Requested timestep between outputs.

Private Static Attributes

static constexpr auto default_run_id = "zzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzzz"#: Default value for run_id. Use ‘z’ because it is not a valid hex character, so this is an invalid UUID

Friends

template<class T> inline friend bool operator==(const VarStr<T> &var, const std::string &name)#: Does var represent field name?

template<class T> friend int local_N_sum(int value, const VarStr<T> &f)#

Helper function for getLocalN: return the number of points to evolve in f, plus the accumulator value

If f.evolve_bndry, includes the boundary (NB: not guard!) points

FIXME: This could be a lambda local to getLocalN with an auto argument in C++14

struct MonitorInfo#

#include <solver.hxx>

Helper struct for Monitor and the name of its time dimension.

Public Members

Monitor *monitor#: Non-owning pointer to a monitor instance.

std::string time_dimension#: Name of the time dimension for writing outputs to file. Monitors with the same period as the outputs have plain “t”, all other monitors have an ascending integer suffix, starting with the front monitor. WARNING! This could change if there are different monitors added on subsequent restarts, or if they are added in a different order.

template<class T> struct VarStr#

#include <solver.hxx>

A structure to hold an evolving variable.

Public Members

bool constraint = {false}#

T *var = {nullptr}#: Does F_var represent a constraint?

T *F_var = {nullptr}#: The evolving variable.

std::unique_ptr<T> MMS_err = {nullptr}#: The time derivative or constraint on var.

CELL_LOC location = {CELL_DEFAULT}#: Error for MMS.

bool covariant = {false}#: For fields and vector components.

bool evolve_bndry = {false}#: For vectors.

std::string name#: Are the boundary regions being evolved?

std::string description = {""}#: Name of the variable.