Actual source code: petscsystypes.h
petsc-3.11.4 2019-09-28
1: #if !defined(PETSCSYSTYPES_H)
2: #define PETSCSYSTYPES_H
4: #include <petscconf.h>
5: #include <petscfix.h>
7: /*MC
8: PetscErrorCode - datatype used for return error code from almost all PETSc functions
10: Level: beginner
12: .seealso: CHKERRQ, SETERRQ
13: M*/
14: typedef int PetscErrorCode;
16: /*MC
18: PetscClassId - A unique id used to identify each PETSc class.
20: Notes:
21: Use PetscClassIdRegister() to obtain a new value for a new class being created. Usually
22: XXXInitializePackage() calls it for each class it defines.
24: Developer Notes:
25: Internal integer stored in the _p_PetscObject data structure.
26: These are all computed by an offset from the lowest one, PETSC_SMALLEST_CLASSID.
28: Level: developer
30: .seealso: PetscClassIdRegister(), PetscLogEventRegister(), PetscHeaderCreate()
31: M*/
32: typedef int PetscClassId;
34: /*MC
35: PetscMPIInt - datatype used to represent 'int' parameters to MPI functions.
37: Level: intermediate
39: Notes:
40: usually this is the same as PetscInt, but if PETSc was built with --with-64-bit-indices but
41: standard C/Fortran integers are 32 bit then this is NOT the same as PetscInt; it remains 32 bit.
43: PetscMPIIntCast(a,&b) checks if the given PetscInt a will fit in a PetscMPIInt, if not it
44: generates a PETSC_ERR_ARG_OUTOFRANGE error.
46: .seealso: PetscBLASInt, PetscInt, PetscMPIIntCast()
48: M*/
49: typedef int PetscMPIInt;
51: /*MC
52: PetscEnum - datatype used to pass enum types within PETSc functions.
54: Level: intermediate
56: .seealso: PetscOptionsGetEnum(), PetscOptionsEnum(), PetscBagRegisterEnum()
57: M*/
58: typedef enum { ENUM_DUMMY } PetscEnum;
60: typedef short PetscShort;
61: typedef char PetscChar;
62: typedef float PetscFloat;
64: /*MC
65: PetscInt - PETSc type that represents an integer, used primarily to
66: 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.
68: Notes:
69: 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.
71: Level: beginner
73: .seealso: PetscBLASInt, PetscMPIInt, PetscReal, PetscScalar, PetscComplex, PetscInt, MPIU_REAL, MPIU_SCALAR, MPIU_COMPLEX, MPIU_INT
74: M*/
76: #if defined(PETSC_HAVE_STDINT_H)
77: # include <stdint.h>
78: #endif
79: #if defined (PETSC_HAVE_INTTYPES_H)
80: # if !defined(__STDC_FORMAT_MACROS)
81: # define __STDC_FORMAT_MACROS /* required for using PRId64 from c++ */
82: # endif
83: # include <inttypes.h>
84: # if !defined(PRId64)
85: # define PRId64 "ld"
86: # endif
87: #endif
89: #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 */
90: typedef int64_t PetscInt64;
91: #elif (PETSC_SIZEOF_LONG_LONG == 8)
92: typedef long long PetscInt64;
93: #elif defined(PETSC_HAVE___INT64)
94: typedef __int64 PetscInt64;
95: #else
96: # error "cannot determine PetscInt64 type"
97: #endif
99: #if defined(PETSC_USE_64BIT_INDICES)
100: typedef PetscInt64 PetscInt;
101: #else
102: typedef int PetscInt;
103: #endif
105: /*MC
106: PetscBLASInt - datatype used to represent 'int' parameters to BLAS/LAPACK functions.
108: Notes:
109: Usually this is the same as PetscInt, but if PETSc was built with --with-64-bit-indices but
110: standard C/Fortran integers are 32 bit then this is NOT the same as PetscInt it remains 32 bit
111: (except on very rare BLAS/LAPACK implementations that support 64 bit integers see the note below).
113: PetscErrorCode PetscBLASIntCast(a,&b) checks if the given PetscInt a will fit in a PetscBLASInt, if not it
114: generates a PETSC_ERR_ARG_OUTOFRANGE error
116: Installation Notes:
117: The 64bit versions of MATLAB ship with BLAS and LAPACK that use 64 bit integers for sizes etc,
118: if you run ./configure with the option
119: --with-blaslapack-lib=[/Applications/MATLAB_R2010b.app/bin/maci64/libmwblas.dylib,/Applications/MATLAB_R2010b.app/bin/maci64/libmwlapack.dylib]
120: but you need to also use --known-64-bit-blas-indices.
122: MKL also ships with 64 bit integer versions of the BLAS and LAPACK, if you select those you must also ./configure with
123: --known-64-bit-blas-indices
125: OpenBLAS can also be built to use 64 bit integers. The ./configure options --download-openblas -download-openblas-64-bit-blas-indices
126: will build a 64 bit integer version
128: Developer Notes:
129: Eventually ./configure should automatically determine the size of the integers used by BLAS/LAPACK.
131: External packages such as hypre, ML, SuperLU etc do not provide any support for passing 64 bit integers to BLAS/LAPACK so cannot
132: be used with PETSc if you have set PetscBLASInt to long int.
134: Level: intermediate
136: .seealso: PetscMPIInt, PetscInt, PetscBLASIntCast()
138: M*/
139: #if defined(PETSC_HAVE_64BIT_BLAS_INDICES)
140: typedef PetscInt64 PetscBLASInt;
141: #else
142: typedef int PetscBLASInt;
143: #endif
145: /*E
146: PetscBool - Logical variable. Actually an int in C and a logical in Fortran.
148: Level: beginner
150: Developer Note:
151: 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
152: boolean values. It is not easy to have a simple macro that that will work properly in all circumstances with all three mechanisms.
154: .seealso: PETSC_TRUE, PETSC_FALSE, PetscNot()
155: E*/
156: typedef enum { PETSC_FALSE,PETSC_TRUE } PetscBool;
158: /*MC
159: PetscReal - PETSc type that represents a real number version of PetscScalar
162: Notes:
163: For MPI calls that require datatypes, use MPIU_REAL as the datatype for PetscScalar and MPIU_SUM, MPIU_MAX, etc. for operations.
164: They will automatically work correctly regardless of the size of PetscReal.
166: See PetscScalar for details on how to ./configure the size of PetscReal.
168: Level: beginner
170: .seealso: PetscScalar, PetscComplex, PetscInt, MPIU_REAL, MPIU_SCALAR, MPIU_COMPLEX, MPIU_INT
171: M*/
173: #if defined(PETSC_USE_REAL_SINGLE)
174: typedef float PetscReal;
175: #elif defined(PETSC_USE_REAL_DOUBLE)
176: typedef double PetscReal;
177: #elif defined(PETSC_USE_REAL___FLOAT128)
178: # if defined(__cplusplus)
179: extern "C" {
180: # endif
181: # include <quadmath.h>
182: # if defined(__cplusplus)
183: }
184: # endif
185: typedef __float128 PetscReal;
186: #elif defined(PETSC_USE_REAL___FP16)
187: typedef __fp16 PetscReal;
188: #endif /* PETSC_USE_REAL_* */
190: /*MC
191: PetscComplex - PETSc type that represents a complex number with precision matching that of PetscReal.
193: Synopsis:
194: #include <petscsys.h>
195: PetscComplex number = 1. + 2.*PETSC_i;
197: Notes:
198: For MPI calls that require datatypes, use MPIU_COMPLEX as the datatype for PetscComplex and MPIU_SUM etc for operations.
199: They will automatically work correctly regardless of the size of PetscComplex.
201: See PetscScalar for details on how to ./configure the size of PetscReal
203: Complex numbers are automatically available if PETSc was able to find a working complex implementation
205: Level: beginner
207: .seealso: PetscReal, PetscScalar, PetscComplex, PetscInt, MPIU_REAL, MPIU_SCALAR, MPIU_COMPLEX, MPIU_INT, PETSC_i
208: M*/
210: #if defined(__cplusplus) && defined(PETSC_HAVE_CXX_COMPLEX) && !defined(PETSC_USE_REAL___FLOAT128)
211: # if !defined(PETSC_SKIP_COMPLEX)
212: /* C++ support of complex number */
213: # define PETSC_HAVE_COMPLEX 1
214: # if defined(PETSC_HAVE_CUDA) && __CUDACC_VER_MAJOR__ > 6
215: /* complex headers in thrust only available in CUDA 7.0 and above */
216: # define petsccomplexlib thrust
217: # include <thrust/complex.h>
218: # else
219: # define petsccomplexlib std
220: # include <complex>
221: # endif
222: # if defined(PETSC_USE_REAL_SINGLE)
223: typedef petsccomplexlib::complex<float> PetscComplex;
224: # elif defined(PETSC_USE_REAL_DOUBLE)
225: typedef petsccomplexlib::complex<double> PetscComplex;
226: # elif defined(PETSC_USE_REAL___FLOAT128)
227: typedef petsccomplexlib::complex<__float128> PetscComplex; /* Notstandard and not expected to work, use __complex128 */
228: # endif /* PETSC_USE_REAL_ */
229: # endif /* ! PETSC_SKIP_COMPLEX */
230: #elif defined(PETSC_HAVE_C99_COMPLEX) && !defined(PETSC_USE_REAL___FP16)
231: # if !defined(PETSC_SKIP_COMPLEX)
232: # define PETSC_HAVE_COMPLEX 1
233: # include <complex.h>
234: # if defined(PETSC_USE_REAL_SINGLE) || defined(PETSC_USE_REAL___FP16)
235: typedef float _Complex PetscComplex;
236: # elif defined(PETSC_USE_REAL_DOUBLE)
237: typedef double _Complex PetscComplex;
238: # elif defined(PETSC_USE_REAL___FLOAT128)
239: typedef __complex128 PetscComplex;
240: # endif /* PETSC_USE_REAL_* */
241: # endif /* !PETSC_SKIP_COMPLEX */
242: #elif (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX))
243: # error "PETSc was configured --with-scalar-type=complex, but a language-appropriate complex library is not available"
244: #endif /* !PETSC_SKIP_COMPLEX */
246: /*MC
247: PetscScalar - PETSc type that represents either a double precision real number, a double precision
248: complex number, a single precision real number, a __float128 real or complex or a __fp16 real - if the code is configured
249: with --with-scalar-type=real,complex --with-precision=single,double,__float128,__fp16
251: Notes:
252: For MPI calls that require datatypes, use MPIU_SCALAR as the datatype for PetscScalar and MPIU_SUM, MPIU_MAX etc for operations. They will automatically work correctly regardless of the size of PetscScalar.
254: Level: beginner
256: .seealso: PetscReal, PetscComplex, PetscInt, MPIU_REAL, MPIU_SCALAR, MPIU_COMPLEX, MPIU_INT, PetscRealPart(), PetscImaginaryPart()
257: M*/
259: #if (defined(PETSC_USE_COMPLEX) && !defined(PETSC_SKIP_COMPLEX))
260: typedef PetscComplex PetscScalar;
261: #else /* PETSC_USE_COMPLEX */
262: typedef PetscReal PetscScalar;
263: #endif /* PETSC_USE_COMPLEX */
265: /*E
266: PetscCopyMode - Determines how an array passed to certain functions is copied or retained
268: Level: beginner
270: $ PETSC_COPY_VALUES - the array values are copied into new space, the user is free to reuse or delete the passed in array
271: $ 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
272: $ delete the array. The array MUST have been obtained with PetscMalloc(). Hence this mode cannot be used in Fortran.
273: $ 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
274: the array but the user must delete the array after the object is destroyed.
276: E*/
277: typedef enum {PETSC_COPY_VALUES, PETSC_OWN_POINTER, PETSC_USE_POINTER} PetscCopyMode;
279: /*MC
280: PETSC_FALSE - False value of PetscBool
282: Level: beginner
284: Note:
285: Zero integer
287: .seealso: PetscBool, PETSC_TRUE
288: M*/
290: /*MC
291: PETSC_TRUE - True value of PetscBool
293: Level: beginner
295: Note:
296: Nonzero integer
298: .seealso: PetscBool, PETSC_FALSE
299: M*/
301: /*MC
302: PetscLogDouble - Used for logging times
304: Notes:
305: Contains double precision numbers that are not used in the numerical computations, but rather in logging, timing etc.
307: Level: developer
309: M*/
310: typedef double PetscLogDouble;
312: /*E
313: PetscDataType - Used for handling different basic data types.
315: Level: beginner
317: Notes:
318: Use of this should be avoided if one can directly use MPI_Datatype instead.
320: Developer comment:
321: It would be nice if we could always just use MPI Datatypes, why can we not?
323: .seealso: PetscBinaryRead(), PetscBinaryWrite(), PetscDataTypeToMPIDataType(),
324: PetscDataTypeGetSize()
326: E*/
327: typedef enum {PETSC_DATATYPE_UNKNOWN = 0,
328: PETSC_DOUBLE = 1, PETSC_COMPLEX = 2, PETSC_LONG = 3, PETSC_SHORT = 4, PETSC_FLOAT = 5,
329: PETSC_CHAR = 6, PETSC_BIT_LOGICAL = 7, PETSC_ENUM = 8, PETSC_BOOL = 9, PETSC___FLOAT128 = 10,
330: PETSC_OBJECT = 11, PETSC_FUNCTION = 12, PETSC_STRING = 13, PETSC___FP16 = 14, PETSC_STRUCT = 15,
331: PETSC_INT = 16, PETSC_INT64 = 17} PetscDataType;
333: #if defined(PETSC_USE_COMPLEX)
334: # define PETSC_SCALAR PETSC_COMPLEX
335: #else
336: # if defined(PETSC_USE_REAL_SINGLE)
337: # define PETSC_SCALAR PETSC_FLOAT
338: # elif defined(PETSC_USE_REAL___FLOAT128)
339: # define PETSC_SCALAR PETSC___FLOAT128
340: # elif defined(PETSC_USE_REAL___FP16)
341: # define PETSC_SCALAR PETSC___FP16
342: # else
343: # define PETSC_SCALAR PETSC_DOUBLE
344: # endif
345: #endif
346: #if defined(PETSC_USE_REAL_SINGLE)
347: # define PETSC_REAL PETSC_FLOAT
348: #elif defined(PETSC_USE_REAL___FLOAT128)
349: # define PETSC_REAL PETSC___FLOAT128
350: #elif defined(PETSC_USE_REAL___FP16)
351: # define PETSC_REAL PETSC___FP16
352: #else
353: # define PETSC_REAL PETSC_DOUBLE
354: #endif
355: #define PETSC_FORTRANADDR PETSC_LONG
357: /*S
358: PetscToken - 'Token' used for managing tokenizing strings
360: Level: intermediate
362: .seealso: PetscTokenCreate(), PetscTokenFind(), PetscTokenDestroy()
363: S*/
364: typedef struct _p_PetscToken* PetscToken;
366: /*S
367: PetscObject - any PETSc object, PetscViewer, Mat, Vec, KSP etc
369: Level: beginner
371: Note:
372: This is the base class from which all PETSc objects are derived from.
374: .seealso: PetscObjectDestroy(), PetscObjectView(), PetscObjectGetName(), PetscObjectSetName(), PetscObjectReference(), PetscObjectDereference()
375: S*/
376: typedef struct _p_PetscObject* PetscObject;
378: /*MC
379: PetscObjectId - unique integer Id for a PetscObject
381: Level: developer
383: Notes:
384: Unlike pointer values, object ids are never reused.
386: .seealso: PetscObjectState, PetscObjectGetId()
387: M*/
388: #if defined(PETSC_USING_F90) && !defined(PETSC_USE_FORTRANKIND) /* compaq F90 */
389: typedef int PetscObjectId;
390: #else
391: typedef PetscInt64 PetscObjectId;
392: #endif
394: /*MC
395: PetscObjectState - integer state for a PetscObject
397: Level: developer
399: Notes:
400: Object state is always-increasing and (for objects that track state) can be used to determine if an object has
401: changed since the last time you interacted with it. It is 64-bit so that it will not overflow for a very long time.
403: .seealso: PetscObjectId, PetscObjectStateGet(), PetscObjectStateIncrease(), PetscObjectStateSet()
404: M*/
405: #if defined(PETSC_USING_F90) && !defined(PETSC_USE_FORTRANKIND) /* compaq F90 */
406: typedef int PetscObjectState;
407: #else
408: typedef PetscInt64 PetscObjectState;
409: #endif
411: /*S
412: PetscFunctionList - Linked list of functions, possibly stored in dynamic libraries, accessed
413: by string name
415: Level: advanced
417: .seealso: PetscFunctionListAdd(), PetscFunctionListDestroy(), PetscOpFlist
418: S*/
419: typedef struct _n_PetscFunctionList *PetscFunctionList;
421: /*E
422: PetscFileMode - Access mode for a file.
424: Level: beginner
426: $ FILE_MODE_READ - open a file at its beginning for reading
427: $ FILE_MODE_WRITE - open a file at its beginning for writing (will create if the file does not exist)
428: $ FILE_MODE_APPEND - open a file at end for writing
429: $ FILE_MODE_UPDATE - open a file for updating, meaning for reading and writing
430: $ FILE_MODE_APPEND_UPDATE - open a file for updating, meaning for reading and writing, at the end
432: .seealso: PetscViewerFileSetMode()
433: E*/
434: typedef enum {FILE_MODE_READ, FILE_MODE_WRITE, FILE_MODE_APPEND, FILE_MODE_UPDATE, FILE_MODE_APPEND_UPDATE} PetscFileMode;
436: typedef void* PetscDLHandle;
437: typedef enum {PETSC_DL_DECIDE=0,PETSC_DL_NOW=1,PETSC_DL_LOCAL=2} PetscDLMode;
439: /*S
440: PetscObjectList - Linked list of PETSc objects, each accessable by string name
442: Level: developer
444: Notes:
445: Used by PetscObjectCompose() and PetscObjectQuery()
447: .seealso: PetscObjectListAdd(), PetscObjectListDestroy(), PetscObjectListFind(), PetscObjectCompose(), PetscObjectQuery(), PetscFunctionList
448: S*/
449: typedef struct _n_PetscObjectList *PetscObjectList;
451: /*S
452: PetscDLLibrary - Linked list of dynamics libraries to search for functions
454: Level: advanced
456: .seealso: PetscDLLibraryOpen()
457: S*/
458: typedef struct _n_PetscDLLibrary *PetscDLLibrary;
460: /*S
461: PetscContainer - Simple PETSc object that contains a pointer to any required data
463: Level: advanced
465: .seealso: PetscObject, PetscContainerCreate()
466: S*/
467: typedef struct _p_PetscContainer* PetscContainer;
469: /*S
470: PetscRandom - Abstract PETSc object that manages generating random numbers
472: Level: intermediate
474: Concepts: random numbers
476: .seealso: PetscRandomCreate(), PetscRandomGetValue(), PetscRandomType
477: S*/
478: typedef struct _p_PetscRandom* PetscRandom;
480: /*
481: In binary files variables are stored using the following lengths,
482: regardless of how they are stored in memory on any one particular
483: machine. Use these rather then sizeof() in computing sizes for
484: PetscBinarySeek().
485: */
486: #define PETSC_BINARY_INT_SIZE (32/8)
487: #define PETSC_BINARY_FLOAT_SIZE (32/8)
488: #define PETSC_BINARY_CHAR_SIZE (8/8)
489: #define PETSC_BINARY_SHORT_SIZE (16/8)
490: #define PETSC_BINARY_DOUBLE_SIZE (64/8)
491: #define PETSC_BINARY_SCALAR_SIZE sizeof(PetscScalar)
493: /*E
494: PetscBinarySeekType - argument to PetscBinarySeek()
496: Level: advanced
498: .seealso: PetscBinarySeek(), PetscBinarySynchronizedSeek()
499: E*/
500: typedef enum {PETSC_BINARY_SEEK_SET = 0,PETSC_BINARY_SEEK_CUR = 1,PETSC_BINARY_SEEK_END = 2} PetscBinarySeekType;
502: /*E
503: PetscBuildTwoSidedType - algorithm for setting up two-sided communication
505: $ PETSC_BUILDTWOSIDED_ALLREDUCE - classical algorithm using an MPI_Allreduce with
506: $ a buffer of length equal to the communicator size. Not memory-scalable due to
507: $ the large reduction size. Requires only MPI-1.
508: $ PETSC_BUILDTWOSIDED_IBARRIER - nonblocking algorithm based on MPI_Issend and MPI_Ibarrier.
509: $ Proved communication-optimal in Hoefler, Siebert, and Lumsdaine (2010). Requires MPI-3.
510: $ PETSC_BUILDTWOSIDED_REDSCATTER - similar to above, but use more optimized function
511: $ that only communicates the part of the reduction that is necessary. Requires MPI-2.
513: Level: developer
515: .seealso: PetscCommBuildTwoSided(), PetscCommBuildTwoSidedSetType(), PetscCommBuildTwoSidedGetType()
516: E*/
517: typedef enum {
518: PETSC_BUILDTWOSIDED_NOTSET = -1,
519: PETSC_BUILDTWOSIDED_ALLREDUCE = 0,
520: PETSC_BUILDTWOSIDED_IBARRIER = 1,
521: PETSC_BUILDTWOSIDED_REDSCATTER = 2
522: /* Updates here must be accompanied by updates in finclude/petscsys.h and the string array in mpits.c */
523: } PetscBuildTwoSidedType;
525: /*E
526: InsertMode - Whether entries are inserted or added into vectors or matrices
528: Level: beginner
530: .seealso: VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(),
531: VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(),
532: MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd()
533: E*/
534: typedef enum {NOT_SET_VALUES, INSERT_VALUES, ADD_VALUES, MAX_VALUES, INSERT_ALL_VALUES, ADD_ALL_VALUES, INSERT_BC_VALUES, ADD_BC_VALUES} InsertMode;
536: /*MC
537: INSERT_VALUES - Put a value into a vector or matrix, overwrites any previous value
539: Level: beginner
541: .seealso: InsertMode, VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(),
542: VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), ADD_VALUES,
543: MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd(), MAX_VALUES
545: M*/
547: /*MC
548: ADD_VALUES - Adds a value into a vector or matrix, if there previously was no value, just puts the
549: value into that location
551: Level: beginner
553: .seealso: InsertMode, VecSetValues(), MatSetValues(), VecSetValue(), VecSetValuesBlocked(),
554: VecSetValuesLocal(), VecSetValuesBlockedLocal(), MatSetValuesBlocked(), INSERT_VALUES,
555: MatSetValuesBlockedLocal(), MatSetValuesLocal(), VecScatterBegin(), VecScatterEnd(), MAX_VALUES
557: M*/
559: /*MC
560: MAX_VALUES - Puts the maximum of the scattered/gathered value and the current value into each location
562: Level: beginner
564: .seealso: InsertMode, VecScatterBegin(), VecScatterEnd(), ADD_VALUES, INSERT_VALUES
566: M*/
568: /*S
569: PetscSubcomm - A decomposition of an MPI communicator into subcommunicators
571: Notes:
572: After a call to PetscSubcommSetType(), PetscSubcommSetTypeGeneral(), or PetscSubcommSetFromOptions() one may call
573: $ PetscSubcommChild() returns the associated subcommunicator on this process
574: $ PetscSubcommContiguousParent() returns a parent communitor but with all child of the same subcommunicator having contiguous rank
576: Sample Usage:
577: PetscSubcommCreate()
578: PetscSubcommSetNumber()
579: PetscSubcommSetType(PETSC_SUBCOMM_INTERLACED);
580: ccomm = PetscSubcommChild()
581: PetscSubcommDestroy()
583: Level: advanced
585: Concepts: communicator, create
587: Notes:
588: $ 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
589: $ PETSC_SUBCOMM_CONTIGUOUS - each new communicator contains a set of process with contiguous ranks in the original MPI communicator
590: $ PETSC_SUBCOMM_INTERLACED - each new communictor contains a set of processes equally far apart in rank from the others in that new communicator
592: Example: Consider a communicator with six processes split into 3 subcommunicators.
593: $ 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
594: $ PETSC_SUBCOMM_INTERLACED - the first communicator contains rank 0,3, the second 1,4 and the third 2,5
596: Developer Notes:
597: This is used in objects such as PCREDUNDANT to manage the subcommunicators on which the redundant computations
598: are performed.
601: .seealso: PetscSubcommCreate(), PetscSubcommSetNumber(), PetscSubcommSetType(), PetscSubcommView(), PetscSubcommSetFromOptions()
603: S*/
604: typedef struct _n_PetscSubcomm* PetscSubcomm;
605: typedef enum {PETSC_SUBCOMM_GENERAL=0,PETSC_SUBCOMM_CONTIGUOUS=1,PETSC_SUBCOMM_INTERLACED=2} PetscSubcommType;
607: /*S
608: PetscHeap - A simple class for managing heaps
610: Level: intermediate
612: Concepts: random numbers
614: .seealso: PetscHeapCreate(), PetscHeapAdd(), PetscHeapPop(), PetscHeapPeek(), PetscHeapStash(), PetscHeapUnstash(), PetscHeapView(), PetscHeapDestroy()
615: S*/
616: typedef struct _PetscHeap *PetscHeap;
618: typedef struct _n_PetscShmComm* PetscShmComm;
619: typedef struct _n_PetscOmpCtrl* PetscOmpCtrl;
621: /*S
622: PetscSegBuffer - a segmented extendable buffer
624: Level: developer
626: .seealso: PetscSegBufferCreate(), PetscSegBufferGet(), PetscSegBufferExtract(), PetscSegBufferDestroy()
627: S*/
628: typedef struct _n_PetscSegBuffer *PetscSegBuffer;
630: typedef struct _n_PetscOptionsHelpPrinted *PetscOptionsHelpPrinted;
632: #endif