File physicsmodel.hxx#

Base class for Physics Models.

Changelog:

2013-08 Ben Dudson benjamin.dudson@york.ac.uk

  • Initial version

Copyright 2013 B.D.Dudson

Contact: Ben Dudson, bd512@york.ac.uk

This file is part of BOUT++.

BOUT++ is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

BOUT++ is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with BOUT++. If not, see http://www.gnu.org/licenses/.

Defines

BOUTMAIN(ModelClass)#

Macro to define a simple main() which creates the given model and runs it. This should be sufficient for most use cases, but a user can define their own main() function if needed.

Example

class MyModel : public PhysicsModel { .. };

BOUTMAIN(MyModel);

SOLVE_FOR1(var)#

Macro to replace solver->add, passing variable name.

SOLVE_FOR2(var1, var2)#
SOLVE_FOR3(var1, var2, var3)#
SOLVE_FOR4(var1, var2, var3, var4)#
SOLVE_FOR5(var1, var2, var3, var4, var5)#
SOLVE_FOR6(var1, var2, var3, var4, var5, var6)#
SOLVE_FOR(...)#

Add fields to the solver. This should accept up to ten arguments

SAVE_ONCE1(var)#

Write this variable once to the grid file.

SAVE_ONCE2(var1, var2)#
SAVE_ONCE3(var1, var2, var3)#
SAVE_ONCE4(var1, var2, var3, var4)#
SAVE_ONCE5(var1, var2, var3, var4, var5)#
SAVE_ONCE6(var1, var2, var3, var4, var5, var6)#
SAVE_ONCE(...)#
SAVE_REPEAT1(var)#

Write this variable every timestep.

SAVE_REPEAT2(var1, var2)#
SAVE_REPEAT3(var1, var2, var3)#
SAVE_REPEAT4(var1, var2, var3, var4)#
SAVE_REPEAT5(var1, var2, var3, var4, var5)#
SAVE_REPEAT6(var1, var2, var3, var4, var5, var6)#
SAVE_REPEAT(...)#
class PhysicsModel#
#include <physicsmodel.hxx>

Base class for physics models

Public Types

using preconfunc = int (PhysicsModel::*)(BoutReal t, BoutReal gamma, BoutReal delta)#
using jacobianfunc = int (PhysicsModel::*)(BoutReal t)#
template<class Model, typename = std::enable_if_t<std::is_base_of_v<PhysicsModel, Model>>>
using ModelPreconFunc = int (Model::*)(BoutReal t, BoutReal gamma, BoutReal delta)#
template<class Model, typename = std::enable_if_t<std::is_base_of_v<PhysicsModel, Model>>>
using ModelJacobianFunc = int (Model::*)(BoutReal t)#

Public Functions

PhysicsModel()#
inline PhysicsModel(Mesh *mesh_, bool output_enabled_, bool restart_enabled_)#
virtual ~PhysicsModel() = default#
void initialise(Solver *s)#

Initialise the model, calling the init() and postInit() methods

Note: this is usually only called by the Solver

int runRHS(BoutReal time, bool linear = false)#

Run the RHS function, to calculate the time derivatives

Input

The system state should be in the evolving variables

The time derivatives will be put in the ddt() variables

Returns a flag: 0 indicates success, non-zero an error flag

Parameters:
  • time[in] The simulation time

  • linear[in] True if the function can be a linearised form

bool splitOperator()#

True if this model uses split operators

int runConvective(BoutReal time, bool linear = false)#

Run the convective (usually explicit) part of the model

Parameters:
  • time[in] The simulation time

  • linear[in] True if the function can be a linearised form

int runDiffusive(BoutReal time, bool linear = false)#

Run the diffusive (usually implicit) part of the model

Parameters:
  • time[in] The simulation time

  • linear[in] True if the function can be a linearised form

bool hasPrecon()#

True if a preconditioner has been defined

int runPrecon(BoutReal t, BoutReal gamma, BoutReal delta)#

Run the preconditioner. The system state should be in the evolving variables, and the vector to be solved in the ddt() variables. The result will be put in the ddt() variables.

Note: this is usually only called by the Solver

bool hasJacobian()#

True if a Jacobian function has been defined

int runJacobian(BoutReal t)#

Run the Jacobian-vector multiplication function

Note: this is usually only called by the Solver

inline int runTimestepMonitor(BoutReal simtime, BoutReal dt)#
void writeOutputFile(const Options &options)#

Write options to output_file

void writeOutputFile(const Options &options, const std::string &time_dimension)#

Write variables with time_dimension from options to output_file

void finishOutputTimestep() const#

Finish the output for this timestep, verifying all evolving variables have the correct length

Public Members

Mesh *mesh = {nullptr}#
bout::DataFileFacade dump = {}#
bout::DataFileFacade restart = {}#

Protected Functions

virtual int init(bool restarting) = 0#

This function is called once by the solver at the start of a simulation.

A valid PhysicsModel must implement this function

Variables should be read from the inputs, and the variables to be evolved should be specified.

virtual int postInit(bool restarting)#

Post-initialise. This reads the restart file

Parameters:

restarting[in] If true, will load state from restart file

inline virtual int rhs(BoutReal t)#

This function is called by the time integration solver at least once per time step.

Variables being evolved will be set by the solver before the call, and this function must calculate and set the time-derivatives.

By default this function just returns an error, which will stop the simulation.

An optional second argument can be used, linear, which is set to true when the rhs() function can be linearised. This is used in e.g. linear iterative solves.

inline virtual int rhs(BoutReal t, bool linear)#
virtual void outputVars(Options &options)#

Output additional variables other than the evolving variables.

virtual void restartVars(Options &options)#

Add additional variables other than the evolving variables to the restart files.

inline virtual int convective(BoutReal t)#
inline virtual int convective(BoutReal t, bool linear)#
inline virtual int diffusive(BoutReal t)#
inline virtual int diffusive(BoutReal t, bool linear)#
inline virtual int outputMonitor(BoutReal simtime, int iter, int NOUT)#

Implemented by user code to monitor solution at output times

inline virtual int timestepMonitor(BoutReal simtime, BoutReal dt)#

Timestep monitor. If enabled by setting solver:monitor_timestep=true then this function is called every internal timestep.

inline void setSplitOperator(bool split = true)#

Specify that this model is split into a convective and diffusive part.

inline void setPrecon(preconfunc pset)#

Specify a preconditioner function.

template<class Model>
inline void setPrecon(ModelPreconFunc<Model> preconditioner)#
inline void setJacobian(jacobianfunc jset)#

Specify a Jacobian-vector multiply function.

template<class Model>
inline void setJacobian(ModelJacobianFunc<Model> jacobian)#
void bout_solve(Field2D &var, const char *name, const std::string &description = "")#

Specify a variable for the solver to evolve

Note that the variable must not be destroyed (e.g. go out of scope) after this call, since a pointer to var is stored in the solver.

To evolve the state, the solver will set var, and the user-supplied rhs() function should calculate ddt(var).

Parameters:
  • var[in] The variable to evolve

  • name[in] The name to use for variable initialisation and output

void bout_solve(Field3D &var, const char *name, const std::string &description = "")#
void bout_solve(Vector2D &var, const char *name, const std::string &description = "")#
void bout_solve(Vector3D &var, const char *name, const std::string &description = "")#
inline Options &readFromRestartFile(const std::string &name)#

Helper function for reading from restart_options.

void writeRestartFile()#

Write the restart file to disk now.

void writeOutputFile()#

Write the output file to disk now.

bool bout_constrain(Field3D &var, Field3D &F_var, const char *name)#

Specify a constrained variable var, which will be adjusted to make F_var equal to zero. If the solver does not support constraints then this will throw an exception

Parameters:
  • var[in] The variable the solver should modify

  • F_var[in] The control variable, which the user will set

  • name[in] The name to use for initialisation and output

Protected Attributes

Solver *solver = {nullptr}#

This is set by a call to initialise, and can be used by models to specify evolving variables.

Private Members

Options output_options#

State for outputs.

std::unique_ptr<bout::OptionsIO> output_file#

File to write the outputs to.

bool output_enabled = {true}#

Should we write output files.

Options restart_options#

Stores the state for restarting.

std::unique_ptr<bout::OptionsIO> restart_file#

File to write the restart-state to.

bool restart_enabled = {true}#

Should we write restart files.

bool splitop = {false}#

Split operator model?

preconfunc userprecon = {nullptr}#

Pointer to user-supplied preconditioner function.

jacobianfunc userjacobian = {nullptr}#

Pointer to user-supplied Jacobian-vector multiply function.

bool initialised = {false}#

True if model already initialised.

PhysicsModelMonitor modelMonitor = {this}#

write restarts and pass outputMonitor method inside a Monitor subclass

class PhysicsModelMonitor : public Monitor#
#include <physicsmodel.hxx>

Monitor class for PhysicsModel

Public Functions

PhysicsModelMonitor() = delete#
inline PhysicsModelMonitor(PhysicsModel *model)#
virtual int call(Solver *solver, BoutReal simtime, int iter, int nout) override#

Callback function for the solver, called after timestep_ has passed

Parameters:
  • solver[in] The solver calling this monitor

  • time[in] The current simulation time

  • iter[in] The current simulation iteration

  • nout[in] The total number of iterations for this simulation

Returns:

non-zero if simulation should be stopped

Private Members

PhysicsModel *model#
namespace bout

Provides access to the Hypre library, handling initialisation and finalisation.

Usage

#include <bout/hyprelib.hxx>

class MyClass { public:

private: HypreLib lib; };

This will then automatically initialise Hypre the first time an object is created, and finalise it when the last object is destroyed.

Copyright 2012 B.D.Dudson, S.Farley, M.V.Umansky, X.Q.Xu

Contact: Ben Dudson, bd512@york.ac.uk

This file is part of BOUT++.

BOUT++ is free software: you can redistribute it and/or modify it under the terms of the GNU Lesser General Public License as published by the Free Software Foundation, either version 3 of the License, or (at your option) any later version.

BOUT++ is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License along with BOUT++. If not, see http://www.gnu.org/licenses/.

Information about the version of BOUT++

The build system will update this file on every commit, which may result in files that include it getting rebuilt. Therefore it should be included in as few places as possible

Information about the version of BOUT++

The build system will update this file at configure-time

SNB model

class DataFileFacade#
#include <physicsmodel.hxx>

Stop-gap shim for DataFile: allows PhysicsModels to still save outputs using the old DataFile interface without upgrading to the new OptionsNetCDF.

Stores pointers to objects in a vector, using a variant, which can then be pushed into an Options and written to file using the default implementaion of PhysicsModel::outputVars

Public Types

using ValueType = bout::utils::variant<bool*, int*, BoutReal*, std::string*, Field2D*, Field3D*, FieldPerp*, Array<BoutReal>*, Matrix<BoutReal>*, Tensor<BoutReal>*>#

This is Options::ValueType but using pointers rather than values.

Public Functions

template<class T>
inline void addRepeat(T &value, const std::string &name)#
template<class T>
inline void addOnce(T &value, const std::string &name)#
template<class T>
inline void add(T &value, const std::string &name, bool save_repeat = false)#
inline void addRepeat(ValueType value, const std::string &name)#
inline void addOnce(ValueType value, const std::string &name)#
void add(ValueType value, const std::string &name, bool save_repeat = false)#
void add(Vector2D *value, const std::string &name, bool save_repeat = false)#
void add(Vector3D *value, const std::string &name, bool save_repeat = false)#
inline const std::vector<StoredValue> &getData()#

Private Members

std::vector<StoredValue> data#
struct StoredValue#

Helper struct to save enough information so that we can save an object to file later

Public Functions

inline StoredValue(std::string name_, ValueType value_, bool repeat_)#

Public Members

std::string name#
ValueType value#
bool repeat#
struct OptionsConversionVisitor#

Public Functions

inline OptionsConversionVisitor(Options &options_, std::string name_)#
template<class T>
inline void operator()(T *value)#

Private Members

Options &options#
std::string name#