1#ifndef DUNE_FEM_PETSCCOMMON_HH
2#define DUNE_FEM_PETSCCOMMON_HH
13#include <dune/common/parallel/communication.hh>
14#include <dune/common/stdstreams.hh>
15#include <dune/common/exceptions.hh>
23#define PETSC_HAVE_BROKEN_RECURSIVE_MACRO
42 class Exception :
public Dune::Exception {};
47 typedef ::PetscErrorCode ErrorCode;
55 inline static ErrorCode ErrorCheckHelper ( ErrorCode errorCode ) { CHKERRQ( errorCode );
return 0; }
57 inline ErrorCode ErrorHandler ( MPI_Comm comm,
int line,
const char *function,
const char *file,
58#
if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 5
61 ErrorCode errorCode, PetscErrorType p,
const char *message,
void *context )
63 std::ostringstream msgout;
64 msgout <<
"PETSc Error in the PETSc function '" << function <<
"' at " << file <<
":" << line <<
":";
67 PetscErrorMessage( errorCode, &text, PETSC_NULL );
69 msgout <<
" '" << text <<
"'. Error message: '" << message <<
"'";
71 msgout <<
" Unable to retrieve error text";
73 DUNE_THROW( Exception, msgout.str() );
78 inline static void ErrorCheck ( ErrorCode errorCode )
82 DUNE_THROW( Exception,
"Generic PETSC exception" );
90 inline void initialize(
const bool verbose,
int &argc,
char **&args,
const char file[] = 0 ,
const char help[] = 0,
bool ownHandler =
true )
92 ::PetscInitialize( &argc, &args, file, help );
99 dverb <<
"INFORMATION: Setting up an own error handler for PETSc errors. If you want the default PETSc handler,\n"
100 <<
"INFORMATION: set the last argument of Dune::Petsc::initialize(...) to 'false'.\n";
102 ::PetscPushErrorHandler( &ErrorHandler, 0 );
109 inline void finalize ()
112 ::PetscPopErrorHandler();
117 template <
class Comm>
118 inline auto castToPetscComm(
const Comm& comm )
122 if constexpr ( std::is_same< Dune::No_Comm, Comm > :: value ||
123 std::is_same< Dune::Communication< No_Comm >, Comm >::value )
125 return PETSC_COMM_SELF;
139 template <
class Comm>
140 inline static void KSPCreate (
const Comm& comm, KSP *inksp )
142 ErrorCheck( ::KSPCreate( castToPetscComm( comm ), inksp ) );
145 inline static void KSPDestroy ( KSP *ksp ) {
146#if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
147 ErrorCheck( ::KSPDestroy( *ksp ) );
149 ErrorCheck( ::KSPDestroy( ksp ) );
152 inline static void KSPGetPC ( KSP ksp, PC *pc ) { ErrorCheck( ::KSPGetPC( ksp, pc ) ); }
153 inline static void KSPSetFromOptions ( KSP ksp ) { ErrorCheck( ::KSPSetFromOptions( ksp ) ); }
154 inline static void KSPSetUp ( KSP ksp ) { ErrorCheck( ::KSPSetUp( ksp ) ); }
155 inline static void KSPSetType ( KSP ksp,
const KSPType type ) { ErrorCheck( ::KSPSetType( ksp, type ) ); }
156 inline static void KSPGMRESSetRestart ( KSP ksp, PetscInt restart ) { ErrorCheck( ::KSPGMRESSetRestart( ksp, restart ) ); }
158 template <
class Comm>
159 inline static void KSPView (
const Comm& comm,
163 ErrorCheck( ::KSPView( ksp, viewer ) );
166 template <
class Comm>
167 inline static void KSPView (
const Comm& comm,
170 KSPView( comm, ksp, PETSC_VIEWER_STDOUT_( castToPetscComm(comm) ) );
173 inline static void KSPMonitorSet (KSP ksp, PetscErrorCode (*monitor)(KSP,PetscInt,PetscReal,
void*),
174#
if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
175 void *mctx,PetscErrorCode (*monitordestroy)(
void*)
177 void *mctx,PetscErrorCode (*monitordestroy)(
void**)
181 ErrorCheck( ::KSPMonitorSet( ksp, monitor, mctx, monitordestroy ) );
183 inline static void KSPGetIterationNumber( KSP ksp, PetscInt* its )
184 { ErrorCheck( ::KSPGetIterationNumber( ksp, its ) ); }
185 inline static void KSPGetConvergedReason(KSP ksp, KSPConvergedReason *reason)
186 { ErrorCheck( ::KSPGetConvergedReason( ksp, reason ) ); }
189#if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 5
190 inline static void KSPSetOperators (KSP ksp, Mat Amat, Mat Pmat, MatStructure flag ) { ErrorCheck( ::KSPSetOperators( ksp, Amat, Pmat, flag ) ); }
192 inline static void KSPSetOperators (KSP ksp, Mat Amat, Mat Pmat ) { ErrorCheck( ::KSPSetOperators( ksp, Amat, Pmat ) ); }
194 inline static void KSPSetTolerances ( KSP ksp, PetscReal rtol, PetscReal abstol, PetscReal dtol, PetscInt maxits )
195 { ErrorCheck( ::KSPSetTolerances( ksp, rtol, abstol, dtol, maxits ) ); }
196 inline static void KSPSetInitialGuessNonzero( KSP ksp, PetscBool flg ) { ErrorCheck( ::KSPSetInitialGuessNonzero( ksp, flg ) ); };
197 inline static void KSPSolve ( KSP ksp, Vec b, Vec x ) { ErrorCheck( ::KSPSolve( ksp, b, x ) ); }
198 inline static void KSPSetPC ( KSP ksp, PC pc ) { ErrorCheck( ::KSPSetPC( ksp, pc ) ); }
201 template <
class Comm>
202 inline static void PCCreate (
const Comm& comm, PC* pc)
204 ErrorCheck( ::PCCreate( castToPetscComm( comm ), pc ) );
207 inline static void PCDestroy ( PC* pc) {
208#if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
209 ErrorCheck( ::PCDestroy( *pc ) );
211 ErrorCheck( ::PCDestroy( pc ) );
214 inline static void PCSetType ( PC pc,
const PCType type ) { ErrorCheck( ::PCSetType( pc, type ) ); }
215 inline static void PCSetFromOptions ( PC pc ) { ErrorCheck( ::PCSetFromOptions( pc ) ); }
216 inline static void PCSetUp ( PC pc ) { ErrorCheck( ::PCSetUp( pc ) ); }
217 inline static void PCFactorSetLevels( PC pc, PetscInt level ) { ErrorCheck( ::PCFactorSetLevels( pc, level ) ); }
218 inline static void PCSORSetOmega( PC pc, PetscReal omega ) { ErrorCheck( ::PCSORSetOmega( pc, omega ) ); }
219#if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 9
220 inline static void PCFactorSetMatSolverPackage( PC pc,
const MatSolverPackage type )
222 ErrorCheck( ::PCFactorSetMatSolverPackage( pc, type ) );
225 inline static void PCFactorSetMatSolverPackage( PC pc,
const MatSolverType type )
227 ErrorCheck( ::PCFactorSetMatSolverType( pc, type ) );
230 inline static void PCHYPRESetType( PC pc,
const char* type )
232 ErrorCheck( ::PCHYPRESetType( pc, type ) );
236 inline static void MatAssemblyBegin ( Mat mat, MatAssemblyType type ) { ErrorCheck( ::MatAssemblyBegin( mat, type ) ); }
237 inline static void MatAssemblyEnd ( Mat mat, MatAssemblyType type ) { ErrorCheck( ::MatAssemblyEnd( mat, type ) ); }
238 inline static void MatAssembled( Mat mat, PetscBool* assembled ) { ErrorCheck( ::MatAssembled( mat, assembled ) ); }
240 template <
class Comm>
241 inline static void MatCreate (
const Comm& comm, Mat *A )
243 ErrorCheck( ::MatCreate( castToPetscComm( comm ), A) );
246 template <
class Comm>
247 inline static void MatCreateBlockMat (
const Comm& comm,
249 PetscInt m, PetscInt n, PetscInt bs, PetscInt nz, PetscInt* nnz )
251 ErrorCheck( ::MatCreateBlockMat( castToPetscComm( comm ), n, m, bs, nz, nnz, A) );
253 inline static void MatDestroy ( Mat *A ) {
254 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
255 ErrorCheck( ::MatDestroy( *A ) );
257 ErrorCheck( ::MatDestroy( A ) );
260 inline static void MatSetUp( Mat mat )
262 ErrorCheck( ::MatSetUp(mat));
264 inline static void MatSetUp( Mat mat, PetscInt bs, PetscInt nz )
268 ErrorCheck( ::MatSeqAIJSetPreallocation(mat,nz,PETSC_NULL) );
269 ErrorCheck( ::MatMPIAIJSetPreallocation(mat,nz,PETSC_NULL,nz/2,PETSC_NULL) );
273 ErrorCheck( ::MatSeqBAIJSetPreallocation(mat,bs,nz,PETSC_NULL) );
274 ErrorCheck( ::MatMPIBAIJSetPreallocation(mat,bs,nz,PETSC_NULL,nz/2,PETSC_NULL) );
277 ErrorCheck( ::MatSetOption(mat, MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE) );
281 ErrorCheck( ::MatSetUp(mat));
283 inline static void MatSetUp( Mat mat, PetscInt bs,
const PetscInt *d_nnz,
const PetscInt *o_nnz )
287 ErrorCheck( ::MatSeqAIJSetPreallocation(mat,0,d_nnz ) );
288 ErrorCheck( ::MatMPIAIJSetPreallocation(mat,0,d_nnz,5,o_nnz) );
292 ErrorCheck( ::MatSeqBAIJSetPreallocation(mat,bs,0,d_nnz ) );
293 ErrorCheck( ::MatMPIBAIJSetPreallocation(mat,bs,0,d_nnz,5,PETSC_NULL) );
296 ErrorCheck( ::MatSetOption(mat, MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE) );
297 ErrorCheck( ::MatSetUp(mat));
300 inline static void MatGetOwnershipRange ( Mat mat, PetscInt *m, PetscInt* n ) { ErrorCheck( ::MatGetOwnershipRange( mat, m, n ) ); }
301 inline static void MatGetSize ( Mat mat, PetscInt *m, PetscInt* n ) { ErrorCheck( ::MatGetSize( mat, m, n ) ); }
302 inline static void MatMult ( Mat mat, Vec x, Vec y ) { ErrorCheck( ::MatMult( mat, x, y ) ); }
303 inline static void MatSetBlockSize ( Mat A, PetscInt bs ) { ErrorCheck( ::MatSetBlockSize( A, bs ) ); }
304 inline static void MatSetSizes ( Mat A, PetscInt m, PetscInt n, PetscInt M, PetscInt N ) { ErrorCheck( ::MatSetSizes( A, m, n, M, N ) ); }
305 inline static void MatSetFromOptions ( Mat B ) { ErrorCheck( ::MatSetFromOptions( B ) ); }
306 inline static void MatSetType ( Mat mat,
const MatType matype ) { ErrorCheck( ::MatSetType( mat, matype ) ); }
308 inline static void MatSetValue ( Mat v, PetscInt i, PetscInt j, PetscScalar va, InsertMode mode )
314 ErrorCheck( ::MatSetValue( v, i, j, va, mode ) );
318 inline static void MatSetValues ( Mat mat, PetscInt m,
const PetscInt idxm[], PetscInt n,
const PetscInt idxn[],
const PetscScalar v[], InsertMode addv )
324 ErrorCheck( ::MatSetValues( mat, m, idxm, n, idxn, v, addv ) );
328 inline static void MatSetValuesBlocked ( Mat mat, PetscInt m,
const PetscInt idxm[], PetscInt n,
const PetscInt idxn[],
const PetscScalar v[], InsertMode addv )
334 ErrorCheck( ::MatSetValuesBlocked( mat, m, idxm, n, idxn, v, addv ) );
338 inline static void MatGetValues ( Mat mat, PetscInt m,
const PetscInt idxm[], PetscInt n,
const PetscInt idxn[], PetscScalar v[] )
339 { ErrorCheck( ::MatGetValues( mat, m, idxm, n, idxn, v ) ); }
341 inline static void MatZeroRows ( Mat mat, PetscInt m,
const PetscInt idxm[],
const PetscScalar v )
347 ErrorCheck( ::MatZeroRows( mat, m, idxm, v, 0, 0 ) );
351 inline static void MatView ( Mat mat, PetscViewer viewer ) { ErrorCheck( ::MatView( mat, viewer ) ); }
352 inline static void MatZeroEntries ( Mat mat )
358 ErrorCheck( ::MatZeroEntries( mat ) );
362 inline static void PetscBarrier ( PetscObject obj ) { ErrorCheck( ::PetscBarrier( obj ) ); }
363 inline static void PetscFinalize () { ErrorCheck( ::PetscFinalize() ); }
364 inline static void PetscInitialize(
int *argc,
char ***args,
const char file[],
const char help[] ) { ErrorCheck( ::PetscInitialize( argc, args, file, help ) ); }
365 inline static void PetscViewerASCIIOpen ( MPI_Comm comm,
const char name[], PetscViewer *lab ) { ErrorCheck( ::PetscViewerASCIIOpen( comm, name, lab ) ); }
366 inline static void PetscViewerDestroy ( PetscViewer *viewer ) {
367 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
368 ErrorCheck( ::PetscViewerDestroy( *viewer ) );
370 ErrorCheck( ::PetscViewerDestroy( viewer ) );
373 inline static void PetscViewerSetFormat ( PetscViewer viewer, PetscViewerFormat format )
375 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 5
376 ErrorCheck( ::PetscViewerSetFormat( viewer, format ) );
378 ErrorCheck( ::PetscViewerPushFormat( viewer, format ) );
381 inline static void VecAssemblyBegin ( Vec vec ) { ErrorCheck( ::VecAssemblyBegin( vec ) ); }
382 inline static void VecAssemblyEnd ( Vec vec ) { ErrorCheck( ::VecAssemblyEnd( vec ) ); }
383 inline static void VecAXPY ( Vec y, PetscScalar alpha, Vec x) { ErrorCheck( ::VecAXPY( y, alpha, x ) ); }
384 inline static void VecCopy ( Vec x, Vec y ) { ErrorCheck( ::VecCopy( x, y ) ); }
386 template <
class Comm>
387 inline static void VecCreate (
const Comm& comm, Vec *vec )
389 ErrorCheck( ::VecCreate( castToPetscComm( comm ), vec ) );
392 template <
class Comm>
393 inline static void VecCreateGhost (
const Comm& comm, PetscInt n, PetscInt N, PetscInt nghost,
const PetscInt ghosts[], Vec *vv )
394 { ErrorCheck( ::VecCreateGhost( castToPetscComm( comm ), n, N, nghost, ghosts, vv ) ); }
396 template <
class Comm>
397 inline static void VecCreateGhostBlock (
const Comm& comm, PetscInt bs, PetscInt n, PetscInt N, PetscInt nghost,
const PetscInt ghosts[], Vec *vv )
399#if PETSC_VERSION_MAJOR < 3 || (PETSC_VERSION_MAJOR == 3 && PETSC_VERSION_MINOR <= 2)
400 std::unique_ptr< PetscInt[] > ghostsCopy(
new PetscInt[ nghost * bs ] );
401 for( PetscInt i = 0; i < nghost; ++i )
403 for( PetscInt j = 0; j < bs; ++j )
404 ghostsCopy[ i*bs + j ] = ghosts[ i ]*bs + j;
406 VecCreateGhost( n, N, nghost * bs, ghostsCopy.get(), vv );
408 ErrorCheck( ::VecCreateGhostBlock( castToPetscComm( comm ), bs, n, N, nghost, ghosts, vv ) );
412 inline static void VecDestroy ( Vec *v ) {
413 #if PETSC_VERSION_MAJOR <= 3 && PETSC_VERSION_MINOR < 2
414 ErrorCheck( ::VecDestroy( *v ) );
416 ErrorCheck( ::VecDestroy( v ) );
419 inline static void VecDot ( Vec x, Vec y, PetscScalar *val ) { ErrorCheck( ::VecDot( x, y, val ) ); }
420 inline static void VecDuplicate ( Vec v, Vec *newv ) { ErrorCheck( ::VecDuplicate( v, newv ) ); }
421 inline static void VecGetBlockSize ( Vec x, PetscInt* bs ) { ErrorCheck( ::VecGetBlockSize( x, bs ) ); }
422 inline static void VecGetLocalSize ( Vec x, PetscInt *size ) { ErrorCheck( ::VecGetLocalSize( x, size ) ); }
423 inline static void VecGetOwnershipRange ( Vec x, PetscInt *low, PetscInt *high ) { ErrorCheck( ::VecGetOwnershipRange( x, low, high ) ); }
424 inline static void VecGetSize ( Vec x, PetscInt *size ) { ErrorCheck( ::VecGetSize( x, size ) ); }
425 inline static void VecGetValues ( Vec x, PetscInt ni,
const PetscInt ix[], PetscScalar y[] ) { ErrorCheck( ::VecGetValues( x, ni, ix, y ) ); }
426 inline static void VecGhostGetLocalForm ( Vec g, Vec *l ) { ErrorCheck( ::VecGhostGetLocalForm( g, l ) ); }
427 inline static void VecGhostRestoreLocalForm ( Vec g, Vec *l ) { ErrorCheck( ::VecGhostRestoreLocalForm( g, l ) ); }
428 inline static void VecGhostUpdateBegin ( Vec g, InsertMode insertmode, ScatterMode scattermode ) { ErrorCheck( ::VecGhostUpdateBegin( g, insertmode, scattermode ) ); }
429 inline static void VecGhostUpdateEnd ( Vec g, InsertMode insertmode, ScatterMode scattermode ) { ErrorCheck( ::VecGhostUpdateEnd( g, insertmode, scattermode ) ); }
430 inline static void VecNorm ( Vec x, NormType type, PetscReal *val ) { ErrorCheck( ::VecNorm( x, type, val ) ); }
431 inline static void VecScale ( Vec x, PetscScalar alpha ) { ErrorCheck( ::VecScale( x, alpha ) ); }
432 inline static void VecSet ( Vec x, PetscScalar alpha ) { ErrorCheck( ::VecSet( x, alpha ) ); }
433 inline static void VecSetBlockSize ( Vec x, PetscInt bs ) { ErrorCheck( ::VecSetBlockSize( x, bs ) ); }
434 inline static void VecSetFromOptions ( Vec vec ) { ErrorCheck( ::VecSetFromOptions( vec ) ); }
435 inline static void VecSetType ( Vec vec,
const VecType method ) { ErrorCheck( ::VecSetType( vec, method ) ); }
436 inline static void VecSetSizes ( Vec v, PetscInt n, PetscInt N ) { ErrorCheck( ::VecSetSizes( v, n, N ) ); }
437 inline static void VecSetValue ( Vec v,
int row, PetscScalar value, InsertMode mode ) { ErrorCheck( ::VecSetValue( v, row, value, mode ) ); }
438 inline static void VecSetValuesBlocked ( Vec v, PetscInt ni,
const PetscInt xi[],
const PetscScalar values[], InsertMode mode ) { ErrorCheck( ::VecSetValuesBlocked( v, ni, xi, values, mode ) ); }
439 inline static void VecView ( Vec vec, PetscViewer viewer ) { ErrorCheck( ::VecView( vec, viewer ) ); }
Definition: bindguard.hh:11