PETScKrylovSolver.h

Note

The documentation on this page was automatically extracted from the DOLFIN C++ code and may need to be edited or expanded.

class PETScKrylovSolver

Parent class(es)

This class implements Krylov methods for linear systems of the form Ax = b. It is a wrapper for the Krylov solvers of PETSc.

enum class norm_type

Norm types used in convergence testing. Not all solvers types support all norm types (see http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/KSP/KSPSetNormType.html). Note that ‘default’ is a reserved keyword, so we use ‘default_norm’

PETScKrylovSolver(MPI_Comm comm, std::string method = "default", std::string preconditioner = "default")

Create Krylov solver for a particular method and named preconditioner

PETScKrylovSolver(std::string method = "default", std::string preconditioner = "default")

Create Krylov solver for a particular method and named preconditioner

PETScKrylovSolver(MPI_Comm comm, std::string method, std::shared_ptr<PETScPreconditioner> preconditioner)

Create Krylov solver for a particular method and PETScPreconditioner (shared_ptr version)

PETScKrylovSolver(std::string method, std::shared_ptr<PETScPreconditioner> preconditioner)

Create Krylov solver for a particular method and PETScPreconditioner (shared_ptr version)

explicit PETScKrylovSolver(KSP ksp)

Create solver wrapper of a PETSc KSP object

void set_operator(std::shared_ptr<const GenericLinearOperator> A)

Set operator (matrix)

void set_operator(const PETScBaseMatrix &A)

Set operator (PETScMatrix). This is memory-safe as PETSc will increase the reference count to the underlying PETSc object.

void set_operators(std::shared_ptr<const GenericLinearOperator> A, std::shared_ptr<const GenericLinearOperator> P)

Set operator (matrix) and preconditioner matrix

void set_operators(const PETScBaseMatrix &A, const PETScBaseMatrix &P)

Set operator and preconditioner matrix (PETScMatrix). This is memory-safe as PETSc will increase the reference count to the underlying PETSc objects.

std::size_t solve(GenericVector &x, const GenericVector &b)

Solve linear system Ax = b and return number of iterations

std::size_t solve(GenericVector &x, const GenericVector &b, bool transpose)

Solve linear system Ax = b and return number of iterations (A^t x = b if transpose is true)

std::size_t solve(PETScVector &x, const PETScVector &b, bool transpose = false)

Solve linear system Ax = b and return number of iterations (A^t x = b if transpose is true)

std::size_t solve(const GenericLinearOperator &A, GenericVector &x, const GenericVector &b)

Solve linear system Ax = b and return number of iterations

void set_nonzero_guess(bool nonzero_guess)

Use nonzero initial guess for solution function (nonzero_guess=true, the solution vector x will not be zeroed before the solver starts)

void set_reuse_preconditioner(bool reuse_pc)

Reuse preconditioner if true, even if matrix operator changes (by default preconditioner will be re-built if the matrix changes)

void set_tolerances(double relative, double absolute, double diverged, int max_iter)

Set tolerances (relative residual, alsolute residial, maximum number of iterations)

void set_norm_type(norm_type type)

Set norm type used in convergence testing - not all solvers types support all norm types

norm_type get_norm_type() const

Get norm type used in convergence testing

void monitor(bool monitor_convergence)

Monitor residual at each iteration

void set_options_prefix(std::string options_prefix)

Sets the prefix used by PETSc when searching the PETSc options database

std::string get_options_prefix() const

Returns the prefix used by PETSc when searching the PETSc options database

void set_from_options() const

Set options from PETSc options database

std::string str(bool verbose) const

Return informal string representation (pretty-print)

MPI_Comm mpi_comm() const

Return MPI communicator

KSP ksp() const

Return PETSc KSP pointer

static std::map<std::string, std::string> methods()

Return a list of available solver methods

static std::map<std::string, std::string> preconditioners()

Return a list of available named preconditioners

static Parameters default_parameters()

Default parameter values

std::string parameter_type() const

Return parameter type: “krylov_solver” or “lu_solver”

void set_dm(DM dm)

Set the DM

void set_dm_active(bool val)

Activate/deactivate DM