4 #include <petsc/private/kspimpl.h> 6 #include <amdis/Output.hpp> 15 inline PetscErrorCode KSPMonitorCyclic(KSP ksp, PetscInt n, PetscReal rnorm,
void* ctx)
18 msg(
"iteration {:>4}: resid {:<.6e}", n, rnorm);
23 inline PetscErrorCode KSPMonitorNoisy(KSP ksp, PetscInt n, PetscReal rnorm,
void* ctx)
26 KSPBuildResidual(ksp, NULL, NULL, &resid);
28 PetscReal truenorm = 0.0;
29 VecNorm(resid, NORM_2, &truenorm);
32 PetscReal bnorm = 0.0;
33 VecNorm(ksp->vec_rhs, NORM_2, &bnorm);
35 msg(
"iteration {:>4}: resid {:<.6e}, true resid {:<.2e}, rel. resid {:<.2e}", n, rnorm, truenorm, truenorm/bnorm);
40 inline PetscErrorCode KSPMonitorTrueResidualNorm(KSP ksp, PetscInt n, PetscReal rnorm,
void* ctx)
42 assert(ctx !=
nullptr);
43 #if (PETSC_VERSION_MINOR >= 7) 44 return ::KSPMonitorTrueResidualNorm(ksp, n, rnorm, (PetscViewerAndFormat*)(ctx));
46 return ::KSPMonitorTrueResidualNorm(ksp, n, rnorm, ctx);
51 inline PetscErrorCode KSPMonitorDefault(KSP ksp, PetscInt n, PetscReal rnorm,
void* ctx)
53 assert(ctx !=
nullptr);
54 #if (PETSC_VERSION_MINOR >= 7) 55 return ::KSPMonitorDefault(ksp, n, rnorm, (PetscViewerAndFormat*)(ctx));
57 return ::KSPMonitorDefault(ksp, n, rnorm, ctx);
62 template <
class Monitor>
63 inline PetscErrorCode KSPMonitorSet(KSP ksp, Monitor monitor)
65 #if (PETSC_VERSION_MINOR >= 7) 66 PetscViewerAndFormat *vf;
68 ierr = ::PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_DEFAULT,&vf);CHKERRQ(ierr);
69 ierr = ::KSPMonitorSet(ksp,(PetscErrorCode (*)(KSP,PetscInt,PetscReal,
void*))monitor,vf,(PetscErrorCode (*)(
void**))PetscViewerAndFormatDestroy);CHKERRQ(ierr);
72 return ::KSPMonitorSet(ksp, monitor, PETSC_NULL, PETSC_NULL);
77 inline PetscErrorCode KSPGetOperators(KSP ksp, Mat *Amat, Mat *Pmat)
79 #if (PETSC_VERSION_MINOR >= 5) 80 return ::KSPGetOperators(ksp, Amat, Pmat);
82 return ::KSPGetOperators(ksp, Amat, Pmat, PETSC_NULL);
87 inline PetscErrorCode KSPSetOperators(KSP ksp, Mat Amat, Mat Pmat)
89 #if (PETSC_VERSION_MINOR >= 5) 90 return ::KSPSetOperators(ksp, Amat, Pmat);
92 return ::KSPSetOperators(ksp, Amat, Pmat, SAME_NONZERO_PATTERN);
99 inline PetscErrorCode MatCreateVecs(Mat mat, Vec *right, Vec *left)
101 #if (PETSC_VERSION_MINOR >= 6) 102 return ::MatCreateVecs(mat, right, left);
104 return ::MatGetVecs(mat, right, left);
109 inline PetscErrorCode MatCreateSubMatrix(Mat mat, IS isrow, IS iscol, MatReuse cll, Mat *newmat)
111 #if (PETSC_VERSION_MINOR >= 8) 112 return ::MatCreateSubMatrix(mat, isrow, iscol, cll, newmat);
114 return ::MatGetSubMatrix(mat, isrow, iscol, cll, newmat);
119 inline PetscErrorCode MatNullSpaceRemove(MatNullSpace sp, Vec vec)
121 #if (PETSC_VERSION_MINOR >= 5) 122 return ::MatNullSpaceRemove(sp, vec);
124 return ::MatNullSpaceRemove(sp, vec, PETSC_NULL);
130 #if (PETSC_VERSION_MINOR >= 9) 131 inline PetscErrorCode PCFactorSetMatSolverType(PC pc, MatSolverType stype)
134 return ::PCFactorSetMatSolverType(pc, stype);
137 template <
class MatSolverType>
139 inline PetscErrorCode PCFactorSetMatSolverType(PC pc, MatSolverType stype)
141 return ::PCFactorSetMatSolverPackage(pc, stype);
146 inline PetscErrorCode PetscOptionsView(PetscViewer viewer)
148 #if (PETSC_VERSION_MINOR >= 7) 149 return ::PetscOptionsView(PETSC_NULL, viewer);
151 return ::PetscOptionsView(viewer);
158 inline PetscErrorCode PetscOptionsInsertString(
const char in_str[])
160 #if (PETSC_VERSION_MINOR >= 7) 161 return ::PetscOptionsInsertString(PETSC_NULL, in_str);
163 return ::PetscOptionsInsertString(in_str);
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
void msg(std::string const &str, Args &&... args)
print a message
Definition: Output.hpp:98