AMDiS  0.3
The Adaptive Multi-Dimensional Simulation Toolbox
Interface.hpp
1 #pragma once
2 
3 #include <petscksp.h>
4 #include <petsc/private/kspimpl.h>
5 
6 #include <amdis/Output.hpp>
7 
8 namespace AMDiS
9 {
10  namespace PETSc
11  {
12  // KSP routines
13 
15  inline PetscErrorCode KSPMonitorCyclic(KSP ksp, PetscInt n, PetscReal rnorm, void* ctx)
16  {
17  if (n % 10 == 0)
18  msg("iteration {:>4}: resid {:<.6e}", n, rnorm);
19  return 0;
20  }
21 
23  inline PetscErrorCode KSPMonitorNoisy(KSP ksp, PetscInt n, PetscReal rnorm, void* ctx)
24  {
25  Vec resid;
26  KSPBuildResidual(ksp, NULL, NULL, &resid);
27 
28  PetscReal truenorm = 0.0;
29  VecNorm(resid, NORM_2, &truenorm);
30  VecDestroy(&resid);
31 
32  PetscReal bnorm = 0.0;
33  VecNorm(ksp->vec_rhs, NORM_2, &bnorm);
34 
35  msg("iteration {:>4}: resid {:<.6e}, true resid {:<.2e}, rel. resid {:<.2e}", n, rnorm, truenorm, truenorm/bnorm);
36  return 0;
37  }
38 
40  inline PetscErrorCode KSPMonitorTrueResidualNorm(KSP ksp, PetscInt n, PetscReal rnorm, void* ctx)
41  {
42  assert(ctx != nullptr);
43 #if (PETSC_VERSION_MINOR >= 7)
44  return ::KSPMonitorTrueResidualNorm(ksp, n, rnorm, (PetscViewerAndFormat*)(ctx));
45 #else
46  return ::KSPMonitorTrueResidualNorm(ksp, n, rnorm, ctx);
47 #endif
48  }
49 
51  inline PetscErrorCode KSPMonitorDefault(KSP ksp, PetscInt n, PetscReal rnorm, void* ctx)
52  {
53  assert(ctx != nullptr);
54 #if (PETSC_VERSION_MINOR >= 7)
55  return ::KSPMonitorDefault(ksp, n, rnorm, (PetscViewerAndFormat*)(ctx));
56 #else
57  return ::KSPMonitorDefault(ksp, n, rnorm, ctx);
58 #endif
59  }
60 
62  template <class Monitor>
63  inline PetscErrorCode KSPMonitorSet(KSP ksp, Monitor monitor)
64  {
65 #if (PETSC_VERSION_MINOR >= 7)
66  PetscViewerAndFormat *vf;
67  PetscErrorCode ierr;
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);
70  return ierr;
71 #else
72  return ::KSPMonitorSet(ksp, monitor, PETSC_NULL, PETSC_NULL);
73 #endif
74  }
75 
77  inline PetscErrorCode KSPGetOperators(KSP ksp, Mat *Amat, Mat *Pmat)
78  {
79 #if (PETSC_VERSION_MINOR >= 5)
80  return ::KSPGetOperators(ksp, Amat, Pmat);
81 #else
82  return ::KSPGetOperators(ksp, Amat, Pmat, PETSC_NULL);
83 #endif
84  }
85 
87  inline PetscErrorCode KSPSetOperators(KSP ksp, Mat Amat, Mat Pmat)
88  {
89 #if (PETSC_VERSION_MINOR >= 5)
90  return ::KSPSetOperators(ksp, Amat, Pmat);
91 #else
92  return ::KSPSetOperators(ksp, Amat, Pmat, SAME_NONZERO_PATTERN);
93 #endif
94  }
95 
96  // Mat routines
97 
99  inline PetscErrorCode MatCreateVecs(Mat mat, Vec *right, Vec *left)
100  {
101 #if (PETSC_VERSION_MINOR >= 6)
102  return ::MatCreateVecs(mat, right, left);
103 #else
104  return ::MatGetVecs(mat, right, left);
105 #endif
106  }
107 
109  inline PetscErrorCode MatCreateSubMatrix(Mat mat, IS isrow, IS iscol, MatReuse cll, Mat *newmat)
110  {
111 #if (PETSC_VERSION_MINOR >= 8)
112  return ::MatCreateSubMatrix(mat, isrow, iscol, cll, newmat);
113 #else
114  return ::MatGetSubMatrix(mat, isrow, iscol, cll, newmat);
115 #endif
116  }
117 
119  inline PetscErrorCode MatNullSpaceRemove(MatNullSpace sp, Vec vec)
120  {
121 #if (PETSC_VERSION_MINOR >= 5)
122  return ::MatNullSpaceRemove(sp, vec);
123 #else
124  return ::MatNullSpaceRemove(sp, vec, PETSC_NULL);
125 #endif
126  }
127 
128  // PC routines
129 
130 #if (PETSC_VERSION_MINOR >= 9)
131  inline PetscErrorCode PCFactorSetMatSolverType(PC pc, MatSolverType stype)
133  {
134  return ::PCFactorSetMatSolverType(pc, stype);
135  }
136 #else
137  template <class MatSolverType>
139  inline PetscErrorCode PCFactorSetMatSolverType(PC pc, MatSolverType stype)
140  {
141  return ::PCFactorSetMatSolverPackage(pc, stype);
142  }
143 #endif
144 
146  inline PetscErrorCode PetscOptionsView(PetscViewer viewer)
147  {
148 #if (PETSC_VERSION_MINOR >= 7)
149  return ::PetscOptionsView(PETSC_NULL, viewer);
150 #else
151  return ::PetscOptionsView(viewer);
152 #endif
153  }
154 
155  // PETSc routine
156 
158  inline PetscErrorCode PetscOptionsInsertString(const char in_str[])
159  {
160 #if (PETSC_VERSION_MINOR >= 7)
161  return ::PetscOptionsInsertString(PETSC_NULL, in_str);
162 #else
163  return ::PetscOptionsInsertString(in_str);
164 #endif
165  }
166 
167  } // end namespace PETSc
168 } // end namespace AMDiS
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