File naulin_laplace.hxx#

class LaplaceNaulin : public Laplacian#
#include <naulin_laplace.hxx>

Solves the 2D Laplacian equation.

Public Functions

LaplaceNaulin(Options *opt = NULL, const CELL_LOC loc = CELL_CENTRE, Mesh *mesh_in = nullptr, Solver *solver = nullptr)#
~LaplaceNaulin() = default#
inline virtual void setCoefA(const Field2D &val) override#

Set coefficients for inversion. Re-builds matrices if necessary.

inline virtual void setCoefA(const Field3D &val) override#
inline virtual void setCoefC(const Field2D &val) override#
inline virtual void setCoefC(const Field3D &val) override#
inline virtual void setCoefC1(const Field3D &val) override#
inline virtual void setCoefC1(const Field2D &val) override#
inline virtual void setCoefC2(const Field3D &val) override#
inline virtual void setCoefC2(const Field2D &val) override#
inline virtual void setCoefD(const Field3D &val) override#
inline virtual void setCoefD(const Field2D &val) override#
inline virtual void setCoefEx(const Field2D &val) override#
inline virtual void setCoefEz(const Field2D &val) override#
inline virtual bool uses3DCoefs() const override#

Does this solver use Field3D coefficients (true) or only their DC component (false)

inline virtual FieldPerp solve(const FieldPerp &b) override#
inline virtual FieldPerp solve(const FieldPerp &b, const FieldPerp &x0) override#
virtual Field3D solve(const Field3D &b, const Field3D &x0) override#

Performs the laplacian inversion y-slice by y-slice

Parameters:
  • b[in] All the y-slices of b_slice, which is the right hand side of the equation A*x_slice = b_slice

  • x0[in] All the y-slices of the variable eventually used to set BC

Returns:

x All the y-slices of x_slice in the equation A*x_slice = b_slice

inline virtual Field3D solve(const Field3D &b) override#
inline virtual void setGlobalFlags(int f) override#
inline virtual void setInnerBoundaryFlags(int f) override#
inline virtual void setOuterBoundaryFlags(int f) override#
inline BoutReal getMeanIterations() const#
inline void resetMeanIterations()#
void outputVars(Options &output_options, const std::string &time_dimension) const override#

Private Functions

LaplaceNaulin(const LaplaceNaulin&)#
LaplaceNaulin &operator=(const LaplaceNaulin&)#
void copy_x_boundaries(Field3D &x, const Field3D &x0, Mesh *mesh)#

Copy the boundary guard cells from the input ‘initial guess’ x0 into x. These may be used to set non-zero-value boundary conditions

Private Members

Field3D Acoef#
Field3D C1coef#
Field3D C2coef#
Field3D Dcoef#
std::unique_ptr<Laplacian> delp2solver = {nullptr}#

Laplacian solver used to solve the equation with constant-in-z coefficients.

BoutReal rtol#

Solver tolerances.

BoutReal atol#
int maxits#

Maximum number of iterations.

BoutReal initial_underrelax_factor = {1.}#

Initial choice for under-relaxation factor, should be greater than 0 and less than or equal to 1. Value of 1 means no underrelaxation

BoutReal naulinsolver_mean_its#

Mean number of iterations taken by the solver.

BoutReal naulinsolver_mean_underrelax_counts = {0.}#

Mean number of times the underrelaxation factor is reduced.

int ncalls#

Counter for the number of times the solver has been called.