![]() |
Reference documentation for deal.II version GIT 9c182271f7 2023-03-28 14:30:01+00:00
|
#include <deal.II/lac/petsc_solver.h>
Classes | |
struct | AdditionalData |
struct | SolverDataMUMPS |
Public Member Functions | |
SparseDirectMUMPS (SolverControl &cn, const AdditionalData &data=AdditionalData()) | |
SparseDirectMUMPS (SolverControl &cn, const MPI_Comm &mpi_communicator, const AdditionalData &data=AdditionalData()) | |
void | solve (const MatrixBase &A, VectorBase &x, const VectorBase &b) |
void | set_symmetric_mode (const bool flag) |
void | solve (const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionBase &preconditioner) |
virtual void | reset () |
void | set_prefix (const std::string &prefix) |
SolverControl & | control () const |
void | initialize (const PreconditionBase &preconditioner) |
Protected Member Functions | |
virtual void | set_solver_type (KSP &ksp) const override |
Protected Attributes | |
const AdditionalData | additional_data |
SolverControl & | solver_control |
std::string | prefix_name |
Static Private Member Functions | |
static PetscErrorCode | convergence_test (KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control) |
Private Attributes | |
std::unique_ptr< SolverDataMUMPS > | solver_data |
bool | symmetric_mode |
An implementation of the solver interface using the sparse direct MUMPS solver through PETSc. This class has the usual interface of all other solver classes but it is of course different in that it doesn't implement an iterative solver. As a consequence, things like the SolverControl object have no particular meaning here.
MUMPS allows to make use of symmetry in this matrix. In this class this is made possible by the set_symmetric_mode() function. If your matrix is symmetric, you can use this class as follows:
Definition at line 899 of file petsc_solver.h.
PETScWrappers::SparseDirectMUMPS::SparseDirectMUMPS | ( | SolverControl & | cn, |
const AdditionalData & | data = AdditionalData() |
||
) |
Constructor.
Definition at line 654 of file petsc_solver.cc.
PETScWrappers::SparseDirectMUMPS::SparseDirectMUMPS | ( | SolverControl & | cn, |
const MPI_Comm & | mpi_communicator, | ||
const AdditionalData & | data = AdditionalData() |
||
) |
Constructor. This constructor is deprecated and ignores the MPI communicator argument. Use the other constructor instead.
Definition at line 663 of file petsc_solver.cc.
void PETScWrappers::SparseDirectMUMPS::solve | ( | const MatrixBase & | A, |
VectorBase & | x, | ||
const VectorBase & | b | ||
) |
The method to solve the linear system.
Definition at line 707 of file petsc_solver.cc.
void PETScWrappers::SparseDirectMUMPS::set_symmetric_mode | ( | const bool | flag | ) |
The method allows to take advantage if the system matrix is symmetric by using LDL^T decomposition instead of more expensive LU. The argument indicates whether the matrix is symmetric or not.
Definition at line 920 of file petsc_solver.cc.
|
overrideprotectedvirtual |
Function that takes a Krylov Subspace Solver context object, and sets the type of solver that is requested by the derived class.
Implements PETScWrappers::SolverBase.
Definition at line 680 of file petsc_solver.cc.
|
staticprivate |
A function that is used in PETSc as a callback to check convergence. It takes the information provided from PETSc and checks it against deal.II's own SolverControl objects to see if convergence has been reached.
Definition at line 881 of file petsc_solver.cc.
|
inherited |
Solve the linear system Ax=b
. Depending on the information provided by derived classes and the object passed as a preconditioner, one of the linear solvers and preconditioners of PETSc is chosen. Repeated calls to solve() do not reconstruct the preconditioner for performance reasons. See class Documentation.
Definition at line 52 of file petsc_solver.cc.
|
virtualinherited |
Resets the contained preconditioner and solver object. See class description for more details.
Definition at line 149 of file petsc_solver.cc.
|
inherited |
Sets a prefix name for the solver object. Useful when customizing the PETSc KSP object with command-line options.
Definition at line 142 of file petsc_solver.cc.
|
inherited |
Access to object that controls convergence.
Definition at line 156 of file petsc_solver.cc.
|
inherited |
initialize the solver with the preconditioner. This function is intended for use with SLEPc spectral transformation class.
Definition at line 203 of file petsc_solver.cc.
|
protected |
Store a copy of flags for this particular solver.
Definition at line 943 of file petsc_solver.h.
|
private |
Definition at line 978 of file petsc_solver.h.
|
private |
Flag specifies whether matrix being factorized is symmetric or not. It influences the type of the used preconditioner (PCLU or PCCHOLESKY)
Definition at line 984 of file petsc_solver.h.
|
protectedinherited |
Reference to the object that controls convergence of the iterative solver. In fact, for these PETSc wrappers, PETSc does so itself, but we copy the data from this object before starting the solution process, and copy the data back into it afterwards.
Definition at line 151 of file petsc_solver.h.
|
protectedinherited |
Solver prefix name to qualify options specific to the PETSc KSP object in the current context. Note: A hyphen (-) must NOT be given at the beginning of the prefix name. The first character of all runtime options is AUTOMATICALLY the hyphen.
Definition at line 166 of file petsc_solver.h.