File invert_laplace.hxx#
Perpendicular Laplacian inversion using FFT and Tridiagonal solver
Equation solved is: \(d*\nabla^2_\perp x + (1/c)\nabla_perp c\cdot\nabla_\perp x + a x = b\)
Where a, c and d are functions of x and y only (not z)
Copyright 2010 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/.
Defines
-
PVEC_REAL_MPI_TYPE
Typedefs
-
template<class DerivedType>
using RegisterLaplace = LaplaceFactory::RegisterInFactory<DerivedType># Simpler name for Factory registration helper class
Usage:
#include <bout/laplacefactory.hxx> namespace { RegisterLaplace<MyLaplace> registerlaplacemine("mylaplace"); }
Functions
Variables
-
constexpr auto LAPLACE_SPT = "spt"#
-
constexpr auto LAPLACE_TRI = "tri"#
-
constexpr auto LAPLACE_BAND = "band"#
-
constexpr auto LAPLACE_PETSC = "petsc"#
-
constexpr auto LAPLACE_PETSCAMG = "petscamg"#
-
constexpr auto LAPLACE_PETSC3DAMG = "petsc3damg"#
-
constexpr auto LAPLACE_HYPRE3D = "hypre3d"#
-
constexpr auto LAPLACE_CYCLIC = "cyclic"#
-
constexpr auto LAPLACE_MULTIGRID = "multigrid"#
-
constexpr auto LAPLACE_NAULIN = "naulin"#
-
constexpr auto LAPLACE_IPT = "ipt"#
-
constexpr auto LAPLACE_PCR = "pcr"#
-
constexpr auto LAPLACE_PCR_THOMAS = "pcr_thomas"#
-
constexpr int INVERT_DC_GRAD = 1#
Zero-gradient for DC (constant in Z) component. Default is zero value.
-
constexpr int INVERT_AC_GRAD = 2#
Zero-gradient for AC (non-constant in Z) component. Default is zero value.
-
constexpr int INVERT_AC_LAP = 4#
Use zero-laplacian (decaying solution) to AC component.
-
constexpr int INVERT_SYM = 8#
Use symmetry to enforce either zero-value or zero-gradient.
-
constexpr int INVERT_RHS = 32#
Use input value in RHS boundary.
-
constexpr int INVERT_DC_LAP = 64#
Use zero-laplacian solution for DC component.
-
constexpr int INVERT_BNDRY_ONE = 128#
Only use one boundary point.
-
constexpr int INVERT_DC_GRADPAR = 256#
-
constexpr int INVERT_DC_GRADPARINV = 512#
-
constexpr int INVERT_IN_CYLINDER = 1024#
For use in cylindrical coordiate system.
-
constexpr int INVERT_ZERO_DC = 1#
Zero the DC (constant in Z) component of the solution.
-
constexpr int INVERT_START_NEW = 2#
Iterative method start from solution=0. Has no effect for direct solvers.
-
constexpr int INVERT_BOTH_BNDRY_ONE = 4#
Sets the width of the boundaries to 1.
-
constexpr int INVERT_4TH_ORDER = 8#
Use band solver for 4th order in x.
-
constexpr int INVERT_KX_ZERO = 16#
Zero the kx=0, n = 0 component.
-
class LaplaceFactory : public Factory<Laplacian, LaplaceFactory, Options*, CELL_LOC, Mesh*, Solver*>#
Public Functions
Public Static Attributes
-
static constexpr auto type_name = "Laplacian"#
-
static constexpr auto section_name = "laplace"#
-
static constexpr auto option_name = "type"#
-
static constexpr auto default_type = LAPLACE_CYCLIC#
-
static constexpr auto type_name = "Laplacian"#
-
class Laplacian#
- #include <invert_laplace.hxx>
Base class for Laplacian inversion.
Subclassed by LaplaceHypre3d, LaplaceNaulin, LaplacePCR, LaplacePCR_THOMAS, LaplacePetsc, LaplacePetsc3dAmg, LaplaceSPT, LaplaceSerialTri
Public Functions
-
Laplacian(Options *options = nullptr, const CELL_LOC loc = CELL_CENTRE, Mesh *mesh_in = nullptr, Solver *solver = nullptr)#
Laplacian inversion initialisation. Called once at the start to get settings.
-
virtual ~Laplacian() = default#
-
virtual void setCoefA(const Field2D &val) = 0#
Set coefficients for inversion. Re-builds matrices if necessary.
-
inline virtual void setGlobalFlags(int f)#
-
inline virtual void setInnerBoundaryFlags(int f)#
-
inline virtual void setOuterBoundaryFlags(int f)#
-
inline virtual int getGlobalFlags() const#
-
inline virtual int getInnerBoundaryFlags() const#
-
inline virtual int getOuterBoundaryFlags() const#
-
inline virtual bool uses3DCoefs() const#
Does this solver use Field3D coefficients (true) or only their DC component (false)
-
virtual Field3D solve(const Field3D &b, const Field3D &x0)#
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
-
void tridagCoefs(int jx, int jy, int jz, dcomplex &a, dcomplex &b, dcomplex &c, const Field2D *ccoef = nullptr, const Field2D *d = nullptr, CELL_LOC loc = CELL_DEFAULT)#
Coefficients in tridiagonal inversion.
-
inline void outputVars(Options &output_options) const#
Add any output variables to
output_options
, for example, performance information, with optional name for the time dimension
-
inline virtual void outputVars([[maybe_unused]] Options &output_options, [[maybe_unused]] const std::string &time_dimension) const#
-
inline void savePerformance(Solver &solver)#
Register performance monitor with
solver
, prefix output withOptions
section name
-
void savePerformance(Solver &solver, const std::string &name)#
Register performance monitor that is call every timestep with
solver
, prefix output withname
. Call this function from yourPhysicsModel::init
to get time-dependent performance information, or calloutputVars
directly in non-model code to get the information at that point in time.To use this for a Laplacian implementation, override
outputVars(Options&, const std::string&)
, and add whatever information you would like to save to theOptions
argument. This will then be called automatically byLaplacianMonitor
.See
LaplaceNaulin::outputVars
for an example.
Public Static Functions
-
static inline std::unique_ptr<Laplacian> create(Options *opts = nullptr, const CELL_LOC location = CELL_CENTRE, Mesh *mesh_in = nullptr, Solver *solver = nullptr)#
Create a new Laplacian solver
- Parameters:
opt – [in] The options section to use. By default “laplace” will be used
-
static void cleanup()#
Frees all memory.
Protected Functions
-
inline bool isGlobalFlagSet(int flag) const#
Return true if global/default
flag
is set.
-
inline bool isInnerBoundaryFlagSet(int flag) const#
Return true if
flag
is set for the inner boundary condition.
-
inline bool isOuterBoundaryFlagSet(int flag) const#
Return true if
flag
is set for the outer boundary condition.
-
bool isInnerBoundaryFlagSetOnFirstX(int flag) const#
Return true if
flag
is set for the inner boundary condition and this is the first proc in X direction
-
bool isOuterBoundaryFlagSetOnLastX(int flag) const#
Return true if
flag
is set for the outer boundary condition and this the last proc in X direction
-
inline void tridagCoefs(int jx, int jy, BoutReal kwave, dcomplex &a, dcomplex &b, dcomplex &c, const Field2D *ccoef = nullptr, const Field2D *d = nullptr, CELL_LOC loc = CELL_DEFAULT)#
-
void tridagCoefs(int jx, int jy, BoutReal kwave, dcomplex &a, dcomplex &b, dcomplex &c, const Field2D *c1coef, const Field2D *c2coef, const Field2D *d, CELL_LOC loc = CELL_DEFAULT)#
-
inline void tridagMatrix(dcomplex *avec, dcomplex *bvec, dcomplex *cvec, dcomplex *bk, int jy, int kz, BoutReal kwave, const Field2D *a, const Field2D *ccoef, const Field2D *d, bool includeguards = true, bool zperiodic = true)#
-
void tridagMatrix(dcomplex *avec, dcomplex *bvec, dcomplex *cvec, dcomplex *bk, int jy, int kz, BoutReal kwave, const Field2D *a, const Field2D *c1coef, const Field2D *c2coef, const Field2D *d, bool includeguards = true, bool zperiodic = true)#
Set the matrix components of A in Ax=b
This function will
Calling tridagCoef, solving
D*Laplace_perp(x) + (1/C1)Grad_perp(C2)*Grad_perp(x) + Ax = B
for each fourier component
Set the boundary conditions by setting the first and last rows properly
- Parameters:
avec – [in] Lower diagonal of the tridiagonal matrix. DO NOT CONFUSE WITH “A”
bvec – [in] The main diagonal
cvec – [in] The upper diagonal. DO NOT CONFUSE WITH “C” (called ccoef here)
bk – [in] The b in Ax = b
jy – [in] Index of the current y-slice
kz – [in] The mode number index
kwave – [in] The mode number (different from kz only if we are taking a part of the z-domain [and not from 0 to 2*pi])
global_flags – [in] Global flags of the inversion
inner_boundary_flags – [in] Flags used to set the inner boundary
outer_boundary_flags – [in] Flags used to set the outer boundary
a – [in] A in the equation above. DO NOT CONFUSE WITH avec
c1coef – [in] C1 in the equation above. DO NOT CONFUSE WITH cvec
c2coef – [in] C2 in the equation above. DO NOT CONFUSE WITH cvec
d – [in] D in the equation above
includeguards – [in] Whether or not the guard points in x should be used
zperiodic – [in] Whether a special case for kx=kz=0 is needed
avec – [out] Lower diagonal of the tridiagonal matrix. DO NOT CONFUSE WITH “A”
bvec – [out] The main diagonal
cvec – [out] The upper diagonal. DO NOT CONFUSE WITH “C” (called ccoef here)
Protected Attributes
-
bool async_send#
If true, use asyncronous send in parallel algorithms.
-
int maxmode#
The maximum Z mode to solve for.
-
bool low_mem#
If true, reduce the amount of memory used.
-
bool all_terms#
applies to Delp2 operator and laplacian inversion
-
bool nonuniform#
Non-uniform mesh correction.
-
bool include_yguards#
solve in y-guard cells, default true.
-
int extra_yguards_lower#
exclude some number of points at the lower boundary, useful for staggered grids or when boundary conditions make inversion redundant
-
int extra_yguards_upper#
exclude some number of points at the upper boundary, useful for staggered grids or when boundary conditions make inversion redundant
-
Coordinates *coords#
Coordinates object, so we only have to call localmesh->getCoordinates(location) once
Private Members
-
int global_flags#
Default flags.
-
int inner_boundary_flags#
Flags to set inner boundary condition.
-
int outer_boundary_flags#
Flags to set outer boundary condition.
-
std::string performance_name#
Name for writing performance infomation; default taken from constructing
Options
section
-
LaplacianMonitor monitor = {*this}#
-
Laplacian(Options *options = nullptr, const CELL_LOC loc = CELL_CENTRE, Mesh *mesh_in = nullptr, Solver *solver = nullptr)#