Actual source code: petschypre.h
1: #pragma once
3: #include <petscsys.h>
4: #include <petscpkg_version.h>
5: #include <HYPRE_config.h>
6: #include <HYPRE_utilities.h>
8: /* from version 2.16 on, HYPRE_BigInt is 64-bit for 64-bit pointer installations
9: and 32-bit for 32-bit installations -> not the best name for a variable */
10: #if PETSC_PKG_HYPRE_VERSION_LT(2, 16, 0)
11: typedef PetscInt HYPRE_BigInt;
12: #endif
14: #if defined(HYPRE_BIGINT) || defined(HYPRE_MIXEDINT)
15: #define PetscHYPRE_BigInt_FMT "lld"
16: #ifdef __cplusplus /* make sure our format specifiers line up */
17: #include <type_traits>
18: static_assert(std::is_same<HYPRE_BigInt, long long int>::value, "");
19: #endif
20: #else
21: #define PetscHYPRE_BigInt_FMT "d"
22: #ifdef __cplusplus /* make sure our format specifiers line up */
23: #include <type_traits>
24: static_assert(std::is_same<HYPRE_BigInt, int>::value, "");
25: #endif
26: #endif
28: /*
29: With scalar type == real, HYPRE_Complex == PetscScalar;
30: With scalar type == complex, HYPRE_Complex is double __complex__ while PetscScalar may be std::complex<double>
31: */
32: static inline PetscErrorCode PetscHYPREScalarCast(PetscScalar a, HYPRE_Complex *b)
33: {
34: PetscFunctionBegin;
35: #if defined(HYPRE_COMPLEX)
36: ((PetscReal *)b)[0] = PetscRealPart(a);
37: ((PetscReal *)b)[1] = PetscImaginaryPart(a);
38: #else
39: *b = a;
40: #endif
41: PetscFunctionReturn(PETSC_SUCCESS);
42: }
44: #if PETSC_PKG_HYPRE_VERSION_GT(2, 28, 0) || (PETSC_PKG_HYPRE_VERSION_EQ(2, 28, 0) && defined(HYPRE_DEVELOP_NUMBER) && HYPRE_DEVELOP_NUMBER >= 22)
45: static inline PetscErrorCode PetscHYPREFinalize_Private(void)
46: {
47: if (HYPRE_Initialized() && !HYPRE_Finalized()) PetscCallExternal(HYPRE_Finalize, );
48: return PETSC_SUCCESS;
49: }
50: #define PetscHYPREInitialize() \
51: do { \
52: if (!HYPRE_Initialized()) { \
53: PetscCallExternal(HYPRE_Initialize, ); \
54: PetscCall(PetscRegisterFinalize(PetscHYPREFinalize_Private)); \
55: } \
56: } while (0)
57: #else
58: #define PetscHYPREInitialize() (void)0
59: #endif
61: #if PETSC_PKG_HYPRE_VERSION_LT(2, 19, 0)
62: typedef int HYPRE_MemoryLocation;
63: #define hypre_IJVectorMemoryLocation(a) 0
64: #define hypre_IJMatrixMemoryLocation(a) 0
65: #endif