2 #ifndef DUNE_FEM_PETSCCOMMON_HH 3 #define DUNE_FEM_PETSCCOMMON_HH 14 #include <dune/common/stdstreams.hh> 15 #include <dune/common/exceptions.hh> 23 #define PETSC_HAVE_BROKEN_RECURSIVE_MACRO 41 #define FEM_PETSC_COMM_DEFAULT PETSC_COMM_WORLD 43 #define FEM_PETSC_COMM_DEFAULT PETSC_COMM_SELF 49 class Exception :
public Dune::Exception {};
54 typedef ::PetscErrorCode ErrorCode;
62 inline ErrorCode ErrorCheckHelper ( ErrorCode errorCode ) { CHKERRQ( errorCode );
return 0; }
64 inline ErrorCode ErrorHandler (
MPI_Comm comm,
int line,
const char *
function,
const char *file,
65 #
if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 5
68 ErrorCode errorCode, PetscErrorType p,
const char *message,
void *context )
70 std::ostringstream msgout;
71 msgout <<
"PETSc Error in the PETSc function '" <<
function <<
"' at " << file <<
":" << line <<
":";
74 PetscErrorMessage( errorCode, &text, PETSC_NULL );
76 msgout <<
" '" << text <<
"'. Error message: '" << message <<
"'";
78 msgout <<
" Unable to retrieve error text";
80 DUNE_THROW( Exception, msgout.str() );
85 inline void ErrorCheck ( ErrorCode errorCode )
89 DUNE_THROW( Exception,
"Generic PETSC exception" );
97 inline void initialize(
const bool verbose,
int &argc,
char **&args,
const char file[] = 0 ,
const char help[] = 0,
bool ownHandler =
true )
99 ::PetscInitialize( &argc, &args, file, help );
106 dverb <<
"INFORMATION: Setting up an own error handler for PETSc errors. If you want the default PETSc handler,\n" 107 <<
"INFORMATION: set the last argument of Dune::Petsc::initialize(...) to 'false'.\n";
109 ::PetscPushErrorHandler( &ErrorHandler, 0 );
118 inline void initialize(
int *argc,
char ***args,
const char file[],
const char help[],
bool ownHandler =
true )
120 ::PetscInitialize( argc, args, file, help );
127 std::cerr <<
"INFORMATION: Setting up an own error handler for PETSc errors. If you want the default PETSc handler,\n" 128 <<
"INFORMATION: set the last argument of Dune::Petsc::initialize(...) to 'false'.\n";
130 ::PetscPushErrorHandler( &ErrorHandler, 0 );
138 inline void finalize ()
141 ::PetscPopErrorHandler();
152 inline void KSPCreate ( KSP *inksp ) { ErrorCheck( ::KSPCreate( FEM_PETSC_COMM_DEFAULT, inksp ) ); }
153 inline void KSPDestroy ( KSP *ksp ) {
154 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 155 ErrorCheck( ::KSPDestroy( *ksp ) );
157 ErrorCheck( ::KSPDestroy( ksp ) );
158 #endif // PETSC_PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 160 inline void KSPGetPC ( KSP ksp, PC *pc ) { ErrorCheck( ::KSPGetPC( ksp, pc ) ); }
161 inline void KSPSetFromOptions ( KSP ksp ) { ErrorCheck( ::KSPSetFromOptions( ksp ) ); }
162 inline void KSPSetType ( KSP ksp,
const KSPType type ) { ErrorCheck( ::KSPSetType( ksp, type ) ); }
163 inline void KSPGMRESSetRestart ( KSP ksp, PetscInt restart ) { ErrorCheck( ::KSPGMRESSetRestart( ksp, restart ) ); }
164 inline void KSPView ( KSP ksp, PetscViewer viewer = PETSC_VIEWER_STDOUT_(FEM_PETSC_COMM_DEFAULT) )
166 ErrorCheck( ::KSPView( ksp, viewer ) );
168 inline void KSPMonitorSet (KSP ksp, PetscErrorCode (*monitor)(KSP,PetscInt,PetscReal,
void*),
169 #
if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
170 void *mctx,PetscErrorCode (*monitordestroy)(
void*)
172 void *mctx,PetscErrorCode (*monitordestroy)(
void**)
176 ErrorCheck( ::KSPMonitorSet( ksp, monitor, mctx, monitordestroy ) );
178 inline void KSPGetIterationNumber( KSP ksp, PetscInt* its )
179 { ErrorCheck( ::KSPGetIterationNumber( ksp, its ) ); }
182 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 5 183 inline void KSPSetOperators (KSP ksp, Mat Amat, Mat Pmat, MatStructure flag ) { ErrorCheck( ::KSPSetOperators( ksp, Amat, Pmat, flag ) ); }
185 inline void KSPSetOperators (KSP ksp, Mat Amat, Mat Pmat ) { ErrorCheck( ::KSPSetOperators( ksp, Amat, Pmat ) ); }
187 inline void KSPSetTolerances ( KSP ksp, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt maxits )
188 { ErrorCheck( ::KSPSetTolerances( ksp, rtol, abstol, dtol, maxits ) ); }
189 inline void KSPSolve ( KSP ksp, Vec b, Vec x ) { ErrorCheck( ::KSPSolve( ksp, b, x ) ); }
190 inline void KSPSetPC ( KSP ksp, PC pc ) { ErrorCheck( ::KSPSetPC( ksp, pc ) ); }
193 inline void PCCreate ( PC* pc) { ErrorCheck( ::PCCreate( FEM_PETSC_COMM_DEFAULT, pc ) ); }
194 inline void PCDestroy ( PC* pc) {
195 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 196 ErrorCheck( ::PCDestroy( *pc ) );
198 ErrorCheck( ::PCDestroy( pc ) );
201 inline void PCSetType ( PC pc,
const PCType type ) { ErrorCheck( ::PCSetType( pc, type ) ); }
202 inline void PCFactorSetLevels( PC pc, PetscInt level ) { ErrorCheck( ::PCFactorSetLevels( pc, level ) ); }
203 inline void PCFactorSetMatSolverPackage( PC pc,
const MatSolverPackage type )
205 ErrorCheck( ::PCFactorSetMatSolverPackage( pc, type ) );
207 inline void PCHYPRESetType( PC pc,
const char* type )
210 ErrorCheck( ::PCHYPRESetType( pc, type ) );
212 DUNE_THROW( InvalidStateException,
"HYPRE not found within PETSc, please re-configure!");
217 inline void MatAssemblyBegin ( Mat mat, MatAssemblyType type ) { ErrorCheck( ::MatAssemblyBegin( mat, type ) ); }
218 inline void MatAssemblyEnd ( Mat mat, MatAssemblyType type ) { ErrorCheck( ::MatAssemblyEnd( mat, type ) ); }
219 inline void MatCreate ( Mat *A ) { ErrorCheck( ::MatCreate( FEM_PETSC_COMM_DEFAULT, A) ); }
220 inline void MatCreateBlockMat ( Mat *A, PetscInt m, PetscInt n, PetscInt bs, PetscInt nz, PetscInt* nnz )
222 ErrorCheck( ::MatCreateBlockMat( FEM_PETSC_COMM_DEFAULT, n, m, bs, nz, nnz, A) );
224 inline void MatDestroy ( Mat *A ) {
225 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 226 ErrorCheck( ::MatDestroy( *A ) );
228 ErrorCheck( ::MatDestroy( A ) );
229 #endif // PETSC_PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 231 inline void MatSetUp( Mat mat )
233 ErrorCheck( ::MatSetUp(mat));
235 inline void MatSetUp( Mat mat, PetscInt bs,
int nz )
239 ErrorCheck( ::MatSeqAIJSetPreallocation(mat,nz,PETSC_NULL) );
240 ErrorCheck( ::MatMPIAIJSetPreallocation(mat,nz,PETSC_NULL,nz/2,PETSC_NULL) );
244 ErrorCheck( ::MatSeqBAIJSetPreallocation(mat,bs,nz,PETSC_NULL) );
245 ErrorCheck( ::MatMPIBAIJSetPreallocation(mat,bs,nz,PETSC_NULL,nz/2,PETSC_NULL) );
248 ErrorCheck( ::MatSetOption(mat, MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE) );
252 ErrorCheck( ::MatSetUp(mat));
254 inline void MatSetUp( Mat mat, PetscInt bs,
const int *d_nnz,
const int *o_nnz )
258 ErrorCheck( ::MatSeqAIJSetPreallocation(mat,0,d_nnz ) );
259 ErrorCheck( ::MatMPIAIJSetPreallocation(mat,0,d_nnz,5,o_nnz) );
263 ErrorCheck( ::MatSeqBAIJSetPreallocation(mat,bs,0,d_nnz ) );
264 ErrorCheck( ::MatMPIBAIJSetPreallocation(mat,bs,0,d_nnz,5,PETSC_NULL) );
267 ErrorCheck( ::MatSetOption(mat, MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE) );
268 ErrorCheck( ::MatSetUp(mat));
271 inline void MatGetOwnershipRange ( Mat mat, PetscInt *m, PetscInt* n ) { ErrorCheck( ::MatGetOwnershipRange( mat, m, n ) ); }
272 inline void MatGetSize ( Mat mat, PetscInt *m, PetscInt* n ) { ErrorCheck( ::MatGetSize( mat, m, n ) ); }
273 inline void MatMult ( Mat mat, Vec x, Vec y ) { ErrorCheck( ::MatMult( mat, x, y ) ); }
274 inline void MatSetBlockSize ( Mat A, PetscInt bs ) { ErrorCheck( ::MatSetBlockSize( A, bs ) ); }
275 inline void MatSetSizes ( Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N ) { ErrorCheck( ::MatSetSizes( A, m, n, M, N ) ); }
276 inline void MatSetFromOptions ( Mat B ) { ErrorCheck( ::MatSetFromOptions( B ) ); }
277 inline void MatSetType ( Mat mat,
const MatType matype ) { ErrorCheck( ::MatSetType( mat, matype ) ); }
278 inline void MatSetValue ( Mat v, PetscInt i, PetscInt j, PetscScalar va, InsertMode mode ) { ErrorCheck( ::MatSetValue( v, i, j, va, mode ) ); }
279 inline void MatSetValues ( Mat mat, PetscInt m,
const PetscInt idxm[], PetscInt n,
const PetscInt idxn[],
const PetscScalar v[], InsertMode addv )
280 { ErrorCheck( ::MatSetValues( mat, m, idxm, n, idxn, v, addv ) ); }
281 inline void MatGetValues ( Mat mat, PetscInt m,
const PetscInt idxm[], PetscInt n,
const PetscInt idxn[], PetscScalar v[] )
282 { ErrorCheck( ::MatGetValues( mat, m, idxm, n, idxn, v ) ); }
283 inline void MatView ( Mat mat, PetscViewer viewer ) { ErrorCheck( ::MatView( mat, viewer ) ); }
284 inline void MatZeroEntries ( Mat mat ) { ErrorCheck( ::MatZeroEntries( mat ) ); }
285 inline void PetscBarrier ( PetscObject obj ) { ErrorCheck( ::PetscBarrier( obj ) ); }
286 inline void PetscFinalize () { ErrorCheck( ::PetscFinalize() ); }
287 inline void PetscInitialize(
int *argc,
char ***args,
const char file[],
const char help[] ) { ErrorCheck( ::PetscInitialize( argc, args, file, help ) ); }
288 inline void PetscViewerASCIIOpen (
MPI_Comm comm,
const char name[], PetscViewer *lab ) { ErrorCheck( ::PetscViewerASCIIOpen( comm, name, lab ) ); }
289 inline void PetscViewerDestroy ( PetscViewer *viewer ) {
290 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 291 ErrorCheck( ::PetscViewerDestroy( *viewer ) );
293 ErrorCheck( ::PetscViewerDestroy( viewer ) );
296 inline void PetscViewerSetFormat ( PetscViewer viewer, PetscViewerFormat format ) { ErrorCheck( ::PetscViewerSetFormat( viewer, format ) ); }
297 inline void VecAssemblyBegin ( Vec vec ) { ErrorCheck( ::VecAssemblyBegin( vec ) ); }
298 inline void VecAssemblyEnd ( Vec vec ) { ErrorCheck( ::VecAssemblyEnd( vec ) ); }
299 inline void VecAXPY ( Vec y, PetscScalar alpha, Vec x) { ErrorCheck( ::VecAXPY( y, alpha, x ) ); }
300 inline void VecCopy ( Vec x, Vec y ) { ErrorCheck( ::VecCopy( x, y ) ); }
301 inline void VecCreate ( Vec *vec ) { ErrorCheck( ::VecCreate( FEM_PETSC_COMM_DEFAULT, vec ) ); }
302 inline void VecCreateGhost ( PetscInt n, PetscInt N, PetscInt nghost,
const PetscInt ghosts[], Vec *vv )
303 { ErrorCheck( ::VecCreateGhost( FEM_PETSC_COMM_DEFAULT, n, N, nghost, ghosts, vv ) ); }
304 inline void VecDestroy ( Vec *v ) {
305 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 306 ErrorCheck( ::VecDestroy( *v ) );
308 ErrorCheck( ::VecDestroy( v ) );
309 #endif // PETSC_PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2 311 inline void VecDot ( Vec x, Vec y, PetscScalar *val ) { ErrorCheck( ::VecDot( x, y, val ) ); }
312 inline void VecDuplicate ( Vec v, Vec *newv ) { ErrorCheck( ::VecDuplicate( v, newv ) ); }
313 inline void VecGetLocalSize ( Vec x, PetscInt *size ) { ErrorCheck( ::VecGetLocalSize( x, size ) ); }
314 inline void VecGetOwnershipRange ( Vec x, PetscInt *low, PetscInt *high ) { ErrorCheck( ::VecGetOwnershipRange( x, low, high ) ); }
315 inline void VecGetSize ( Vec x, PetscInt *size ) { ErrorCheck( ::VecGetSize( x, size ) ); }
316 inline void VecGetValues ( Vec x, PetscInt ni,
const PetscInt ix[], PetscScalar y[] ) { ErrorCheck( ::VecGetValues( x, ni, ix, y ) ); }
317 inline void VecGhostGetLocalForm ( Vec g, Vec *l ) { ErrorCheck( ::VecGhostGetLocalForm( g, l ) ); }
318 inline void VecGhostRestoreLocalForm ( Vec g, Vec *l ) { ErrorCheck( ::VecGhostRestoreLocalForm( g, l ) ); }
319 inline void VecGhostUpdateBegin ( Vec g, InsertMode insertmode, ScatterMode scattermode ) { ErrorCheck( ::VecGhostUpdateBegin( g, insertmode, scattermode ) ); }
320 inline void VecGhostUpdateEnd ( Vec g, InsertMode insertmode, ScatterMode scattermode ) { ErrorCheck( ::VecGhostUpdateEnd( g, insertmode, scattermode ) ); }
321 inline void VecNorm ( Vec x, NormType type, PetscReal *val ) { ErrorCheck( ::VecNorm( x, type, val ) ); }
322 inline void VecScale ( Vec x, PetscScalar alpha ) { ErrorCheck( ::VecScale( x, alpha ) ); }
323 inline void VecSet ( Vec x, PetscScalar alpha ) { ErrorCheck( ::VecSet( x, alpha ) ); }
324 inline void VecSetFromOptions ( Vec vec ) { ErrorCheck( ::VecSetFromOptions( vec ) ); }
325 inline void VecSetType ( Vec vec,
const VecType method ) { ErrorCheck( ::VecSetType( vec, method ) ); }
326 inline void VecSetSizes ( Vec v, PetscInt n, PetscInt N ) { ErrorCheck( ::VecSetSizes( v, n, N ) ); }
327 inline void VecSetValue ( Vec v,
int row, PetscScalar value, InsertMode mode ) { ErrorCheck( ::VecSetValue( v, row, value, mode ) ); }
328 inline void VecView ( Vec vec, PetscViewer viewer ) { ErrorCheck( ::VecView( vec, viewer ) ); }
334 #endif // #if HAVE_PETSC 336 #endif // DUNE_FEM_PETSCCOMMON_HH static int rank()
Definition: mpimanager.hh:116
int MPI_Comm
Definition: odeobjectfiles.cc:12
Definition: coordinate.hh:4