18 #ifndef LIBMESH_PETSC_MACRO_H 19 #define LIBMESH_PETSC_MACRO_H 22 #include "libmesh/libmesh_config.h" 24 #ifdef LIBMESH_HAVE_PETSC 33 #define PETSC_VERSION_LESS_THAN(major,minor,subminor) \ 34 ((LIBMESH_DETECTED_PETSC_VERSION_MAJOR < (major) || \ 35 (LIBMESH_DETECTED_PETSC_VERSION_MAJOR == (major) && (LIBMESH_DETECTED_PETSC_VERSION_MINOR < (minor) || \ 36 (LIBMESH_DETECTED_PETSC_VERSION_MINOR == (minor) && \ 37 LIBMESH_DETECTED_PETSC_VERSION_SUBMINOR < (subminor)))))) 38 #define PETSC_VERSION_EQUALS(major,minor,subminor) \ 39 ((LIBMESH_DETECTED_PETSC_VERSION_MAJOR == (major)) && \ 40 (LIBMESH_DETECTED_PETSC_VERSION_MINOR == (minor)) && \ 41 (LIBMESH_DETECTED_PETSC_VERSION_SUBMINOR == (subminor))) 47 #define EXTERN_C_FOR_PETSC_BEGIN 48 #define EXTERN_C_FOR_PETSC_END 52 #include <libmesh/ignore_warnings.h> 54 # define LIBMESH_SAW_I 58 # undef I // Avoid complex.h contamination 60 #include <libmesh/restore_warnings.h> 64 #define LibMeshVecDestroy(x) VecDestroy(x) 65 #define LibMeshVecScatterDestroy(x) VecScatterDestroy(x) 66 #define LibMeshMatDestroy(x) MatDestroy(x) 67 #define LibMeshISDestroy(x) ISDestroy(x) 68 #define LibMeshKSPDestroy(x) KSPDestroy(x) 69 #define LibMeshSNESDestroy(x) SNESDestroy(x) 70 #define LibMeshPetscViewerDestroy(x) PetscViewerDestroy(x) 71 #define LibMeshPCDestroy(x) PCDestroy(x) 77 #define LibMeshVecScatterCreate(xin,ix,yin,iy,newctx) VecScatterCreate(xin,ix,yin,iy,newctx) 78 #define ISCreateLibMesh(comm,n,idx,mode,is) ISCreateGeneral((comm),(n),(idx),(mode),(is)) 81 #if PETSC_VERSION_LESS_THAN(3,8,0) 82 # define LibMeshCreateSubMatrix MatGetSubMatrix 84 # define LibMeshCreateSubMatrix MatCreateSubMatrix 89 #if PETSC_VERSION_LESS_THAN(3,17,0) 90 # define LIBMESH_SETERRQ1 SETERRQ1 91 # define LIBMESH_SETERRQ2 SETERRQ2 92 # define LIBMESH_SETERRQ3 SETERRQ3 94 # define LIBMESH_SETERRQ1 SETERRQ 95 # define LIBMESH_SETERRQ2 SETERRQ 96 # define LIBMESH_SETERRQ3 SETERRQ 101 #if PETSC_VERSION_LESS_THAN(3,8,0) 102 # define LIBMESH_PETSCINT_FMT "D" 104 # define LIBMESH_PETSCINT_FMT PetscInt_FMT 108 #if PETSC_VERSION_LESS_THAN(3,19,0) 109 # define LIBMESH_PETSC_NULLPTR PETSC_NULL 111 # define LIBMESH_PETSC_NULLPTR PETSC_NULLPTR 117 #if LIBMESH_DEFAULT_QUADRUPLE_PRECISION 118 # include <boost/multiprecision/float128.hpp> 123 std::ostream &
operator<< (std::ostream & os,
const PetscScalar in)
125 os << (boost::multiprecision::float128(in));
165 template <
typename T>
166 PetscScalar
PS(T val)
168 return val.backend().value();
171 template <
typename T>
172 PetscScalar *
pPS(T * ptr)
174 return &(ptr->backend().value());
177 template <
typename T>
178 const PetscScalar *
pPS(
const T * ptr)
180 return &(ptr->backend().value());
183 template <
typename T>
186 return &(ptr->backend().value());
189 template <
typename T>
190 const PetscReal *
pPR(
const T * ptr)
192 return &(ptr->backend().value());
200 template <
typename T>
201 PetscScalar
PS(T val)
206 template <
typename T>
207 PetscScalar *
pPS(T * ptr)
212 template <
typename T>
213 const PetscScalar *
pPS(
const T * ptr)
218 template <
typename T>
219 PetscReal *
pPR(T * ptr)
224 template <
typename T>
225 const PetscReal *
pPR(
const T * ptr)
231 #endif // LIBMESH_ENABLE_QUADRUPLE_PRECISION 233 #else // LIBMESH_HAVE_PETSC 235 #define PETSC_VERSION_LESS_THAN(major,minor,subminor) 1 236 #define PETSC_VERSION_EQUALS(major,minor,revision) 0 238 #endif // LIBMESH_HAVE_PETSC 242 #ifdef LIBMESH_DETECTED_PETSC_VERSION_RELEASE 244 #define PETSC_RELEASE_LESS_THAN(major, minor, subminor) \ 245 (PETSC_VERSION_LESS_THAN(major, minor, subminor) && LIBMESH_DETECTED_PETSC_VERSION_RELEASE) 246 #define PETSC_RELEASE_EQUALS(major, minor, subminor) \ 247 (PETSC_VERSION_EQUALS(major, minor, subminor) && LIBMESH_DETECTED_PETSC_VERSION_RELEASE) 251 #define PETSC_RELEASE_LESS_THAN(major, minor, subminor) \ 252 (PETSC_VERSION_LESS_THAN(major, minor, subminor)) 253 #define PETSC_RELEASE_EQUALS(major, minor, subminor) (PETSC_VERSION_EQUALS(major, minor, subminor)) 257 #define PETSC_RELEASE_LESS_EQUALS(major, minor, subminor) \ 258 (PETSC_RELEASE_LESS_THAN(major, minor, subminor) || PETSC_RELEASE_EQUALS(major, minor, subminor)) 260 #define PETSC_RELEASE_GREATER_EQUALS(major, minor, subminor) \ 261 (0 == PETSC_RELEASE_LESS_THAN(major, minor, subminor)) 263 #define PETSC_RELEASE_GREATER_THAN(major, minor, subminor) \ 264 (0 == PETSC_RELEASE_LESS_EQUALS(major, minor, subminor)) 266 #endif // LIBMESH_PETSC_MACRO_H
std::ostream & operator<<(std::ostream &os, const PetscScalar in)
PetscScalar * pPS(T *ptr)
The libMesh namespace provides an interface to certain functionality in the library.