File petsc_laplace.hxx#

class LaplacePetsc : public Laplacian#

Public Functions

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

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

inline virtual void setCoefC(const Field2D &val) override#
inline virtual void setCoefC1(const Field2D &val) override#
inline virtual void setCoefC2(const Field2D &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 void setCoefA(const Field3D &val) override#
inline virtual void setCoefC(const Field3D &val) override#
inline virtual void setCoefC1(const Field3D &val) override#
inline virtual void setCoefC2(const Field3D &val) override#
inline virtual void setCoefD(const Field3D &val) override#
inline virtual void setCoefEx(const Field3D &val) override#
inline virtual void setCoefEz(const Field3D &val) override#
virtual FieldPerp solve(const FieldPerp &b) override#
virtual FieldPerp solve(const FieldPerp &b, const FieldPerp &x0) override#

Solves Ax=b for x given a b and an initial guess for x (x0)

This function will:

  1. Set the matrix element of the matrix A, used to solve Ax=b (this includes setting the values for the bounary condition)

  2. Solve the matrix Ax = b

Parameters:
  • b[in] The RHS of the equation Ax=b. This is an y-slice of the original field. The field wil be flattened to an 1D array in order to write the equation on the form Ax=b

  • x0[in] The initial guess for the solver. May also contain the boundary condition if flag 32 - INVERT_SET is set

Returns:

sol The solution x of the problem Ax=b.

int precon(Vec x, Vec y)#

Preconditioner function.

void setCoefA(const Field2D &val) = 0

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

inline void setCoefA(const Field3D &val)
inline void setCoefA(BoutReal r)#
void setCoefC(const Field2D &val) = 0
inline void setCoefC(const Field3D &val)
inline void setCoefC(BoutReal r)#
inline void setCoefC1(const Field2D &val)
inline void setCoefC1(const Field3D &val)
inline void setCoefC1(BoutReal r)#
inline void setCoefC2(const Field2D &val)
inline void setCoefC2(const Field3D &val)
inline void setCoefC2(BoutReal r)#
void setCoefD(const Field2D &val) = 0
inline void setCoefD(const Field3D &val)
inline void setCoefD(BoutReal r)#
void setCoefEx(const Field2D &val) = 0
inline void setCoefEx(const Field3D &val)
inline void setCoefEx(BoutReal r)#
void setCoefEz(const Field2D &val) = 0
inline void setCoefEz(const Field3D &val)
inline void setCoefEz(BoutReal r)#
FieldPerp solve(const FieldPerp &b) = 0
Field3D solve(const Field3D &b)#
Field2D solve(const Field2D &b)#
inline FieldPerp solve(const FieldPerp &b, const FieldPerp &x0)
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

Field2D solve(const Field2D &b, const Field2D &x0)#

Private Functions

void Element(int i, int x, int z, int xshift, int zshift, PetscScalar ele, Mat &MatA)#

Sets the elements of the matrix A, which is used to solve the problem Ax=b.

Parameters:
  • i[in] The row of the PETSc matrix

  • x[in] Local x index of the mesh

  • z[in] Local z index of the mesh

  • xshift[in] The shift in rows from the index x

  • zshift[in] The shift in columns from the index z

  • ele[in] Value of the element

  • MatA[in] The matrix A used in the inversion

  • MatA[out] The matrix A used in the inversion

void Coeffs(int x, int y, int z, BoutReal &A1, BoutReal &A2, BoutReal &A3, BoutReal &A4, BoutReal &A5)#

Set the matrix components of A in Ax=b, solving D*Laplace_perp(x) + (1/C1)Grad_perp(C2)*Grad_perp(x) + Ax = B

Note

“A” in the equation above is not added here. For calculations of the coefficients, please refer to the user manual.

Parameters:
  • x[in] The current x index

  • y[in] The current y index

  • z[in] The current y index

  • coef1[in] Placeholder for convenient variable used to set matrix (see manual for details)

  • coef2[in] Convenient variable used to set matrix (see manual for details)

  • coef3[in] Placeholder for convenient variable used to set matrix (see manual for details)

  • coef4[in] Placeholder for convenient variable used to set matrix (see manual for details)

  • coef5[in] Placeholder for convenient variable used to set matrix (see manual for details)

  • coef1[out] Convenient variable used to set matrix (see manual for details)

  • coef2[out] Convenient variable used to set matrix (see manual for details)

  • coef3[out] Convenient variable used to set matrix (see manual for details)

  • coef4[out] Convenient variable used to set matrix (see manual for details)

  • coef5[out] Convenient variable used to set matrix (see manual for details)

void vecToField(Vec x, FieldPerp &f)#
void fieldToVec(const FieldPerp &f, Vec x)#

Private Members

Field3D A#
Field3D C1#
Field3D C2#
Field3D D#
Field3D Ex#
Field3D Ez#
bool issetD#
bool issetC#
bool issetE#
FieldPerp sol#
int Istart#
int Iend#
int meshx#
int meshz#
int size#
int localN#
MPI_Comm comm#
Mat MatA#
Vec xs#
Vec bs#
KSP ksp#
Options *opts#
std::string ksptype#

KSP solver type.

std::string pctype#

Preconditioner type.

BoutReal richardson_damping_factor#
BoutReal chebyshev_max#
BoutReal chebyshev_min#
int gmres_max_steps#
BoutReal rtol#
BoutReal atol#
BoutReal dtol#
int maxits#
bool direct#
bool fourth_order#
PetscLib lib#
bool rightprec#
std::unique_ptr<Laplacian> pcsolve#
int implemented_flags#
int implemented_boundary_flags#