17 #ifndef __deal2__petsc_solver_h
18 #define __deal2__petsc_solver_h
21 #include <deal.II/base/config.h>
23 #ifdef DEAL_II_WITH_PETSC
25 # include <deal.II/lac/exceptions.h>
26 # include <deal.II/lac/solver_control.h>
27 # include <deal.II/base/std_cxx1x/shared_ptr.h>
29 # include <petscksp.h>
31 DEAL_II_NAMESPACE_OPEN
34 namespace PETScWrappers
39 class PreconditionerBase;
151 virtual void reset();
174 <<
"An error with error number " << arg1
175 <<
" occurred while calling a PETSc function");
228 const PetscInt iteration,
229 const PetscReal residual_norm,
230 KSPConvergedReason *reason,
1220 const PetscInt iteration,
1221 const PetscReal residual_norm,
1222 KSPConvergedReason *reason,
1245 std_cxx1x::shared_ptr<SolverDataMUMPS> solver_data;
1258 DEAL_II_NAMESPACE_CLOSE
1260 #endif // DEAL_II_WITH_PETSC
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b, const PreconditionerBase &preconditioner)
virtual void set_solver_type(KSP &ksp) const =0
const AdditionalData additional_data
DeclException1(ExcPETScError, int,<< "An error with error number "<< arg1<< " occurred while calling a PETSc function")
const AdditionalData additional_data
const AdditionalData additional_data
SolverCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
const AdditionalData additional_data
const AdditionalData additional_data
const AdditionalData additional_data
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
SolverBicgstab(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
SparseDirectMUMPS(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
void set_prefix(const std::string &prefix)
const MPI_Comm mpi_communicator
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
bool right_preconditioning
SolverBase(SolverControl &cn, const MPI_Comm &mpi_communicator)
const AdditionalData additional_data
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
SolverCR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
SolverControl & solver_control
const AdditionalData additional_data
SolverBiCG(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverLSQR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverChebychev(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverCGS(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverRichardson(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
const AdditionalData additional_data
virtual void set_solver_type(KSP &ksp) const
virtual void set_solver_type(KSP &ksp) const
const AdditionalData additional_data
AdditionalData(const unsigned int restart_parameter=30, const bool right_preconditioning=false)
void solve(const MatrixBase &A, VectorBase &x, const VectorBase &b)
SolverTCQMR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
virtual void set_solver_type(KSP &ksp) const
const AdditionalData additional_data
SolverTFQMR(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverGMRES(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
SolverControl & control() const
void set_symmetric_mode(const bool flag)
std_cxx1x::shared_ptr< SolverData > solver_data
SolverPreOnly(SolverControl &cn, const MPI_Comm &mpi_communicator=PETSC_COMM_SELF, const AdditionalData &data=AdditionalData())
static PetscErrorCode convergence_test(KSP ksp, const PetscInt iteration, const PetscReal residual_norm, KSPConvergedReason *reason, void *solver_control)
const AdditionalData additional_data
unsigned int restart_parameter
AdditionalData(const double omega=1)
const AdditionalData additional_data