Actual source code: petscsystypes.h
1: /* Portions of this code are under:
2: Copyright (c) 2022 Advanced Micro Devices, Inc. All rights reserved.
3: */
5: #ifndef PETSCSYSTYPES_H
6: #define PETSCSYSTYPES_H
8: #include <petscconf.h>
9: #include <petscconf_poison.h>
10: #include <petscfix.h>
11: #include <petscmacros.h>
12: #include <stddef.h>
14: /* SUBMANSEC = Sys */
16: #include <limits.h> // INT_MIN, INT_MAX
18: #if defined(__clang__) || (PETSC_CPP_VERSION >= 17)
19: // clang allows both [[nodiscard]] and __attribute__((warn_unused_result)) on type
20: // definitions. GCC, however, does not, so check that we are using C++17 [[nodiscard]]
21: // instead of __attribute__((warn_unused_result))
22: #define PETSC_ERROR_CODE_NODISCARD PETSC_NODISCARD
23: #else
24: #define PETSC_ERROR_CODE_NODISCARD
25: #endif
27: #ifdef PETSC_CLANG_STATIC_ANALYZER
28: #undef PETSC_USE_STRICT_PETSCERRORCODE
29: #endif
31: #ifdef PETSC_USE_STRICT_PETSCERRORCODE
32: #define PETSC_ERROR_CODE_TYPEDEF typedef
33: #define PETSC_ERROR_CODE_ENUM_NAME PetscErrorCode
34: #else
35: #define PETSC_ERROR_CODE_TYPEDEF
36: #define PETSC_ERROR_CODE_ENUM_NAME
37: #endif
39: /*E
40: PetscErrorCode - Datatype used to return PETSc error codes.
42: Level: beginner
44: Notes:
45: Virtually all PETSc functions return an error code. It is the callers responsibility to check
46: the value of the returned error code after each PETSc call to determine if any errors
47: occurred. A set of convenience macros (e.g. `PetscCall()`, `PetscCallVoid()`) are provided
48: for this purpose. Failing to properly check for errors is not supported, as errors may leave
49: PETSc in an undetermined state.
51: One can retrieve the error string corresponding to a particular error code using
52: `PetscErrorMessage()`.
54: The user can also configure PETSc with the `--with-strict-petscerrorcode` option to enable
55: compiler warnings when the returned error codes are not captured and checked. Users are
56: *heavily* encouraged to opt-in to this option, as it will become enabled by default in a
57: future release.
59: Developer Notes:
61: These are the generic error codes. These error codes are used in many different places in the
62: PETSc source code. The C-string versions are at defined in `PetscErrorStrings[]` in
63: `src/sys/error/err.c`, while the fortran versions are defined in
64: `src/sys/f90-mod/petscerror.h`. Any changes here must also be made in both locations.
66: .seealso: `PetscErrorMessage()`, `PetscCall()`, `SETERRQ()`
67: E*/
68: PETSC_ERROR_CODE_TYPEDEF enum PETSC_ERROR_CODE_NODISCARD {
69: PETSC_SUCCESS = 0,
70: PETSC_ERR_BOOLEAN_MACRO_FAILURE = 1, /* do not use */
72: PETSC_ERR_MIN_VALUE = 54, /* should always be one less then the smallest value */
74: PETSC_ERR_MEM = 55, /* unable to allocate requested memory */
75: PETSC_ERR_SUP = 56, /* no support for requested operation */
76: PETSC_ERR_SUP_SYS = 57, /* no support for requested operation on this computer system */
77: PETSC_ERR_ORDER = 58, /* operation done in wrong order */
78: PETSC_ERR_SIG = 59, /* signal received */
79: PETSC_ERR_FP = 72, /* floating point exception */
80: PETSC_ERR_COR = 74, /* corrupted PETSc object */
81: PETSC_ERR_LIB = 76, /* error in library called by PETSc */
82: PETSC_ERR_PLIB = 77, /* PETSc library generated inconsistent data */
83: PETSC_ERR_MEMC = 78, /* memory corruption */
84: PETSC_ERR_CONV_FAILED = 82, /* iterative method (KSP or SNES) failed */
85: PETSC_ERR_USER = 83, /* user has not provided needed function */
86: PETSC_ERR_SYS = 88, /* error in system call */
87: PETSC_ERR_POINTER = 70, /* pointer does not point to valid address */
88: PETSC_ERR_MPI_LIB_INCOMP = 87, /* MPI library at runtime is not compatible with MPI user compiled with */
90: PETSC_ERR_ARG_SIZ = 60, /* nonconforming object sizes used in operation */
91: PETSC_ERR_ARG_IDN = 61, /* two arguments not allowed to be the same */
92: PETSC_ERR_ARG_WRONG = 62, /* wrong argument (but object probably ok) */
93: PETSC_ERR_ARG_CORRUPT = 64, /* null or corrupted PETSc object as argument */
94: PETSC_ERR_ARG_OUTOFRANGE = 63, /* input argument, out of range */
95: PETSC_ERR_ARG_BADPTR = 68, /* invalid pointer argument */
96: PETSC_ERR_ARG_NOTSAMETYPE = 69, /* two args must be same object type */
97: PETSC_ERR_ARG_NOTSAMECOMM = 80, /* two args must be same communicators */
98: PETSC_ERR_ARG_WRONGSTATE = 73, /* object in argument is in wrong state, e.g. unassembled mat */
99: PETSC_ERR_ARG_TYPENOTSET = 89, /* the type of the object has not yet been set */
100: PETSC_ERR_ARG_INCOMP = 75, /* two arguments are incompatible */
101: PETSC_ERR_ARG_NULL = 85, /* argument is null that should not be */
102: PETSC_ERR_ARG_UNKNOWN_TYPE = 86, /* type name doesn't match any registered type */
104: PETSC_ERR_FILE_OPEN = 65, /* unable to open file */
105: PETSC_ERR_FILE_READ = 66, /* unable to read from file */
106: PETSC_ERR_FILE_WRITE = 67, /* unable to write to file */
107: PETSC_ERR_FILE_UNEXPECTED = 79, /* unexpected data in file */
109: PETSC_ERR_MAT_LU_ZRPVT = 71, /* detected a zero pivot during LU factorization */
110: PETSC_ERR_MAT_CH_ZRPVT = 81, /* detected a zero pivot during Cholesky factorization */
112: PETSC_ERR_INT_OVERFLOW = 84,
113: PETSC_ERR_FLOP_COUNT = 90,
114: PETSC_ERR_NOT_CONVERGED = 91, /* solver did not converge */
115: PETSC_ERR_MISSING_FACTOR = 92, /* MatGetFactor() failed */
116: PETSC_ERR_OPT_OVERWRITE = 93, /* attempted to over write options which should not be changed */
117: PETSC_ERR_WRONG_MPI_SIZE = 94, /* example/application run with number of MPI ranks it does not support */
118: PETSC_ERR_USER_INPUT = 95, /* missing or incorrect user input */
119: PETSC_ERR_GPU_RESOURCE = 96, /* unable to load a GPU resource, for example cuBLAS */
120: PETSC_ERR_GPU = 97, /* An error from a GPU call, this may be due to lack of resources on the GPU or a true error in the call */
121: PETSC_ERR_MPI = 98, /* general MPI error */
122: PETSC_ERR_RETURN = 99, /* PetscError() incorrectly returned an error code of 0 */
123: PETSC_ERR_MAX_VALUE = 100, /* this is always the one more than the largest error code */
125: /*
126: do not use, exist purely to make the enum bounds equal that of a regular int (so conversion
127: to int in main() is not undefined behavior)
128: */
129: PETSC_ERR_MIN_SIGNED_BOUND_DO_NOT_USE = INT_MIN,
130: PETSC_ERR_MAX_SIGNED_BOUND_DO_NOT_USE = INT_MAX
131: } PETSC_ERROR_CODE_ENUM_NAME;
133: #ifndef PETSC_USE_STRICT_PETSCERRORCODE
134: typedef int PetscErrorCode;
136: /*
137: Needed so that C++ lambdas can deduce the return type as PetscErrorCode from
138: PetscFunctionReturn(PETSC_SUCCESS). Otherwise we get
140: error: return type '(unnamed enum at include/petscsystypes.h:50:1)' must match previous
141: return type 'int' when lambda expression has unspecified explicit return type
142: PetscFunctionReturn(PETSC_SUCCESS);
143: ^
144: */
145: #define PETSC_SUCCESS ((PetscErrorCode)0)
146: #endif
148: #undef PETSC_ERROR_CODE_NODISCARD
149: #undef PETSC_ERROR_CODE_TYPEDEF
150: #undef PETSC_ERROR_CODE_ENUM_NAME
152: /*MC
154: PetscClassId - A unique id used to identify each PETSc class.
156: Notes:
157: Use `PetscClassIdRegister()` to obtain a new value for a new class being created. Usually
158: XXXInitializePackage() calls it for each class it defines.
160: Developer Notes:
161: Internal integer stored in the `_p_PetscObject` data structure. These are all computed by an offset from the lowest one, `PETSC_SMALLEST_CLASSID`.
163: Level: developer
165: .seealso: `PetscClassIdRegister()`, `PetscLogEventRegister()`, `PetscHeaderCreate()`
166: M*/
167: typedef int PetscClassId;
169: /*MC
170: PetscMPIInt - datatype used to represent 'int' parameters to MPI functions.
172: Level: intermediate
174: Notes:
175: This is always a 32 bit integer, sometimes it is the same as `PetscInt`, but if PETSc was built with --with-64-bit-indices but
176: standard C/Fortran integers are 32 bit then this is NOT the same as `PetscInt`; it remains 32 bit.
178: `PetscMPIIntCast`(a,&b) checks if the given `PetscInt` a will fit in a `PetscMPIInt`, if not it
179: generates a `PETSC_ERR_ARG_OUTOFRANGE` error.
181: .seealso: `PetscBLASInt`, `PetscInt`, `PetscMPIIntCast()`
182: M*/
183: typedef int PetscMPIInt;
185: /* Limit MPI to 32-bits */
186: enum {
187: PETSC_MPI_INT_MIN = INT_MIN,
188: PETSC_MPI_INT_MAX = INT_MAX
189: };
191: /*MC
192: PetscSizeT - datatype used to represent sizes in memory (like size_t)
194: Level: intermediate
196: Notes:
197: This is equivalent to size_t, but defined for consistency with Fortran, which lacks a native equivalent of size_t.
199: .seealso: `PetscInt`, `PetscInt64`, `PetscCount`
200: M*/
201: typedef size_t PetscSizeT;
203: /*MC
204: PetscCount - signed datatype used to represent counts
206: Level: intermediate
208: Notes:
209: This is equivalent to ptrdiff_t, but defined for consistency with Fortran, which lacks a native equivalent of ptrdiff_t.
211: Use `PetscCount_FMT` to format with `PetscPrintf()`, `printf()`, and related functions.
213: .seealso: `PetscInt`, `PetscInt64`, `PetscSizeT`
214: M*/
215: typedef ptrdiff_t PetscCount;
216: #define PetscCount_FMT "td"
218: /*MC
219: PetscEnum - datatype used to pass enum types within PETSc functions.
221: Level: intermediate
223: .seealso: `PetscOptionsGetEnum()`, `PetscOptionsEnum()`, `PetscBagRegisterEnum()`
224: M*/
225: typedef enum {
226: ENUM_DUMMY
227: } PetscEnum;
229: typedef short PetscShort;
230: typedef char PetscChar;
231: typedef float PetscFloat;
233: /*MC
234: PetscInt - PETSc type that represents an integer, used primarily to
235: represent size of arrays and indexing into arrays. Its size can be configured with the option --with-64-bit-indices to be either 32-bit (default) or 64-bit.
237: Level: beginner
239: Notes:
240: For MPI calls that require datatypes, use `MPIU_INT` as the datatype for `PetscInt`. It will automatically work correctly regardless of the size of PetscInt.
242: .seealso: `PetscBLASInt`, `PetscMPIInt`, `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`
243: M*/
245: #if defined(PETSC_HAVE_STDINT_H)
246: #include <stdint.h>
247: #endif
248: #if defined(PETSC_HAVE_INTTYPES_H)
251: #endif
252: #include <inttypes.h>
253: #if !defined(PRId64)
254: #define PRId64 "ld"
255: #endif
256: #endif
258: #if defined(PETSC_HAVE_STDINT_H) && defined(PETSC_HAVE_INTTYPES_H) && defined(PETSC_HAVE_MPI_INT64_T) /* MPI_INT64_T is not guaranteed to be a macro */
259: typedef int64_t PetscInt64;
261: #define PETSC_INT64_MIN INT64_MIN
262: #define PETSC_INT64_MAX INT64_MAX
264: #elif (PETSC_SIZEOF_LONG_LONG == 8)
265: typedef long long PetscInt64;
267: #define PETSC_INT64_MIN LLONG_MIN
268: #define PETSC_INT64_MAX LLONG_MAX
270: #elif defined(PETSC_HAVE___INT64)
271: typedef __int64 PetscInt64;
273: #define PETSC_INT64_MIN INT64_MIN
274: #define PETSC_INT64_MAX INT64_MAX
276: #else
277: #error "cannot determine PetscInt64 type"
278: #endif
280: #if defined(PETSC_USE_64BIT_INDICES)
281: typedef PetscInt64 PetscInt;
283: #define PETSC_INT_MIN PETSC_INT64_MIN
284: #define PETSC_INT_MAX PETSC_INT64_MAX
285: #define PetscInt_FMT PetscInt64_FMT
286: #else
287: typedef int PetscInt;
289: enum {
290: PETSC_INT_MIN = INT_MIN,
291: PETSC_INT_MAX = INT_MAX
292: };
294: #define PetscInt_FMT "d"
295: #endif
297: #define PETSC_MIN_INT PETSC_INT_MIN
298: #define PETSC_MAX_INT PETSC_INT_MAX
299: #define PETSC_MAX_UINT16 65535
301: #if defined(PETSC_HAVE_STDINT_H) && defined(PETSC_HAVE_INTTYPES_H) && defined(PETSC_HAVE_MPI_INT64_T) /* MPI_INT64_T is not guaranteed to be a macro */
302: #define MPIU_INT64 MPI_INT64_T
303: #define PetscInt64_FMT PRId64
304: #elif (PETSC_SIZEOF_LONG_LONG == 8)
305: #define MPIU_INT64 MPI_LONG_LONG_INT
306: #define PetscInt64_FMT "lld"
307: #elif defined(PETSC_HAVE___INT64)
308: #define MPIU_INT64 MPI_INT64_T
309: #define PetscInt64_FMT "ld"
310: #else
311: #error "cannot determine PetscInt64 type"
312: #endif
314: /*MC
315: PetscBLASInt - datatype used to represent 'int' parameters to BLAS/LAPACK functions.
317: Level: intermediate
319: Notes:
320: Usually this is the same as `PetscInt`, but if PETSc was built with --with-64-bit-indices but
321: standard C/Fortran integers are 32 bit then this may not be the same as `PetscInt`,
322: except on some BLAS/LAPACK implementations that support 64 bit integers see the notes below.
324: `PetscErrorCode` `PetscBLASIntCast`(a,&b) checks if the given `PetscInt` a will fit in a `PetscBLASInt`, if not it
325: generates a `PETSC_ERR_ARG_OUTOFRANGE` error
327: Installation Notes:
328: ./configure automatically determines the size of the integers used by BLAS/LAPACK except when --with-batch is used
329: in that situation one must know (by some other means) if the integers used by BLAS/LAPACK are 64 bit and if so pass the flag --known-64-bit-blas-indice
331: MATLAB ships with BLAS and LAPACK that use 64 bit integers, for example if you run ./configure with, the option
332: --with-blaslapack-lib=[/Applications/MATLAB_R2010b.app/bin/maci64/libmwblas.dylib,/Applications/MATLAB_R2010b.app/bin/maci64/libmwlapack.dylib]
334: MKL ships with both 32 and 64 bit integer versions of the BLAS and LAPACK. If you pass the flag -with-64-bit-blas-indices PETSc will link
335: against the 64 bit version, otherwise it use the 32 bit version
337: OpenBLAS can be built to use 64 bit integers. The ./configure options --download-openblas -with-64-bit-blas-indices will build a 64 bit integer version
339: External packages such as hypre, ML, SuperLU etc do not provide any support for passing 64 bit integers to BLAS/LAPACK so cannot
340: be used with PETSc when PETSc links against 64 bit integer BLAS/LAPACK. ./configure will generate an error if you attempt to link PETSc against any of
341: these external libraries while using 64 bit integer BLAS/LAPACK.
343: .seealso: `PetscMPIInt`, `PetscInt`, `PetscBLASIntCast()`
344: M*/
345: #if defined(PETSC_HAVE_64BIT_BLAS_INDICES)
346: typedef PetscInt64 PetscBLASInt;
348: #define PETSC_BLAS_INT_MIN PETSC_INT64_MIN
349: #define PETSC_BLAS_INT_MAX PETSC_INT64_MAX
350: #define PetscBLASInt_FMT PetscInt64_FMT
351: #else
352: typedef int PetscBLASInt;
354: enum {
355: PETSC_BLAS_INT_MIN = INT_MIN,
356: PETSC_BLAS_INT_MAX = INT_MAX
357: };
359: #define PetscBLASInt_FMT "d"
360: #endif
362: /*MC
363: PetscCuBLASInt - datatype used to represent 'int' parameters to cuBLAS/cuSOLVER functions.
365: Level: intermediate
367: Notes:
368: As of this writing `PetscCuBLASInt` is always the system `int`.
370: `PetscErrorCode` `PetscCuBLASIntCast`(a,&b) checks if the given `PetscInt` a will fit in a `PetscCuBLASInt`, if not it
371: generates a `PETSC_ERR_ARG_OUTOFRANGE` error
373: .seealso: `PetscBLASInt`, `PetscMPIInt`, `PetscInt`, `PetscCuBLASIntCast()`
374: M*/
375: typedef int PetscCuBLASInt;
377: enum {
378: PETSC_CUBLAS_INT_MIN = INT_MIN,
379: PETSC_CUBLAS_INT_MAX = INT_MAX
380: };
382: /*MC
383: PetscHipBLASInt - datatype used to represent 'int' parameters to hipBLAS/hipSOLVER functions.
385: Level: intermediate
387: Notes:
388: As of this writing `PetscHipBLASInt` is always the system `int`.
390: `PetscErrorCode` `PetscHipBLASIntCast`(a,&b) checks if the given `PetscInt` a will fit in a `PetscHipBLASInt`, if not it
391: generates a `PETSC_ERR_ARG_OUTOFRANGE` error
393: .seealso: PetscBLASInt, PetscMPIInt, PetscInt, PetscHipBLASIntCast()
394: M*/
395: typedef int PetscHipBLASInt;
397: enum {
398: PETSC_HIPBLAS_INT_MIN = INT_MIN,
399: PETSC_HIPBLAS_INT_MAX = INT_MAX
400: };
402: /*E
403: PetscBool - Logical variable. Actually an enum in C and a logical in Fortran.
405: Level: beginner
407: Developer Note:
408: Why have `PetscBool`, why not use bool in C? The problem is that K and R C, C99 and C++ all have different mechanisms for
409: boolean values. It is not easy to have a simple macro that that will work properly in all circumstances with all three mechanisms.
411: .seealso: `PETSC_TRUE`, `PETSC_FALSE`, `PetscNot()`, `PetscBool3`
412: E*/
413: typedef enum {
414: PETSC_FALSE,
415: PETSC_TRUE
416: } PetscBool;
417: PETSC_EXTERN const char *const PetscBools[];
419: /*E
420: PetscBool3 - Ternary logical variable. Actually an enum in C and a 4 byte integer in Fortran.
422: Level: beginner
424: Note:
425: Should not be used with the if (flg) or if (!flg) syntax.
427: .seealso: `PETSC_TRUE`, `PETSC_FALSE`, `PetscNot()`, `PETSC_BOOL3_TRUE`, `PETSC_BOOL3_FALSE`, `PETSC_BOOL3_UNKNOWN`
428: E*/
429: typedef enum {
430: PETSC_BOOL3_FALSE,
431: PETSC_BOOL3_TRUE,
432: PETSC_BOOL3_UNKNOWN = -1
433: } PetscBool3;
435: #define PetscBool3ToBool(a) ((a) == PETSC_BOOL3_TRUE ? PETSC_TRUE : PETSC_FALSE)
436: #define PetscBoolToBool3(a) ((a) == PETSC_TRUE ? PETSC_BOOL3_TRUE : PETSC_BOOL3_FALSE)
438: /*MC
439: PetscReal - PETSc type that represents a real number version of `PetscScalar`
441: Level: beginner
443: Notes:
444: For MPI calls that require datatypes, use `MPIU_REAL` as the datatype for `PetscReal` and `MPIU_SUM`, `MPIU_MAX`, etc. for operations.
445: They will automatically work correctly regardless of the size of `PetscReal`.
447: See `PetscScalar` for details on how to ./configure the size of `PetscReal`.
449: .seealso: `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`
450: M*/
452: #if defined(PETSC_USE_REAL_SINGLE)
453: typedef float PetscReal;
454: #elif defined(PETSC_USE_REAL_DOUBLE)
455: typedef double PetscReal;
456: #elif defined(PETSC_USE_REAL___FLOAT128)
457: #if defined(__cplusplus)
458: extern "C" {
459: #endif
460: #include <quadmath.h>
461: #if defined(__cplusplus)
462: }
463: #endif
464: typedef __float128 PetscReal;
465: #elif defined(PETSC_USE_REAL___FP16)
466: typedef __fp16 PetscReal;
467: #endif /* PETSC_USE_REAL_* */
469: /*MC
470: PetscComplex - PETSc type that represents a complex number with precision matching that of `PetscReal`.
472: Synopsis:
473: #include <petscsys.h>
474: PetscComplex number = 1. + 2.*PETSC_i;
476: Level: beginner
478: Notes:
479: For MPI calls that require datatypes, use `MPIU_COMPLEX` as the datatype for `PetscComplex` and `MPIU_SUM` etc for operations.
480: They will automatically work correctly regardless of the size of `PetscComplex`.
482: See PetscScalar for details on how to ./configure the size of `PetscReal`
484: Complex numbers are automatically available if PETSc was able to find a working complex implementation
486: Petsc has a 'fix' for complex numbers to support expressions such as std::complex<PetscReal> + `PetscInt`, which are not supported by the standard
487: C++ library, but are convenient for petsc users. If the C++ compiler is able to compile code in petsccxxcomplexfix.h (This is checked by
488: configure), we include petsccxxcomplexfix.h to provide this convenience.
490: If the fix causes conflicts, or one really does not want this fix for a particular C++ file, one can define `PETSC_SKIP_CXX_COMPLEX_FIX`
491: at the beginning of the C++ file to skip the fix.
493: .seealso: `PetscReal`, `PetscScalar`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`, `PETSC_i`
494: M*/
495: #if !defined(PETSC_SKIP_COMPLEX)
496: #if defined(PETSC_CLANGUAGE_CXX)
497: #if !defined(PETSC_USE_REAL___FP16) && !defined(PETSC_USE_REAL___FLOAT128)
498: #if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) /* enable complex for library code */
499: #define PETSC_HAVE_COMPLEX 1
500: #elif !defined(__cplusplus) && defined(PETSC_HAVE_C99_COMPLEX) && defined(PETSC_HAVE_CXX_COMPLEX) /* User code only - conditional on library code complex support */
501: #define PETSC_HAVE_COMPLEX 1
502: #endif
503: #elif defined(PETSC_USE_REAL___FLOAT128) && defined(PETSC_HAVE_C99_COMPLEX)
504: #define PETSC_HAVE_COMPLEX 1
505: #endif
506: #else /* !PETSC_CLANGUAGE_CXX */
507: #if !defined(PETSC_USE_REAL___FP16)
509: #define PETSC_HAVE_COMPLEX 1
510: #elif defined(__cplusplus) && defined(PETSC_HAVE_C99_COMPLEX) && defined(PETSC_HAVE_CXX_COMPLEX) /* User code only - conditional on library code complex support */
511: #define PETSC_HAVE_COMPLEX 1
512: #endif
513: #endif
514: #endif /* PETSC_CLANGUAGE_CXX */
515: #endif /* !PETSC_SKIP_COMPLEX */
517: #if defined(PETSC_HAVE_COMPLEX)
518: #if defined(__cplusplus) /* C++ complex support */
519: /* Locate a C++ complex template library */
520: #if defined(PETSC_DESIRE_KOKKOS_COMPLEX) /* Defined in petscvec_kokkos.hpp for *.kokkos.cxx files */
521: #define petsccomplexlib Kokkos
522: #include <Kokkos_Complex.hpp>
523: #elif defined(__CUDACC__) || defined(__HIPCC__)
524: #define petsccomplexlib thrust
525: #include <thrust/complex.h>
526: #elif defined(PETSC_USE_REAL___FLOAT128)
527: #include <complex.h>
528: #else
529: #define petsccomplexlib std
530: #include <complex>
531: #endif
533: /* Define PetscComplex based on the precision */
534: #if defined(PETSC_USE_REAL_SINGLE)
535: typedef petsccomplexlib::complex<float> PetscComplex;
536: #elif defined(PETSC_USE_REAL_DOUBLE)
537: typedef petsccomplexlib::complex<double> PetscComplex;
538: #elif defined(PETSC_USE_REAL___FLOAT128)
539: typedef __complex128 PetscComplex;
540: #endif
542: /* Include a PETSc C++ complex 'fix'. Check PetscComplex manual page for details */
543: #if defined(PETSC_HAVE_CXX_COMPLEX_FIX) && !defined(PETSC_SKIP_CXX_COMPLEX_FIX)
544: #include <petsccxxcomplexfix.h>
545: #endif
546: #else /* c99 complex support */
547: #include <complex.h>
548: #if defined(PETSC_USE_REAL_SINGLE) || defined(PETSC_USE_REAL___FP16)
549: typedef float _Complex PetscComplex;
550: #elif defined(PETSC_USE_REAL_DOUBLE)
551: typedef double _Complex PetscComplex;
552: #elif defined(PETSC_USE_REAL___FLOAT128)
553: typedef __complex128 PetscComplex;
554: #endif /* PETSC_USE_REAL_* */
555: #endif /* !__cplusplus */
556: #endif /* PETSC_HAVE_COMPLEX */
558: /*MC
559: PetscScalar - PETSc type that represents either a double precision real number, a double precision
560: complex number, a single precision real number, a __float128 real or complex or a __fp16 real - if the code is configured
561: with --with-scalar-type=real,complex --with-precision=single,double,__float128,__fp16
563: Notes:
564: For MPI calls that require datatypes, use `MPIU_SCALAR` as the datatype for `PetscScalar` and `MPIU_SUM`, etc for operations. They will automatically work correctly regardless of the size of `PetscScalar`.
566: Level: beginner
568: .seealso: `PetscReal`, `PetscComplex`, `PetscInt`, `MPIU_REAL`, `MPIU_SCALAR`, `MPIU_COMPLEX`, `MPIU_INT`, `PetscRealPart()`, `PetscImaginaryPart()`
569: M*/
571: #if defined(PETSC_USE_COMPLEX) && defined(PETSC_HAVE_COMPLEX)
572: typedef PetscComplex PetscScalar;
573: #else /* PETSC_USE_COMPLEX */
574: typedef PetscReal PetscScalar;
575: #endif /* PETSC_USE_COMPLEX */
577: /*E
578: PetscCopyMode - Determines how an array or `PetscObject` passed to certain functions is copied or retained by the aggregate `PetscObject`
580: Level: beginner
582: Values for array input:
583: + `PETSC_COPY_VALUES` - the array values are copied into new space, the user is free to reuse or delete the passed in array
584: . `PETSC_OWN_POINTER` - the array values are NOT copied, the object takes ownership of the array and will free it later, the user cannot change or
585: delete the array. The array MUST have been obtained with PetscMalloc(). Hence this mode cannot be used in Fortran.
586: - `PETSC_USE_POINTER` - the array values are NOT copied, the object uses the array but does NOT take ownership of the array. The user cannot use
587: the array but the user must delete the array after the object is destroyed.
589: Values for PetscObject:
590: + `PETSC_COPY_VALUES` - the input `PetscObject` is cloned into the aggregate `PetscObject`; the user is free to reuse/modify the input `PetscObject` without side effects.
591: . `PETSC_OWN_POINTER` - the input `PetscObject` is referenced by pointer (with reference count), thus should not be modified by the user.
592: increases its reference count).
593: - `PETSC_USE_POINTER` - invalid for `PetscObject` inputs.
594: E*/
595: typedef enum {
596: PETSC_COPY_VALUES,
597: PETSC_OWN_POINTER,
598: PETSC_USE_POINTER
599: } PetscCopyMode;
600: PETSC_EXTERN const char *const PetscCopyModes[];
602: /*MC
603: PETSC_FALSE - False value of `PetscBool`
605: Level: beginner
607: Note:
608: Zero integer
610: .seealso: `PetscBool`, `PetscBool3`, `PETSC_TRUE`
611: M*/
613: /*MC
614: PETSC_TRUE - True value of `PetscBool`
616: Level: beginner
618: Note:
619: Nonzero integer
621: .seealso: `PetscBool`, `PetscBool3`, `PETSC_FALSE`
622: M*/
624: /*MC
625: PetscLogDouble - Used for logging times
627: Level: developer
629: Note:
630: Contains double precision numbers that are not used in the numerical computations, but rather in logging, timing etc.
632: M*/
633: typedef double PetscLogDouble;
635: /*E
636: PetscDataType - Used for handling different basic data types.
638: Level: beginner
640: Notes:
641: Use of this should be avoided if one can directly use `MPI_Datatype` instead.
643: `PETSC_INT` is the datatype for a `PetscInt`, regardless of whether it is 4 or 8 bytes.
644: `PETSC_REAL`, `PETSC_COMPLEX` and `PETSC_SCALAR` are the datatypes for `PetscReal`, `PetscComplex` and `PetscScalar`, regardless of their sizes.
646: Developer Notes:
647: It would be nice if we could always just use MPI Datatypes, why can we not?
649: If you change any values in `PetscDatatype` make sure you update their usage in
650: share/petsc/matlab/PetscBagRead.m and share/petsc/matlab/@PetscOpenSocket/read/write.m
652: TODO:
653: Add PETSC_INT32 and remove use of improper `PETSC_ENUM`
655: .seealso: `PetscBinaryRead()`, `PetscBinaryWrite()`, `PetscDataTypeToMPIDataType()`,
656: `PetscDataTypeGetSize()`
657: E*/
658: typedef enum {
659: PETSC_DATATYPE_UNKNOWN = 0,
660: PETSC_DOUBLE = 1,
661: PETSC_COMPLEX = 2,
662: PETSC_LONG = 3,
663: PETSC_SHORT = 4,
664: PETSC_FLOAT = 5,
665: PETSC_CHAR = 6,
666: PETSC_BIT_LOGICAL = 7,
667: PETSC_ENUM = 8,
668: PETSC_BOOL = 9,
669: PETSC___FLOAT128 = 10,
670: PETSC_OBJECT = 11,
671: PETSC_FUNCTION = 12,
672: PETSC_STRING = 13,
673: PETSC___FP16 = 14,
674: PETSC_STRUCT = 15,
675: PETSC_INT = 16,
676: PETSC_INT64 = 17,
677: PETSC_COUNT = 18
678: } PetscDataType;
679: PETSC_EXTERN const char *const PetscDataTypes[];
681: #if defined(PETSC_USE_REAL_SINGLE)
682: #define PETSC_REAL PETSC_FLOAT
683: #elif defined(PETSC_USE_REAL_DOUBLE)
684: #define PETSC_REAL PETSC_DOUBLE
685: #elif defined(PETSC_USE_REAL___FLOAT128)
686: #define PETSC_REAL PETSC___FLOAT128
687: #elif defined(PETSC_USE_REAL___FP16)
688: #define PETSC_REAL PETSC___FP16
689: #else
690: #define PETSC_REAL PETSC_DOUBLE
691: #endif
693: #if defined(PETSC_USE_COMPLEX)
694: #define PETSC_SCALAR PETSC_COMPLEX
695: #else
696: #define PETSC_SCALAR PETSC_REAL
697: #endif
699: #define PETSC_FORTRANADDR PETSC_LONG
701: /*S
702: PetscToken - 'Token' used for managing tokenizing strings
704: Level: intermediate
706: .seealso: `PetscTokenCreate()`, `PetscTokenFind()`, `PetscTokenDestroy()`
707: S*/
708: typedef struct _p_PetscToken *PetscToken;
710: /*S
711: PetscObject - any PETSc object, `PetscViewer`, `Mat`, `Vec`, `KSP` etc
713: Level: beginner
715: Notes:
716: This is the base class from which all PETSc objects are derived from.
718: In certain situations one can cast an object, for example a `Vec`, to a `PetscObject` with (`PetscObject`)vec
720: .seealso: `PetscObjectDestroy()`, `PetscObjectView()`, `PetscObjectGetName()`, `PetscObjectSetName()`, `PetscObjectReference()`, `PetscObjectDereference()`
721: S*/
722: typedef struct _p_PetscObject *PetscObject;
724: /*MC
725: PetscObjectId - unique integer Id for a `PetscObject`
727: Level: developer
729: Note:
730: Unlike pointer values, object ids are never reused so one may save a `PetscObjectId` and compare it to one obtained later from a `PetscObject` to determine
731: if the objects are the same. Never compare two object pointer values.
733: .seealso: `PetscObjectState`, `PetscObjectGetId()`
734: M*/
735: typedef PetscInt64 PetscObjectId;
737: /*MC
738: PetscObjectState - integer state for a `PetscObject`
740: Level: developer
742: Notes:
743: Object state is always-increasing and (for objects that track state) can be used to determine if an object has
744: changed since the last time you interacted with it. It is 64-bit so that it will not overflow for a very long time.
746: .seealso: `PetscObjectId`, `PetscObjectStateGet()`, `PetscObjectStateIncrease()`, `PetscObjectStateSet()`
747: M*/
748: typedef PetscInt64 PetscObjectState;
750: /*S
751: PetscFunctionList - Linked list of functions, possibly stored in dynamic libraries, accessed
752: by string name
754: Level: advanced
756: .seealso: `PetscFunctionListAdd()`, `PetscFunctionListDestroy()`
757: S*/
758: typedef struct _n_PetscFunctionList *PetscFunctionList;
760: /*E
761: PetscFileMode - Access mode for a file.
763: Level: beginner
765: Values:
766: + `FILE_MODE_UNDEFINED` - initial invalid value
767: . `FILE_MODE_READ` - open a file at its beginning for reading
768: . `FILE_MODE_WRITE` - open a file at its beginning for writing (will create if the file does not exist)
769: . `FILE_MODE_APPEND` - open a file at end for writing
770: . `FILE_MODE_UPDATE` - open a file for updating, meaning for reading and writing
771: - `FILE_MODE_APPEND_UPDATE` - open a file for updating, meaning for reading and writing, at the end
773: .seealso: `PetscViewerFileSetMode()`
774: E*/
775: typedef enum {
776: FILE_MODE_UNDEFINED = -1,
777: FILE_MODE_READ = 0,
778: FILE_MODE_WRITE,
779: FILE_MODE_APPEND,
780: FILE_MODE_UPDATE,
781: FILE_MODE_APPEND_UPDATE
782: } PetscFileMode;
783: PETSC_EXTERN const char *const PetscFileModes[];
785: typedef void *PetscDLHandle;
786: typedef enum {
787: PETSC_DL_DECIDE = 0,
788: PETSC_DL_NOW = 1,
789: PETSC_DL_LOCAL = 2
790: } PetscDLMode;
792: /*S
793: PetscObjectList - Linked list of PETSc objects, each accessible by string name
795: Level: developer
797: Note:
798: Used by `PetscObjectCompose()` and `PetscObjectQuery()`
800: .seealso: `PetscObjectListAdd()`, `PetscObjectListDestroy()`, `PetscObjectListFind()`, `PetscObjectCompose()`, `PetscObjectQuery()`, `PetscFunctionList`
801: S*/
802: typedef struct _n_PetscObjectList *PetscObjectList;
804: /*S
805: PetscDLLibrary - Linked list of dynamic libraries to search for functions
807: Level: advanced
809: .seealso: `PetscDLLibraryOpen()`
810: S*/
811: typedef struct _n_PetscDLLibrary *PetscDLLibrary;
813: /*S
814: PetscContainer - Simple PETSc object that contains a pointer to any required data
816: Level: advanced
818: Note:
819: This is useful to attach arbitrary data to a `PetscObject` with `PetscObjectCompose()` and `PetscObjectQuery()`
821: .seealso: `PetscObject`, `PetscContainerCreate()`, `PetscObjectCompose()`, `PetscObjectQuery()`
822: S*/
823: typedef struct _p_PetscContainer *PetscContainer;
825: /*S
826: PetscRandom - Abstract PETSc object that manages generating random numbers
828: Level: intermediate
830: .seealso: `PetscRandomCreate()`, `PetscRandomGetValue()`, `PetscRandomType`
831: S*/
832: typedef struct _p_PetscRandom *PetscRandom;
834: /*
835: In binary files variables are stored using the following lengths,
836: regardless of how they are stored in memory on any one particular
837: machine. Use these rather then sizeof() in computing sizes for
838: PetscBinarySeek().
839: */
840: #define PETSC_BINARY_INT_SIZE (32 / 8)
841: #define PETSC_BINARY_FLOAT_SIZE (32 / 8)
842: #define PETSC_BINARY_CHAR_SIZE (8 / 8)
843: #define PETSC_BINARY_SHORT_SIZE (16 / 8)
844: #define PETSC_BINARY_DOUBLE_SIZE (64 / 8)
845: #define PETSC_BINARY_SCALAR_SIZE sizeof(PetscScalar)
847: /*E
848: PetscBinarySeekType - argument to `PetscBinarySeek()`
850: Level: advanced
852: .seealso: `PetscBinarySeek()`, `PetscBinarySynchronizedSeek()`
853: E*/
854: typedef enum {
855: PETSC_BINARY_SEEK_SET = 0,
856: PETSC_BINARY_SEEK_CUR = 1,
857: PETSC_BINARY_SEEK_END = 2
858: } PetscBinarySeekType;
860: /*E
861: PetscBuildTwoSidedType - algorithm for setting up two-sided communication
863: Values:
864: + `PETSC_BUILDTWOSIDED_ALLREDUCE` - classical algorithm using an `MPI_Allreduce()` with
865: a buffer of length equal to the communicator size. Not memory-scalable due to
866: the large reduction size. Requires only MPI-1.
867: . `PETSC_BUILDTWOSIDED_IBARRIER` - nonblocking algorithm based on `MPI_Issend()` and `MPI_Ibarrier()`.
868: Proved communication-optimal in Hoefler, Siebert, and Lumsdaine (2010). Requires MPI-3.
869: - `PETSC_BUILDTWOSIDED_REDSCATTER` - similar to above, but use more optimized function
870: that only communicates the part of the reduction that is necessary. Requires MPI-2.
872: Level: developer
874: .seealso: `PetscCommBuildTwoSided()`, `PetscCommBuildTwoSidedSetType()`, `PetscCommBuildTwoSidedGetType()`
875: E*/
876: typedef enum {
877: PETSC_BUILDTWOSIDED_NOTSET = -1,
878: PETSC_BUILDTWOSIDED_ALLREDUCE = 0,
879: PETSC_BUILDTWOSIDED_IBARRIER = 1,
880: PETSC_BUILDTWOSIDED_REDSCATTER = 2
881: /* Updates here must be accompanied by updates in finclude/petscsys.h and the string array in mpits.c */
882: } PetscBuildTwoSidedType;
883: PETSC_EXTERN const char *const PetscBuildTwoSidedTypes[];
885: /* NOTE: If you change this, you must also change the values in src/vec/f90-mod/petscvec.h */
886: /*E
887: InsertMode - How the entries are combined with the current values in the vectors or matrices
889: Level: beginner
891: Values:
892: + `NOT_SET_VALUES` - do not actually use the values
893: . `INSERT_VALUES` - replace the current values with the provided values, unless the index is marked as constrained by the `PetscSection`
894: . `ADD_VALUES` - add the values to the current values, unless the index is marked as constrained by the `PetscSection`
895: . `MAX_VALUES` - use the maximum of each current value and provided value
896: . `MIN_VALUES` - use the minimum of each current value and provided value
897: . `INSERT_ALL_VALUES` - insert, even if indices that are not marked as constrained by the `PetscSection`
898: . `ADD_ALL_VALUES` - add, even if indices that are not marked as constrained by the `PetscSection`
899: . `INSERT_BC_VALUES` - insert, but ignore indices that are not marked as constrained by the `PetscSection`
900: - `ADD_BC_VALUES` - add, but ignore indices that are not marked as constrained by the `PetscSection`
902: Note:
903: The `PetscSection` that determines the effects of the `InsertMode` values can be obtained by the `Vec` object with `VecGetDM()`
904: and `DMGetLocalSection()`.
906: Not all options are supported for all operations or PETSc object types.
908: .seealso: `VecSetValues()`, `MatSetValues()`, `VecSetValue()`, `VecSetValuesBlocked()`,
909: `VecSetValuesLocal()`, `VecSetValuesBlockedLocal()`, `MatSetValuesBlocked()`,
910: `MatSetValuesBlockedLocal()`, `MatSetValuesLocal()`, `VecScatterBegin()`, `VecScatterEnd()`
911: E*/
912: typedef enum {
913: NOT_SET_VALUES,
914: INSERT_VALUES,
915: ADD_VALUES,
916: MAX_VALUES,
917: MIN_VALUES,
918: INSERT_ALL_VALUES,
919: ADD_ALL_VALUES,
920: INSERT_BC_VALUES,
921: ADD_BC_VALUES
922: } InsertMode;
924: /*MC
925: INSERT_VALUES - Put a value into a vector or matrix, overwrites any previous value
927: Level: beginner
929: .seealso: `InsertMode`, `VecSetValues()`, `MatSetValues()`, `VecSetValue()`, `VecSetValuesBlocked()`,
930: `VecSetValuesLocal()`, `VecSetValuesBlockedLocal()`, `MatSetValuesBlocked()`, `ADD_VALUES`,
931: `MatSetValuesBlockedLocal()`, `MatSetValuesLocal()`, `VecScatterBegin()`, `VecScatterEnd()`, `MAX_VALUES`
932: M*/
934: /*MC
935: ADD_VALUES - Adds a value into a vector or matrix, if there previously was no value, just puts the
936: value into that location
938: Level: beginner
940: .seealso: `InsertMode`, `VecSetValues()`, `MatSetValues()`, `VecSetValue()`, `VecSetValuesBlocked()`,
941: `VecSetValuesLocal()`, `VecSetValuesBlockedLocal()`, `MatSetValuesBlocked()`, `INSERT_VALUES`,
942: `MatSetValuesBlockedLocal()`, `MatSetValuesLocal()`, `VecScatterBegin()`, `VecScatterEnd()`, `MAX_VALUES`
943: M*/
945: /*MC
946: MAX_VALUES - Puts the maximum of the scattered/gathered value and the current value into each location
948: Level: beginner
950: .seealso: `InsertMode`, `VecScatterBegin()`, `VecScatterEnd()`, `ADD_VALUES`, `INSERT_VALUES`
951: M*/
953: /*MC
954: MIN_VALUES - Puts the minimal of the scattered/gathered value and the current value into each location
956: Level: beginner
958: .seealso: `InsertMode`, `VecScatterBegin()`, `VecScatterEnd()`, `ADD_VALUES`, `INSERT_VALUES`
959: M*/
961: /*S
962: PetscSubcomm - A decomposition of an MPI communicator into subcommunicators
964: Values:
965: + `PETSC_SUBCOMM_GENERAL` - similar to `MPI_Comm_split()` each process sets the new communicator (color) they will belong to and the order within that communicator
966: . `PETSC_SUBCOMM_CONTIGUOUS` - each new communicator contains a set of process with contiguous ranks in the original MPI communicator
967: - `PETSC_SUBCOMM_INTERLACED` - each new communictor contains a set of processes equally far apart in rank from the others in that new communicator
969: Notes:
970: After a call to `PetscSubcommSetType()`, `PetscSubcommSetTypeGeneral()`, or `PetscSubcommSetFromOptions()` one may call
971: .vb
972: `PetscSubcommChild()` returns the associated subcommunicator on this process
973: `PetscSubcommContiguousParent()` returns a parent communitor but with all child of the same subcommunicator having contiguous rank
974: .ve
976: Sample Usage:
977: .vb
978: `PetscSubcommCreate()`
979: `PetscSubcommSetNumber()`
980: `PetscSubcommSetType`(`PETSC_SUBCOMM_INTERLACED`);
981: ccomm = `PetscSubcommChild()`
982: `PetscSubcommDestroy()`
983: .ve
985: Example:
986: Consider a communicator with six processes split into 3 subcommunicators.
987: .vb
988: `PETSC_SUBCOMM_CONTIGUOUS` - the first communicator contains rank 0,1 the second rank 2,3 and the third rank 4,5 in the original ordering of the original communicator
989: `PETSC_SUBCOMM_INTERLACED` - the first communicator contains rank 0,3, the second 1,4 and the third 2,5
990: .ve
992: Level: advanced
994: Developer Note:
995: This is used in objects such as `PCREDUNDANT` to manage the subcommunicators on which the redundant computations
996: are performed.
998: .seealso: `PetscSubcommCreate()`, `PetscSubcommSetNumber()`, `PetscSubcommSetType()`, `PetscSubcommView()`, `PetscSubcommSetFromOptions()`
999: S*/
1000: typedef struct _n_PetscSubcomm *PetscSubcomm;
1001: typedef enum {
1002: PETSC_SUBCOMM_GENERAL = 0,
1003: PETSC_SUBCOMM_CONTIGUOUS = 1,
1004: PETSC_SUBCOMM_INTERLACED = 2
1005: } PetscSubcommType;
1006: PETSC_EXTERN const char *const PetscSubcommTypes[];
1008: /*S
1009: PetscHeap - A simple class for managing heaps
1011: Level: intermediate
1013: .seealso: `PetscHeapCreate()`, `PetscHeapAdd()`, `PetscHeapPop()`, `PetscHeapPeek()`, `PetscHeapStash()`, `PetscHeapUnstash()`, `PetscHeapView()`, `PetscHeapDestroy()`
1014: S*/
1015: typedef struct _PetscHeap *PetscHeap;
1017: typedef struct _n_PetscShmComm *PetscShmComm;
1018: typedef struct _n_PetscOmpCtrl *PetscOmpCtrl;
1020: /*S
1021: PetscSegBuffer - a segmented extendable buffer
1023: Level: developer
1025: .seealso: `PetscSegBufferCreate()`, `PetscSegBufferGet()`, `PetscSegBufferExtract()`, `PetscSegBufferDestroy()`
1026: S*/
1027: typedef struct _n_PetscSegBuffer *PetscSegBuffer;
1029: typedef struct _n_PetscOptionsHelpPrinted *PetscOptionsHelpPrinted;
1030: #endif