Actual source code: pinit.c

  1: #define PETSC_DESIRE_FEATURE_TEST_MACROS
  2: /*
  3:    This file defines the initialization of PETSc, including PetscInitialize()
  4: */
  5: #include <petsc/private/petscimpl.h>
  6: #include <petsc/private/logimpl.h>
  7: #include <petscviewer.h>
  8: #include <petsc/private/garbagecollector.h>

 10: #if !defined(PETSC_HAVE_WINDOWS_COMPILERS)
 11: #include <petsc/private/valgrind/valgrind.h>
 12: #endif

 14: #if defined(PETSC_USE_FORTRAN_BINDINGS)
 15: #include <petsc/private/fortranimpl.h>
 16: #endif

 18: #if PetscDefined(USE_COVERAGE)
 19: EXTERN_C_BEGIN
 20:   #if defined(PETSC_HAVE___GCOV_DUMP)
 22:   #endif
 23: void __gcov_flush(void);
 24: EXTERN_C_END
 25: #endif

 27: #if defined(PETSC_SERIALIZE_FUNCTIONS)
 28: PETSC_INTERN PetscFPT PetscFPTData;
 29: PetscFPT              PetscFPTData = 0;
 30: #endif

 32: #if PetscDefined(HAVE_SAWS)
 33: #include <petscviewersaws.h>
 34: #endif

 36: PETSC_INTERN FILE *petsc_history;

 38: PETSC_INTERN PetscErrorCode PetscInitialize_DynamicLibraries(void);
 39: PETSC_INTERN PetscErrorCode PetscFinalize_DynamicLibraries(void);
 40: PETSC_INTERN PetscErrorCode PetscSequentialPhaseBegin_Private(MPI_Comm, int);
 41: PETSC_INTERN PetscErrorCode PetscSequentialPhaseEnd_Private(MPI_Comm, int);
 42: PETSC_INTERN PetscErrorCode PetscCloseHistoryFile(FILE **);

 44: /* user may set these BEFORE calling PetscInitialize() */
 45: MPI_Comm PETSC_COMM_WORLD = MPI_COMM_NULL;
 46: #if PetscDefined(HAVE_MPI_INIT_THREAD)
 47: PetscMPIInt PETSC_MPI_THREAD_REQUIRED = PETSC_DECIDE;
 48: #else
 49: PetscMPIInt PETSC_MPI_THREAD_REQUIRED = MPI_THREAD_SINGLE;
 50: #endif

 52: PetscMPIInt Petsc_Counter_keyval      = MPI_KEYVAL_INVALID;
 53: PetscMPIInt Petsc_InnerComm_keyval    = MPI_KEYVAL_INVALID;
 54: PetscMPIInt Petsc_OuterComm_keyval    = MPI_KEYVAL_INVALID;
 55: PetscMPIInt Petsc_ShmComm_keyval      = MPI_KEYVAL_INVALID;
 56: PetscMPIInt Petsc_CreationIdx_keyval  = MPI_KEYVAL_INVALID;
 57: PetscMPIInt Petsc_Garbage_HMap_keyval = MPI_KEYVAL_INVALID;

 59: PetscMPIInt Petsc_SharedWD_keyval  = MPI_KEYVAL_INVALID;
 60: PetscMPIInt Petsc_SharedTmp_keyval = MPI_KEYVAL_INVALID;

 62: /*
 63:      Declare and set all the string names of the PETSc enums
 64: */
 65: const char *const PetscBools[]     = {"FALSE", "TRUE", "PetscBool", "PETSC_", NULL};
 66: const char *const PetscCopyModes[] = {"COPY_VALUES", "OWN_POINTER", "USE_POINTER", "PetscCopyMode", "PETSC_", NULL};

 68: PetscBool PetscPreLoadingUsed = PETSC_FALSE;
 69: PetscBool PetscPreLoadingOn   = PETSC_FALSE;

 71: PetscInt PetscHotRegionDepth;

 73: PetscBool PETSC_RUNNING_ON_VALGRIND = PETSC_FALSE;

 75: #if defined(PETSC_HAVE_THREADSAFETY)
 76: PetscSpinlock PetscViewerASCIISpinLockOpen;
 77: PetscSpinlock PetscViewerASCIISpinLockStdout;
 78: PetscSpinlock PetscViewerASCIISpinLockStderr;
 79: PetscSpinlock PetscCommSpinLock;
 80: #endif

 82: /*@C
 83:   PetscInitializeNoPointers - Calls PetscInitialize() from C/C++ without the pointers to argc and args

 85:   Collective

 87:   Input Parameters:
 88: + argc     - number of args
 89: . args     - array of command line arguments
 90: . filename - optional name of the program file, pass `NULL` to ignore
 91: - help     - optional help, pass `NULL` to ignore

 93:   Level: advanced

 95:   Notes:
 96:   this is called only by the PETSc Julia interface. Even though it might start MPI it sets the flag to
 97:   indicate that it did NOT start MPI so that the PetscFinalize() does not end MPI, thus allowing PetscInitialize() to
 98:   be called multiple times from Julia without the problem of trying to initialize MPI more than once.

100:   Developer Notes:
101:   Turns off PETSc signal handling to allow Julia to manage signals

103: .seealso: `PetscInitialize()`, `PetscInitializeFortran()`, `PetscInitializeNoArguments()`
104: */
105: PetscErrorCode PetscInitializeNoPointers(int argc, char **args, const char *filename, const char *help)
106: {
107:   int    myargc = argc;
108:   char **myargs = args;

110:   PetscFunctionBegin;
111:   PetscCall(PetscInitialize(&myargc, &myargs, filename, help));
112:   PetscCall(PetscPopSignalHandler());
113:   PetscBeganMPI = PETSC_FALSE;
114:   PetscFunctionReturn(PETSC_SUCCESS);
115: }

117: /*@C
118:   PetscInitializeNoArguments - Calls `PetscInitialize()` from C/C++ without
119:   the command line arguments.

121:   Collective

123:   Level: advanced

125: .seealso: `PetscInitialize()`, `PetscInitializeFortran()`
126: @*/
127: PetscErrorCode PetscInitializeNoArguments(void)
128: {
129:   int    argc = 0;
130:   char **args = NULL;

132:   PetscFunctionBegin;
133:   PetscCall(PetscInitialize(&argc, &args, NULL, NULL));
134:   PetscFunctionReturn(PETSC_SUCCESS);
135: }

137: /*@
138:   PetscInitialized - Determine whether PETSc is initialized.

140:   Output Parameter:
141: . isInitialized - `PETSC_TRUE` if PETSc is initialized, `PETSC_FALSE` otherwise

143:   Level: beginner

145: .seealso: `PetscInitialize()`, `PetscInitializeNoArguments()`, `PetscInitializeFortran()`
146: @*/
147: PetscErrorCode PetscInitialized(PetscBool *isInitialized)
148: {
149:   PetscFunctionBegin;
150:   if (PetscInitializeCalled) PetscAssertPointer(isInitialized, 1);
151:   *isInitialized = PetscInitializeCalled;
152:   PetscFunctionReturn(PETSC_SUCCESS);
153: }

155: /*@
156:   PetscFinalized - Determine whether `PetscFinalize()` has been called yet

158:   Output Parameter:
159: . isFinalized - `PETSC_TRUE` if PETSc is finalized, `PETSC_FALSE` otherwise

161:   Level: developer

163: .seealso: `PetscInitialize()`, `PetscInitializeNoArguments()`, `PetscInitializeFortran()`
164: @*/
165: PetscErrorCode PetscFinalized(PetscBool *isFinalized)
166: {
167:   PetscFunctionBegin;
168:   if (!PetscFinalizeCalled) PetscAssertPointer(isFinalized, 1);
169:   *isFinalized = PetscFinalizeCalled;
170:   PetscFunctionReturn(PETSC_SUCCESS);
171: }

173: PETSC_INTERN PetscErrorCode PetscOptionsCheckInitial_Private(const char[]);

175: /*
176:        This function is the MPI reduction operation used to compute the sum of the
177:    first half of the datatype and the max of the second half.
178: */
179: MPI_Op MPIU_MAXSUM_OP               = 0;
180: MPI_Op Petsc_Garbage_SetIntersectOp = 0;

182: PETSC_INTERN void MPIAPI MPIU_MaxSum_Local(void *in, void *out, int *cnt, MPI_Datatype *datatype)
183: {
184:   PetscInt *xin = (PetscInt *)in, *xout = (PetscInt *)out, i, count = *cnt;

186:   PetscFunctionBegin;
187:   if (*datatype != MPIU_2INT) {
188:     PetscErrorCode ierr = (*PetscErrorPrintf)("Can only handle MPIU_2INT data types");
189:     (void)ierr;
190:     PETSCABORT(MPI_COMM_SELF, PETSC_ERR_ARG_WRONG);
191:   }

193:   for (i = 0; i < count; i++) {
194:     xout[2 * i] = PetscMax(xout[2 * i], xin[2 * i]);
195:     xout[2 * i + 1] += xin[2 * i + 1];
196:   }
197:   PetscFunctionReturnVoid();
198: }

200: /*
201:     Returns the max of the first entry owned by this processor and the
202: sum of the second entry.

204:     The reason sizes[2*i] contains lengths sizes[2*i+1] contains flag of 1 if length is nonzero
205: is so that the MPIU_MAXSUM_OP() can set TWO values, if we passed in only sizes[i] with lengths
206: there would be no place to store the both needed results.
207: */
208: PetscErrorCode PetscMaxSum(MPI_Comm comm, const PetscInt sizes[], PetscInt *max, PetscInt *sum)
209: {
210:   PetscFunctionBegin;
211: #if defined(PETSC_HAVE_MPI_REDUCE_SCATTER_BLOCK)
212:   {
213:     struct {
214:       PetscInt max, sum;
215:     } work;
216:     PetscCallMPI(MPI_Reduce_scatter_block((void *)sizes, &work, 1, MPIU_2INT, MPIU_MAXSUM_OP, comm));
217:     *max = work.max;
218:     *sum = work.sum;
219:   }
220: #else
221:   {
222:     PetscMPIInt size, rank;
223:     struct {
224:       PetscInt max, sum;
225:     } *work;
226:     PetscCallMPI(MPI_Comm_size(comm, &size));
227:     PetscCallMPI(MPI_Comm_rank(comm, &rank));
228:     PetscCall(PetscMalloc1(size, &work));
229:     PetscCall(MPIU_Allreduce((void *)sizes, work, size, MPIU_2INT, MPIU_MAXSUM_OP, comm));
230:     *max = work[rank].max;
231:     *sum = work[rank].sum;
232:     PetscCall(PetscFree(work));
233:   }
234: #endif
235:   PetscFunctionReturn(PETSC_SUCCESS);
236: }

238: #if defined(PETSC_HAVE_REAL___FLOAT128) || defined(PETSC_HAVE_REAL___FP16)
239:   #if defined(PETSC_HAVE_REAL___FLOAT128)
240:     #include <quadmath.h>
241:   #endif
242: MPI_Op MPIU_SUM___FP16___FLOAT128 = 0;
243:   #if defined(PETSC_USE_REAL___FLOAT128) || defined(PETSC_USE_REAL___FP16)
244: MPI_Op MPIU_SUM = 0;
245:   #endif

247: PETSC_EXTERN void MPIAPI PetscSum_Local(void *in, void *out, PetscMPIInt *cnt, MPI_Datatype *datatype)
248: {
249:   PetscInt i, count = *cnt;

251:   PetscFunctionBegin;
252:   if (*datatype == MPIU_REAL) {
253:     PetscReal *xin = (PetscReal *)in, *xout = (PetscReal *)out;
254:     for (i = 0; i < count; i++) xout[i] += xin[i];
255:   }
256:   #if defined(PETSC_HAVE_COMPLEX)
257:   else if (*datatype == MPIU_COMPLEX) {
258:     PetscComplex *xin = (PetscComplex *)in, *xout = (PetscComplex *)out;
259:     for (i = 0; i < count; i++) xout[i] += xin[i];
260:   }
261:   #endif
262:   #if defined(PETSC_HAVE_REAL___FLOAT128)
263:   else if (*datatype == MPIU___FLOAT128) {
264:     __float128 *xin = (__float128 *)in, *xout = (__float128 *)out;
265:     for (i = 0; i < count; i++) xout[i] += xin[i];
266:     #if defined(PETSC_HAVE_COMPLEX)
267:   } else if (*datatype == MPIU___COMPLEX128) {
268:     __complex128 *xin = (__complex128 *)in, *xout = (__complex128 *)out;
269:     for (i = 0; i < count; i++) xout[i] += xin[i];
270:     #endif
271:   }
272:   #endif
273:   #if defined(PETSC_HAVE_REAL___FP16)
274:   else if (*datatype == MPIU___FP16) {
275:     __fp16 *xin = (__fp16 *)in, *xout = (__fp16 *)out;
276:     for (i = 0; i < count; i++) xout[i] += xin[i];
277:   }
278:   #endif
279:   else {
280:   #if !defined(PETSC_HAVE_REAL___FLOAT128) && !defined(PETSC_HAVE_REAL___FP16)
281:     PetscCallAbort(MPI_COMM_SElF, (*PetscErrorPrintf)("Can only handle MPIU_REAL or MPIU_COMPLEX data types"));
282:   #elif !defined(PETSC_HAVE_REAL___FP16)
283:     PetscCallAbort(MPI_COMM_SELF, (*PetscErrorPrintf)("Can only handle MPIU_REAL, MPIU_COMPLEX, MPIU___FLOAT128, or MPIU___COMPLEX128 data types"));
284:   #elif !defined(PETSC_HAVE_REAL___FLOAT128)
285:     PetscCallAbort(MPI_COMM_SELF, (*PetscErrorPrintf)("Can only handle MPIU_REAL, MPIU_COMPLEX, or MPIU___FP16 data types"));
286:   #else
287:     PetscCallAbort(MPI_COMM_SELF, (*PetscErrorPrintf)("Can only handle MPIU_REAL, MPIU_COMPLEX, MPIU___FLOAT128, MPIU___COMPLEX128, or MPIU___FP16 data types"));
288:   #endif
289:     PETSCABORT(MPI_COMM_SELF, PETSC_ERR_ARG_WRONG);
290:   }
291:   PetscFunctionReturnVoid();
292: }
293: #endif

295: #if defined(PETSC_USE_REAL___FLOAT128) || defined(PETSC_USE_REAL___FP16)
296: MPI_Op MPIU_MAX = 0;
297: MPI_Op MPIU_MIN = 0;

299: PETSC_EXTERN void MPIAPI PetscMax_Local(void *in, void *out, PetscMPIInt *cnt, MPI_Datatype *datatype)
300: {
301:   PetscInt i, count = *cnt;

303:   PetscFunctionBegin;
304:   if (*datatype == MPIU_REAL) {
305:     PetscReal *xin = (PetscReal *)in, *xout = (PetscReal *)out;
306:     for (i = 0; i < count; i++) xout[i] = PetscMax(xout[i], xin[i]);
307:   }
308:   #if defined(PETSC_HAVE_COMPLEX)
309:   else if (*datatype == MPIU_COMPLEX) {
310:     PetscComplex *xin = (PetscComplex *)in, *xout = (PetscComplex *)out;
311:     for (i = 0; i < count; i++) xout[i] = PetscRealPartComplex(xout[i]) < PetscRealPartComplex(xin[i]) ? xin[i] : xout[i];
312:   }
313:   #endif
314:   else {
315:     PetscCallAbort(MPI_COMM_SELF, (*PetscErrorPrintf)("Can only handle MPIU_REAL or MPIU_COMPLEX data types"));
316:     PETSCABORT(MPI_COMM_SELF, PETSC_ERR_ARG_WRONG);
317:   }
318:   PetscFunctionReturnVoid();
319: }

321: PETSC_EXTERN void MPIAPI PetscMin_Local(void *in, void *out, PetscMPIInt *cnt, MPI_Datatype *datatype)
322: {
323:   PetscInt i, count = *cnt;

325:   PetscFunctionBegin;
326:   if (*datatype == MPIU_REAL) {
327:     PetscReal *xin = (PetscReal *)in, *xout = (PetscReal *)out;
328:     for (i = 0; i < count; i++) xout[i] = PetscMin(xout[i], xin[i]);
329:   }
330:   #if defined(PETSC_HAVE_COMPLEX)
331:   else if (*datatype == MPIU_COMPLEX) {
332:     PetscComplex *xin = (PetscComplex *)in, *xout = (PetscComplex *)out;
333:     for (i = 0; i < count; i++) xout[i] = PetscRealPartComplex(xout[i]) > PetscRealPartComplex(xin[i]) ? xin[i] : xout[i];
334:   }
335:   #endif
336:   else {
337:     PetscCallAbort(MPI_COMM_SELF, (*PetscErrorPrintf)("Can only handle MPIU_REAL or MPIU_SCALAR data (i.e. double or complex) types"));
338:     PETSCABORT(MPI_COMM_SELF, PETSC_ERR_ARG_WRONG);
339:   }
340:   PetscFunctionReturnVoid();
341: }
342: #endif

344: /*
345:    Private routine to delete internal tag/name counter storage when a communicator is freed.

347:    This is called by MPI, not by users. This is called by MPI_Comm_free() when the communicator that has this  data as an attribute is freed.

349:    Note: this is declared extern "C" because it is passed to MPI_Comm_create_keyval()

351: */
352: PETSC_EXTERN PetscMPIInt MPIAPI Petsc_Counter_Attr_Delete_Fn(MPI_Comm comm, PetscMPIInt keyval, void *count_val, void *extra_state)
353: {
354:   PetscCommCounter      *counter = (PetscCommCounter *)count_val;
355:   struct PetscCommStash *comms   = counter->comms, *pcomm;

357:   PetscFunctionBegin;
358:   PetscCallMPI(PetscInfo(NULL, "Deleting counter data in an MPI_Comm %ld\n", (long)comm));
359:   PetscCallMPI(PetscFree(counter->iflags));
360:   while (comms) {
361:     PetscCallMPI(MPI_Comm_free(&comms->comm));
362:     pcomm = comms;
363:     comms = comms->next;
364:     PetscCall(PetscFree(pcomm));
365:   }
366:   PetscCallMPI(PetscFree(counter));
367:   PetscFunctionReturn(MPI_SUCCESS);
368: }

370: /*
371:   This is invoked on the outer comm as a result of either PetscCommDestroy() (via MPI_Comm_delete_attr) or when the user
372:   calls MPI_Comm_free().

374:   This is the only entry point for breaking the links between inner and outer comms.

376:   This is called by MPI, not by users. This is called when MPI_Comm_free() is called on the communicator.

378:   Note: this is declared extern "C" because it is passed to MPI_Comm_create_keyval()

380: */
381: PETSC_EXTERN PetscMPIInt MPIAPI Petsc_InnerComm_Attr_Delete_Fn(MPI_Comm comm, PetscMPIInt keyval, void *attr_val, void *extra_state)
382: {
383:   union
384:   {
385:     MPI_Comm comm;
386:     void    *ptr;
387:   } icomm;

389:   PetscFunctionBegin;
390:   if (keyval != Petsc_InnerComm_keyval) SETERRMPI(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "Unexpected keyval");
391:   icomm.ptr = attr_val;
392:   if (PetscDefined(USE_DEBUG)) {
393:     /* Error out if the inner/outer comms are not correctly linked through their Outer/InnterComm attributes */
394:     PetscMPIInt flg;
395:     union
396:     {
397:       MPI_Comm comm;
398:       void    *ptr;
399:     } ocomm;
400:     PetscCallMPI(MPI_Comm_get_attr(icomm.comm, Petsc_OuterComm_keyval, &ocomm, &flg));
401:     if (!flg) SETERRMPI(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "Inner comm does not have OuterComm attribute");
402:     if (ocomm.comm != comm) SETERRMPI(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "Inner comm's OuterComm attribute does not point to outer PETSc comm");
403:   }
404:   PetscCallMPI(MPI_Comm_delete_attr(icomm.comm, Petsc_OuterComm_keyval));
405:   PetscCallMPI(PetscInfo(NULL, "User MPI_Comm %ld is being unlinked from inner PETSc comm %ld\n", (long)comm, (long)icomm.comm));
406:   PetscFunctionReturn(MPI_SUCCESS);
407: }

409: /*
410:  * This is invoked on the inner comm when Petsc_InnerComm_Attr_Delete_Fn calls MPI_Comm_delete_attr().  It should not be reached any other way.
411:  */
412: PETSC_EXTERN PetscMPIInt MPIAPI Petsc_OuterComm_Attr_Delete_Fn(MPI_Comm comm, PetscMPIInt keyval, void *attr_val, void *extra_state)
413: {
414:   PetscFunctionBegin;
415:   PetscCallMPI(PetscInfo(NULL, "Removing reference to PETSc communicator embedded in a user MPI_Comm %ld\n", (long)comm));
416:   PetscFunctionReturn(MPI_SUCCESS);
417: }

419: PETSC_EXTERN PetscMPIInt MPIAPI Petsc_ShmComm_Attr_Delete_Fn(MPI_Comm, PetscMPIInt, void *, void *);

421: #if defined(PETSC_USE_PETSC_MPI_EXTERNAL32)
422: PETSC_EXTERN PetscMPIInt PetscDataRep_extent_fn(MPI_Datatype, MPI_Aint *, void *);
423: PETSC_EXTERN PetscMPIInt PetscDataRep_read_conv_fn(void *, MPI_Datatype, PetscMPIInt, void *, MPI_Offset, void *);
424: PETSC_EXTERN PetscMPIInt PetscDataRep_write_conv_fn(void *, MPI_Datatype, PetscMPIInt, void *, MPI_Offset, void *);
425: #endif

427: PetscMPIInt PETSC_MPI_ERROR_CLASS = MPI_ERR_LASTCODE, PETSC_MPI_ERROR_CODE;

429: PETSC_INTERN int    PetscGlobalArgc;
430: PETSC_INTERN char **PetscGlobalArgs;
431: int                 PetscGlobalArgc = 0;
432: char              **PetscGlobalArgs = NULL;
433: PetscSegBuffer      PetscCitationsList;

435: PetscErrorCode PetscCitationsInitialize(void)
436: {
437:   PetscFunctionBegin;
438:   PetscCall(PetscSegBufferCreate(1, 10000, &PetscCitationsList));

440:   PetscCall(PetscCitationsRegister("@TechReport{petsc-user-ref,\n\
441:   Author = {Satish Balay and Shrirang Abhyankar and Mark~F. Adams and Steven Benson and Jed Brown\n\
442:     and Peter Brune and Kris Buschelman and Emil Constantinescu and Lisandro Dalcin and Alp Dener\n\
443:     and Victor Eijkhout and Jacob Faibussowitsch and William~D. Gropp and V\'{a}clav Hapla and Tobin Isaac and Pierre Jolivet\n\
444:     and Dmitry Karpeev and Dinesh Kaushik and Matthew~G. Knepley and Fande Kong and Scott Kruger\n\
445:     and Dave~A. May and Lois Curfman McInnes and Richard Tran Mills and Lawrence Mitchell and Todd Munson\n\
446:     and Jose~E. Roman and Karl Rupp and Patrick Sanan and Jason Sarich and Barry~F. Smith\n\
447:     and Stefano Zampini and Hong Zhang and Hong Zhang and Junchao Zhang},\n\
448:   Title = {{PETSc/TAO} Users Manual},\n\
449:   Number = {ANL-21/39 - Revision 3.20},\n\
450:   Doi = {10.2172/2205494},\n\
451:   Institution = {Argonne National Laboratory},\n\
452:   Year = {2023}\n}\n",
453:                                    NULL));

455:   PetscCall(PetscCitationsRegister("@InProceedings{petsc-efficient,\n\
456:   Author = {Satish Balay and William D. Gropp and Lois Curfman McInnes and Barry F. Smith},\n\
457:   Title = {Efficient Management of Parallelism in Object Oriented Numerical Software Libraries},\n\
458:   Booktitle = {Modern Software Tools in Scientific Computing},\n\
459:   Editor = {E. Arge and A. M. Bruaset and H. P. Langtangen},\n\
460:   Pages = {163--202},\n\
461:   Publisher = {Birkh{\\\"{a}}user Press},\n\
462:   Year = {1997}\n}\n",
463:                                    NULL));

465:   PetscFunctionReturn(PETSC_SUCCESS);
466: }

468: static char programname[PETSC_MAX_PATH_LEN] = ""; /* HP includes entire path in name */

470: PetscErrorCode PetscSetProgramName(const char name[])
471: {
472:   PetscFunctionBegin;
473:   PetscCall(PetscStrncpy(programname, name, sizeof(programname)));
474:   PetscFunctionReturn(PETSC_SUCCESS);
475: }

477: /*@C
478:   PetscGetProgramName - Gets the name of the running program.

480:   Not Collective

482:   Input Parameter:
483: . len - length of the string name

485:   Output Parameter:
486: . name - the name of the running program, provide a string of length `PETSC_MAX_PATH_LEN`

488:   Level: advanced

490: .seealso: `PetscFinalize()`, `PetscInitializeFortran()`, `PetscGetArguments()`, `PetscInitialize()`
491: @*/
492: PetscErrorCode PetscGetProgramName(char name[], size_t len)
493: {
494:   PetscFunctionBegin;
495:   PetscCall(PetscStrncpy(name, programname, len));
496:   PetscFunctionReturn(PETSC_SUCCESS);
497: }

499: /*@C
500:   PetscGetArgs - Allows you to access the raw command line arguments anywhere
501:   after PetscInitialize() is called but before `PetscFinalize()`.

503:   Not Collective

505:   Output Parameters:
506: + argc - count of number of command line arguments
507: - args - the command line arguments

509:   Level: intermediate

511:   Notes:
512:   This is usually used to pass the command line arguments into other libraries
513:   that are called internally deep in PETSc or the application.

515:   The first argument contains the program name as is normal for C arguments.

517: .seealso: `PetscFinalize()`, `PetscInitializeFortran()`, `PetscGetArguments()`, `PetscInitialize()`
518: @*/
519: PetscErrorCode PetscGetArgs(int *argc, char ***args)
520: {
521:   PetscFunctionBegin;
522:   PetscCheck(PetscInitializeCalled || !PetscFinalizeCalled, PETSC_COMM_SELF, PETSC_ERR_ORDER, "You must call after PetscInitialize() but before PetscFinalize()");
523:   *argc = PetscGlobalArgc;
524:   *args = PetscGlobalArgs;
525:   PetscFunctionReturn(PETSC_SUCCESS);
526: }

528: /*@C
529:   PetscGetArguments - Allows you to access the  command line arguments anywhere
530:   after `PetscInitialize()` is called but before `PetscFinalize()`.

532:   Not Collective

534:   Output Parameter:
535: . args - the command line arguments

537:   Level: intermediate

539:   Note:
540:   This does NOT start with the program name and IS `NULL` terminated (final arg is void)

542: .seealso: `PetscFinalize()`, `PetscInitializeFortran()`, `PetscGetArgs()`, `PetscFreeArguments()`, `PetscInitialize()`
543: @*/
544: PetscErrorCode PetscGetArguments(char ***args)
545: {
546:   PetscInt i, argc = PetscGlobalArgc;

548:   PetscFunctionBegin;
549:   PetscCheck(PetscInitializeCalled || !PetscFinalizeCalled, PETSC_COMM_SELF, PETSC_ERR_ORDER, "You must call after PetscInitialize() but before PetscFinalize()");
550:   if (!argc) {
551:     *args = NULL;
552:     PetscFunctionReturn(PETSC_SUCCESS);
553:   }
554:   PetscCall(PetscMalloc1(argc, args));
555:   for (i = 0; i < argc - 1; i++) PetscCall(PetscStrallocpy(PetscGlobalArgs[i + 1], &(*args)[i]));
556:   (*args)[argc - 1] = NULL;
557:   PetscFunctionReturn(PETSC_SUCCESS);
558: }

560: /*@C
561:   PetscFreeArguments - Frees the memory obtained with `PetscGetArguments()`

563:   Not Collective

565:   Output Parameter:
566: . args - the command line arguments

568:   Level: intermediate

570: .seealso: `PetscFinalize()`, `PetscInitializeFortran()`, `PetscGetArgs()`, `PetscGetArguments()`
571: @*/
572: PetscErrorCode PetscFreeArguments(char **args)
573: {
574:   PetscFunctionBegin;
575:   if (args) {
576:     PetscInt i = 0;

578:     while (args[i]) PetscCall(PetscFree(args[i++]));
579:     PetscCall(PetscFree(args));
580:   }
581:   PetscFunctionReturn(PETSC_SUCCESS);
582: }

584: #if PetscDefined(HAVE_SAWS)
585:   #include <petscconfiginfo.h>

587: PETSC_INTERN PetscErrorCode PetscInitializeSAWs(const char help[])
588: {
589:   PetscFunctionBegin;
590:   if (!PetscGlobalRank) {
591:     char      cert[PETSC_MAX_PATH_LEN], root[PETSC_MAX_PATH_LEN], *intro, programname[64], *appline, *options, version[64];
592:     int       port;
593:     PetscBool flg, rootlocal = PETSC_FALSE, flg2, selectport = PETSC_FALSE;
594:     size_t    applinelen, introlen;
595:     char      sawsurl[256];

597:     PetscCall(PetscOptionsHasName(NULL, NULL, "-saws_log", &flg));
598:     if (flg) {
599:       char sawslog[PETSC_MAX_PATH_LEN];

601:       PetscCall(PetscOptionsGetString(NULL, NULL, "-saws_log", sawslog, sizeof(sawslog), NULL));
602:       if (sawslog[0]) {
603:         PetscCallSAWs(SAWs_Set_Use_Logfile, (sawslog));
604:       } else {
605:         PetscCallSAWs(SAWs_Set_Use_Logfile, (NULL));
606:       }
607:     }
608:     PetscCall(PetscOptionsGetString(NULL, NULL, "-saws_https", cert, sizeof(cert), &flg));
609:     if (flg) PetscCallSAWs(SAWs_Set_Use_HTTPS, (cert));
610:     PetscCall(PetscOptionsGetBool(NULL, NULL, "-saws_port_auto_select", &selectport, NULL));
611:     if (selectport) {
612:       PetscCallSAWs(SAWs_Get_Available_Port, (&port));
613:       PetscCallSAWs(SAWs_Set_Port, (port));
614:     } else {
615:       PetscCall(PetscOptionsGetInt(NULL, NULL, "-saws_port", &port, &flg));
616:       if (flg) PetscCallSAWs(SAWs_Set_Port, (port));
617:     }
618:     PetscCall(PetscOptionsGetString(NULL, NULL, "-saws_root", root, sizeof(root), &flg));
619:     if (flg) {
620:       PetscCallSAWs(SAWs_Set_Document_Root, (root));
621:       PetscCall(PetscStrcmp(root, ".", &rootlocal));
622:     } else {
623:       PetscCall(PetscOptionsHasName(NULL, NULL, "-saws_options", &flg));
624:       if (flg) {
625:         PetscCall(PetscStrreplace(PETSC_COMM_WORLD, "${PETSC_DIR}/share/petsc/saws", root, sizeof(root)));
626:         PetscCallSAWs(SAWs_Set_Document_Root, (root));
627:       }
628:     }
629:     PetscCall(PetscOptionsHasName(NULL, NULL, "-saws_local", &flg2));
630:     if (flg2) {
631:       char jsdir[PETSC_MAX_PATH_LEN];
632:       PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_SUP, "-saws_local option requires -saws_root option");
633:       PetscCall(PetscSNPrintf(jsdir, sizeof(jsdir), "%s/js", root));
634:       PetscCall(PetscTestDirectory(jsdir, 'r', &flg));
635:       PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_FILE_READ, "-saws_local option requires js directory in root directory");
636:       PetscCallSAWs(SAWs_Push_Local_Header, ());
637:     }
638:     PetscCall(PetscGetProgramName(programname, sizeof(programname)));
639:     PetscCall(PetscStrlen(help, &applinelen));
640:     introlen = 4096 + applinelen;
641:     applinelen += 1024;
642:     PetscCall(PetscMalloc(applinelen, &appline));
643:     PetscCall(PetscMalloc(introlen, &intro));

645:     if (rootlocal) {
646:       PetscCall(PetscSNPrintf(appline, applinelen, "%s.c.html", programname));
647:       PetscCall(PetscTestFile(appline, 'r', &rootlocal));
648:     }
649:     PetscCall(PetscOptionsGetAll(NULL, &options));
650:     if (rootlocal && help) {
651:       PetscCall(PetscSNPrintf(appline, applinelen, "<center> Running <a href=\"%s.c.html\">%s</a> %s</center><br><center><pre>%s</pre></center><br>\n", programname, programname, options, help));
652:     } else if (help) {
653:       PetscCall(PetscSNPrintf(appline, applinelen, "<center>Running %s %s</center><br><center><pre>%s</pre></center><br>", programname, options, help));
654:     } else {
655:       PetscCall(PetscSNPrintf(appline, applinelen, "<center> Running %s %s</center><br>\n", programname, options));
656:     }
657:     PetscCall(PetscFree(options));
658:     PetscCall(PetscGetVersion(version, sizeof(version)));
659:     PetscCall(PetscSNPrintf(intro, introlen,
660:                             "<body>\n"
661:                             "<center><h2> <a href=\"https://petsc.org/\">PETSc</a> Application Web server powered by <a href=\"https://bitbucket.org/saws/saws\">SAWs</a> </h2></center>\n"
662:                             "<center>This is the default PETSc application dashboard, from it you can access any published PETSc objects or logging data</center><br><center>%s configured with %s</center><br>\n"
663:                             "%s",
664:                             version, petscconfigureoptions, appline));
665:     PetscCallSAWs(SAWs_Push_Body, ("index.html", 0, intro));
666:     PetscCall(PetscFree(intro));
667:     PetscCall(PetscFree(appline));
668:     if (selectport) {
669:       PetscBool silent;

671:       /* another process may have grabbed the port so keep trying */
672:       while (SAWs_Initialize()) {
673:         PetscCallSAWs(SAWs_Get_Available_Port, (&port));
674:         PetscCallSAWs(SAWs_Set_Port, (port));
675:       }

677:       PetscCall(PetscOptionsGetBool(NULL, NULL, "-saws_port_auto_select_silent", &silent, NULL));
678:       if (!silent) {
679:         PetscCallSAWs(SAWs_Get_FullURL, (sizeof(sawsurl), sawsurl));
680:         PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Point your browser to %s for SAWs\n", sawsurl));
681:       }
682:     } else {
683:       PetscCallSAWs(SAWs_Initialize, ());
684:     }
685:     PetscCall(PetscCitationsRegister("@TechReport{ saws,\n"
686:                                      "  Author = {Matt Otten and Jed Brown and Barry Smith},\n"
687:                                      "  Title  = {Scientific Application Web Server (SAWs) Users Manual},\n"
688:                                      "  Institution = {Argonne National Laboratory},\n"
689:                                      "  Year   = 2013\n}\n",
690:                                      NULL));
691:   }
692:   PetscFunctionReturn(PETSC_SUCCESS);
693: }
694: #endif

696: /* Things must be done before MPI_Init() when MPI is not yet initialized, and can be shared between C init and Fortran init */
697: PETSC_INTERN PetscErrorCode PetscPreMPIInit_Private(void)
698: {
699:   PetscFunctionBegin;
700: #if defined(PETSC_HAVE_HWLOC_SOLARIS_BUG)
701:   /* see MPI.py for details on this bug */
702:   (void)setenv("HWLOC_COMPONENTS", "-x86", 1);
703: #endif
704:   PetscFunctionReturn(PETSC_SUCCESS);
705: }

707: #if PetscDefined(HAVE_ADIOS)
708:   #include <adios.h>
709:   #include <adios_read.h>
710: int64_t Petsc_adios_group;
711: #endif
712: #if PetscDefined(HAVE_OPENMP)
713:   #include <omp.h>
714: PetscInt PetscNumOMPThreads;
715: #endif

717: #include <petsc/private/deviceimpl.h>
718: #if PetscDefined(HAVE_CUDA)
719: #include <petscdevice_cuda.h>
720: // REMOVE ME
721: cudaStream_t PetscDefaultCudaStream = NULL;
722: #endif
723: #if PetscDefined(HAVE_HIP)
724: #include <petscdevice_hip.h>
725: // REMOVE ME
726: hipStream_t PetscDefaultHipStream = NULL;
727: #endif

729: #if PetscDefined(HAVE_DLFCN_H)
730:   #include <dlfcn.h>
731: #endif
732: PETSC_INTERN PetscErrorCode PetscLogInitialize(void);
733: #if PetscDefined(HAVE_VIENNACL)
734: PETSC_EXTERN PetscErrorCode PetscViennaCLInit(void);
735: PetscBool                   PetscViennaCLSynchronize = PETSC_FALSE;
736: #endif

738: PetscBool PetscCIEnabled = PETSC_FALSE, PetscCIEnabledPortableErrorOutput = PETSC_FALSE;

740: /*
741:   PetscInitialize_Common  - shared code between C and Fortran initialization

743:   prog:     program name
744:   file:     optional PETSc database file name. Might be in Fortran string format when 'ftn' is true
745:   help:     program help message
746:   ftn:      is it called from Fortran initialization (petscinitializef_)?
747:   readarguments,len: used when fortran is true
748: */
749: PETSC_INTERN PetscErrorCode PetscInitialize_Common(const char *prog, const char *file, const char *help, PetscBool ftn, PetscBool readarguments, PetscInt len)
750: {
751:   PetscMPIInt size;
752:   PetscBool   flg = PETSC_TRUE;
753:   char        hostname[256];

755:   PetscFunctionBegin;
756:   if (PetscInitializeCalled) PetscFunctionReturn(PETSC_SUCCESS);
757:   /* these must be initialized in a routine, not as a constant declaration */
758:   PETSC_STDOUT = stdout;
759:   PETSC_STDERR = stderr;

761:   /* PetscCall can be used from now */
762:   PetscErrorHandlingInitialized = PETSC_TRUE;

764:   /*
765:       The checking over compatible runtime libraries is complicated by the MPI ABI initiative
766:       https://wiki.mpich.org/mpich/index.php/ABI_Compatibility_Initiative which started with
767:         MPICH v3.1 (Released February 2014)
768:         IBM MPI v2.1 (December 2014)
769:         Intel MPI Library v5.0 (2014)
770:         Cray MPT v7.0.0 (June 2014)
771:       As of July 31, 2017 the ABI number still appears to be 12, that is all of the versions
772:       listed above and since that time are compatible.

774:       Unfortunately the MPI ABI initiative has not defined a way to determine the ABI number
775:       at compile time or runtime. Thus we will need to systematically track the allowed versions
776:       and how they are represented in the mpi.h and MPI_Get_library_version() output in order
777:       to perform the checking.

779:       Currently we only check for pre MPI ABI versions (and packages that do not follow the MPI ABI).

781:       Questions:

783:         Should the checks for ABI incompatibility be only on the major version number below?
784:         Presumably the output to stderr will be removed before a release.
785:   */

787: #if defined(PETSC_HAVE_MPI_GET_LIBRARY_VERSION)
788:   {
789:     char        mpilibraryversion[MPI_MAX_LIBRARY_VERSION_STRING];
790:     PetscMPIInt mpilibraryversionlength;

792:     PetscCallMPI(MPI_Get_library_version(mpilibraryversion, &mpilibraryversionlength));
793:     /* check for MPICH versions before MPI ABI initiative */
794:   #if defined(MPICH_VERSION)
795:     #if MPICH_NUMVERSION < 30100000
796:     {
797:       char     *ver, *lf;
798:       PetscBool flg = PETSC_FALSE;

800:       PetscCall(PetscStrstr(mpilibraryversion, "MPICH Version:", &ver));
801:       if (ver) {
802:         PetscCall(PetscStrchr(ver, '\n', &lf));
803:         if (lf) {
804:           *lf = 0;
805:           PetscCall(PetscStrendswith(ver, MPICH_VERSION, &flg));
806:         }
807:       }
808:       if (!flg) {
809:         PetscCall(PetscInfo(NULL, "PETSc warning --- MPICH library version \n%s does not match what PETSc was compiled with %s.\n", mpilibraryversion, MPICH_VERSION));
810:         flg = PETSC_TRUE;
811:       }
812:     }
813:     #endif
814:       /* check for Open MPI version, it is not part of the MPI ABI initiative (is it part of another initiative that needs to be handled?) */
815:   #elif defined(OMPI_MAJOR_VERSION)
816:     {
817:       char     *ver, bs[MPI_MAX_LIBRARY_VERSION_STRING], *bsf;
818:       PetscBool flg                                              = PETSC_FALSE;
819:     #define PSTRSZ 2
820:       char      ompistr1[PSTRSZ][MPI_MAX_LIBRARY_VERSION_STRING] = {"Open MPI", "FUJITSU MPI"};
821:       char      ompistr2[PSTRSZ][MPI_MAX_LIBRARY_VERSION_STRING] = {"v", "Library "};
822:       int       i;
823:       for (i = 0; i < PSTRSZ; i++) {
824:         PetscCall(PetscStrstr(mpilibraryversion, ompistr1[i], &ver));
825:         if (ver) {
826:           PetscCall(PetscSNPrintf(bs, MPI_MAX_LIBRARY_VERSION_STRING, "%s%d.%d", ompistr2[i], OMPI_MAJOR_VERSION, OMPI_MINOR_VERSION));
827:           PetscCall(PetscStrstr(ver, bs, &bsf));
828:           if (bsf) flg = PETSC_TRUE;
829:           break;
830:         }
831:       }
832:       if (!flg) {
833:         PetscCall(PetscInfo(NULL, "PETSc warning --- Open MPI library version \n%s does not match what PETSc was compiled with %d.%d.\n", mpilibraryversion, OMPI_MAJOR_VERSION, OMPI_MINOR_VERSION));
834:         flg = PETSC_TRUE;
835:       }
836:     }
837:   #endif
838:   }
839: #endif

841: #if defined(PETSC_HAVE_DLADDR) && !(defined(__cray__) && defined(__clang__))
842:   /* These symbols are currently in the Open MPI and MPICH libraries; they may not always be, in that case the test will simply not detect the problem */
843:   PetscCheck(!dlsym(RTLD_DEFAULT, "ompi_mpi_init") || !dlsym(RTLD_DEFAULT, "MPID_Abort"), PETSC_COMM_SELF, PETSC_ERR_MPI_LIB_INCOMP, "Application was linked against both Open MPI and MPICH based MPI libraries and will not run correctly");
844: #endif

846:   /* on Windows - set printf to default to printing 2 digit exponents */
847: #if defined(PETSC_HAVE__SET_OUTPUT_FORMAT)
848:   _set_output_format(_TWO_DIGIT_EXPONENT);
849: #endif

851:   PetscCall(PetscOptionsCreateDefault());

853:   PetscFinalizeCalled = PETSC_FALSE;

855:   PetscCall(PetscSetProgramName(prog));
856:   PetscCall(PetscSpinlockCreate(&PetscViewerASCIISpinLockOpen));
857:   PetscCall(PetscSpinlockCreate(&PetscViewerASCIISpinLockStdout));
858:   PetscCall(PetscSpinlockCreate(&PetscViewerASCIISpinLockStderr));
859:   PetscCall(PetscSpinlockCreate(&PetscCommSpinLock));

861:   if (PETSC_COMM_WORLD == MPI_COMM_NULL) PETSC_COMM_WORLD = MPI_COMM_WORLD;
862:   PetscCallMPI(MPI_Comm_set_errhandler(PETSC_COMM_WORLD, MPI_ERRORS_RETURN));

864:   if (PETSC_MPI_ERROR_CLASS == MPI_ERR_LASTCODE) {
865:     PetscCallMPI(MPI_Add_error_class(&PETSC_MPI_ERROR_CLASS));
866:     PetscCallMPI(MPI_Add_error_code(PETSC_MPI_ERROR_CLASS, &PETSC_MPI_ERROR_CODE));
867:   }

869:   /* Done after init due to a bug in MPICH-GM? */
870:   PetscCall(PetscErrorPrintfInitialize());

872:   PetscCallMPI(MPI_Comm_rank(MPI_COMM_WORLD, &PetscGlobalRank));
873:   PetscCallMPI(MPI_Comm_size(MPI_COMM_WORLD, &PetscGlobalSize));

875:   MPIU_BOOL        = MPI_INT;
876:   MPIU_ENUM        = MPI_INT;
877:   MPIU_FORTRANADDR = (sizeof(void *) == sizeof(int)) ? MPI_INT : MPIU_INT64;
878:   if (sizeof(size_t) == sizeof(unsigned)) MPIU_SIZE_T = MPI_UNSIGNED;
879:   else if (sizeof(size_t) == sizeof(unsigned long)) MPIU_SIZE_T = MPI_UNSIGNED_LONG;
880: #if defined(PETSC_SIZEOF_LONG_LONG)
881:   else if (sizeof(size_t) == sizeof(unsigned long long)) MPIU_SIZE_T = MPI_UNSIGNED_LONG_LONG;
882: #endif
883:   else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "Could not find MPI type for size_t");

885:     /*
886:      Initialized the global complex variable; this is because with
887:      shared libraries the constructors for global variables
888:      are not called; at least on IRIX.
889:   */
890: #if defined(PETSC_HAVE_COMPLEX)
891:   {
892:   #if defined(PETSC_CLANGUAGE_CXX) && !defined(PETSC_USE_REAL___FLOAT128)
893:     PetscComplex ic(0.0, 1.0);
894:     PETSC_i = ic;
895:   #else
896:     PETSC_i = _Complex_I;
897:   #endif
898:   }
899: #endif /* PETSC_HAVE_COMPLEX */

901:   /*
902:      Create the PETSc MPI reduction operator that sums of the first
903:      half of the entries and maxes the second half.
904:   */
905:   PetscCallMPI(MPI_Op_create(MPIU_MaxSum_Local, 1, &MPIU_MAXSUM_OP));

907: #if defined(PETSC_HAVE_REAL___FLOAT128)
908:   PetscCallMPI(MPI_Type_contiguous(2, MPI_DOUBLE, &MPIU___FLOAT128));
909:   PetscCallMPI(MPI_Type_commit(&MPIU___FLOAT128));
910:   #if defined(PETSC_HAVE_COMPLEX)
911:   PetscCallMPI(MPI_Type_contiguous(4, MPI_DOUBLE, &MPIU___COMPLEX128));
912:   PetscCallMPI(MPI_Type_commit(&MPIU___COMPLEX128));
913:   #endif
914: #endif
915: #if defined(PETSC_HAVE_REAL___FP16)
916:   PetscCallMPI(MPI_Type_contiguous(2, MPI_CHAR, &MPIU___FP16));
917:   PetscCallMPI(MPI_Type_commit(&MPIU___FP16));
918: #endif

920: #if defined(PETSC_USE_REAL___FLOAT128) || defined(PETSC_USE_REAL___FP16)
921:   PetscCallMPI(MPI_Op_create(PetscSum_Local, 1, &MPIU_SUM));
922:   PetscCallMPI(MPI_Op_create(PetscMax_Local, 1, &MPIU_MAX));
923:   PetscCallMPI(MPI_Op_create(PetscMin_Local, 1, &MPIU_MIN));
924: #elif defined(PETSC_HAVE_REAL___FLOAT128) || defined(PETSC_HAVE_REAL___FP16)
925:   PetscCallMPI(MPI_Op_create(PetscSum_Local, 1, &MPIU_SUM___FP16___FLOAT128));
926: #endif

928:   PetscCallMPI(MPI_Type_contiguous(2, MPIU_SCALAR, &MPIU_2SCALAR));
929:   PetscCallMPI(MPI_Op_create(PetscGarbageKeySortedIntersect, 1, &Petsc_Garbage_SetIntersectOp));
930:   PetscCallMPI(MPI_Type_commit(&MPIU_2SCALAR));

932:   /* create datatypes used by MPIU_MAXLOC, MPIU_MINLOC and PetscSplitReduction_Op */
933: #if !defined(PETSC_HAVE_MPIUNI)
934:   {
935:     PetscMPIInt  blockSizes[2]   = {1, 1};
936:     MPI_Aint     blockOffsets[2] = {offsetof(struct petsc_mpiu_real_int, v), offsetof(struct petsc_mpiu_real_int, i)};
937:     MPI_Datatype blockTypes[2]   = {MPIU_REAL, MPIU_INT}, tmpStruct;

939:     PetscCallMPI(MPI_Type_create_struct(2, blockSizes, blockOffsets, blockTypes, &tmpStruct));
940:     PetscCallMPI(MPI_Type_create_resized(tmpStruct, 0, sizeof(struct petsc_mpiu_real_int), &MPIU_REAL_INT));
941:     PetscCallMPI(MPI_Type_free(&tmpStruct));
942:     PetscCallMPI(MPI_Type_commit(&MPIU_REAL_INT));
943:   }
944:   {
945:     PetscMPIInt  blockSizes[2]   = {1, 1};
946:     MPI_Aint     blockOffsets[2] = {offsetof(struct petsc_mpiu_scalar_int, v), offsetof(struct petsc_mpiu_scalar_int, i)};
947:     MPI_Datatype blockTypes[2]   = {MPIU_SCALAR, MPIU_INT}, tmpStruct;

949:     PetscCallMPI(MPI_Type_create_struct(2, blockSizes, blockOffsets, blockTypes, &tmpStruct));
950:     PetscCallMPI(MPI_Type_create_resized(tmpStruct, 0, sizeof(struct petsc_mpiu_scalar_int), &MPIU_SCALAR_INT));
951:     PetscCallMPI(MPI_Type_free(&tmpStruct));
952:     PetscCallMPI(MPI_Type_commit(&MPIU_SCALAR_INT));
953:   }
954: #endif

956: #if defined(PETSC_USE_64BIT_INDICES)
957:   PetscCallMPI(MPI_Type_contiguous(2, MPIU_INT, &MPIU_2INT));
958:   PetscCallMPI(MPI_Type_commit(&MPIU_2INT));
959: #endif
960:   PetscCallMPI(MPI_Type_contiguous(4, MPI_INT, &MPI_4INT));
961:   PetscCallMPI(MPI_Type_commit(&MPI_4INT));
962:   PetscCallMPI(MPI_Type_contiguous(4, MPIU_INT, &MPIU_4INT));
963:   PetscCallMPI(MPI_Type_commit(&MPIU_4INT));

965:   /*
966:      Attributes to be set on PETSc communicators
967:   */
968:   PetscCallMPI(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, Petsc_Counter_Attr_Delete_Fn, &Petsc_Counter_keyval, (void *)0));
969:   PetscCallMPI(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, Petsc_InnerComm_Attr_Delete_Fn, &Petsc_InnerComm_keyval, (void *)0));
970:   PetscCallMPI(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, Petsc_OuterComm_Attr_Delete_Fn, &Petsc_OuterComm_keyval, (void *)0));
971:   PetscCallMPI(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, Petsc_ShmComm_Attr_Delete_Fn, &Petsc_ShmComm_keyval, (void *)0));
972:   PetscCallMPI(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, MPI_COMM_NULL_DELETE_FN, &Petsc_CreationIdx_keyval, (void *)0));
973:   PetscCallMPI(MPI_Comm_create_keyval(MPI_COMM_NULL_COPY_FN, MPI_COMM_NULL_DELETE_FN, &Petsc_Garbage_HMap_keyval, (void *)0));

975: #if defined(PETSC_USE_FORTRAN_BINDINGS)
976:   if (ftn) PetscCall(PetscInitFortran_Private(readarguments, file, len));
977:   else
978: #endif
979:     PetscCall(PetscOptionsInsert(NULL, &PetscGlobalArgc, &PetscGlobalArgs, file));

981:   /* call a second time so it can look in the options database */
982:   PetscCall(PetscErrorPrintfInitialize());

984:   /*
985:      Check system options and print help
986:   */
987:   PetscCall(PetscOptionsCheckInitial_Private(help));

989:   /*
990:     Creates the logging data structures; this is enabled even if logging is not turned on
991:     This is the last thing we do before returning to the user code to prevent having the
992:     logging numbers contaminated by any startup time associated with MPI
993:   */
994:   PetscCall(PetscLogInitialize());

996:   /*
997:    Initialize PetscDevice and PetscDeviceContext

999:    Note to any future devs thinking of moving this, proper initialization requires:
1000:    1. MPI initialized
1001:    2. Options DB initialized
1002:    3. Petsc error handling initialized, specifically signal handlers. This expects to set up
1003:       its own SIGSEV handler via the push/pop interface.
1004:    4. Logging initialized
1005:   */
1006:   PetscCall(PetscDeviceInitializeFromOptions_Internal(PETSC_COMM_WORLD));

1008: #if PetscDefined(HAVE_VIENNACL)
1009:   flg = PETSC_FALSE;
1010:   PetscCall(PetscOptionsHasName(NULL, NULL, "-log_view", &flg));
1011:   if (!flg) PetscCall(PetscOptionsGetBool(NULL, NULL, "-viennacl_synchronize", &flg, NULL));
1012:   PetscViennaCLSynchronize = flg;
1013:   PetscCall(PetscViennaCLInit());
1014: #endif

1016:   PetscCall(PetscCitationsInitialize());

1018: #if defined(PETSC_HAVE_SAWS)
1019:   PetscCall(PetscInitializeSAWs(ftn ? NULL : help));
1020:   flg = PETSC_FALSE;
1021:   PetscCall(PetscOptionsHasName(NULL, NULL, "-stack_view", &flg));
1022:   if (flg) PetscCall(PetscStackViewSAWs());
1023: #endif

1025:   /*
1026:      Load the dynamic libraries (on machines that support them), this registers all
1027:      the solvers etc. (On non-dynamic machines this initializes the PetscDraw and PetscViewer classes)
1028:   */
1029:   PetscCall(PetscInitialize_DynamicLibraries());

1031:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
1032:   PetscCall(PetscInfo(NULL, "PETSc successfully started: number of processors = %d\n", size));
1033:   PetscCall(PetscGetHostName(hostname, sizeof(hostname)));
1034:   PetscCall(PetscInfo(NULL, "Running on machine: %s\n", hostname));
1035: #if defined(PETSC_HAVE_OPENMP)
1036:   {
1037:     PetscBool omp_view_flag;
1038:     char     *threads = getenv("OMP_NUM_THREADS");

1040:     if (threads) {
1041:       PetscCall(PetscInfo(NULL, "Number of OpenMP threads %s (as given by OMP_NUM_THREADS)\n", threads));
1042:       (void)sscanf(threads, "%" PetscInt_FMT, &PetscNumOMPThreads);
1043:     } else {
1044:       PetscNumOMPThreads = (PetscInt)omp_get_max_threads();
1045:       PetscCall(PetscInfo(NULL, "Number of OpenMP threads %" PetscInt_FMT " (as given by omp_get_max_threads())\n", PetscNumOMPThreads));
1046:     }
1047:     PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "OpenMP options", "Sys");
1048:     PetscCall(PetscOptionsInt("-omp_num_threads", "Number of OpenMP threads to use (can also use environmental variable OMP_NUM_THREADS", "None", PetscNumOMPThreads, &PetscNumOMPThreads, &flg));
1049:     PetscCall(PetscOptionsName("-omp_view", "Display OpenMP number of threads", NULL, &omp_view_flag));
1050:     PetscOptionsEnd();
1051:     if (flg) {
1052:       PetscCall(PetscInfo(NULL, "Number of OpenMP threads %" PetscInt_FMT " (given by -omp_num_threads)\n", PetscNumOMPThreads));
1053:       omp_set_num_threads((int)PetscNumOMPThreads);
1054:     }
1055:     if (omp_view_flag) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "OpenMP: number of threads %" PetscInt_FMT "\n", PetscNumOMPThreads));
1056:   }
1057: #endif

1059: #if defined(PETSC_USE_PETSC_MPI_EXTERNAL32)
1060:   /*
1061:       Tell MPI about our own data representation converter, this would/should be used if extern32 is not supported by the MPI

1063:       Currently not used because it is not supported by MPICH.
1064:   */
1065:   if (!PetscBinaryBigEndian()) PetscCallMPI(MPI_Register_datarep((char *)"petsc", PetscDataRep_read_conv_fn, PetscDataRep_write_conv_fn, PetscDataRep_extent_fn, NULL));
1066: #endif

1068: #if defined(PETSC_SERIALIZE_FUNCTIONS)
1069:   PetscCall(PetscFPTCreate(10000));
1070: #endif

1072: #if defined(PETSC_HAVE_HWLOC)
1073:   {
1074:     PetscViewer viewer;
1075:     PetscCall(PetscOptionsGetViewer(PETSC_COMM_WORLD, NULL, NULL, "-process_view", &viewer, NULL, &flg));
1076:     if (flg) {
1077:       PetscCall(PetscProcessPlacementView(viewer));
1078:       PetscCall(PetscViewerDestroy(&viewer));
1079:     }
1080:   }
1081: #endif

1083:   flg = PETSC_TRUE;
1084:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-viewfromoptions", &flg, NULL));
1085:   if (!flg) PetscCall(PetscOptionsPushGetViewerOff(PETSC_TRUE));

1087: #if defined(PETSC_HAVE_ADIOS)
1088:   PetscCallExternal(adios_init_noxml, PETSC_COMM_WORLD);
1089:   PetscCallExternal(adios_declare_group, &Petsc_adios_group, "PETSc", "", adios_stat_default);
1090:   PetscCallExternal(adios_select_method, Petsc_adios_group, "MPI", "", "");
1091:   PetscCallExternal(adios_read_init_method, ADIOS_READ_METHOD_BP, PETSC_COMM_WORLD, "");
1092: #endif

1094: #if defined(__VALGRIND_H)
1095:   PETSC_RUNNING_ON_VALGRIND = RUNNING_ON_VALGRIND ? PETSC_TRUE : PETSC_FALSE;
1096:   #if defined(PETSC_USING_DARWIN) && defined(PETSC_BLASLAPACK_SDOT_RETURNS_DOUBLE)
1097:   if (PETSC_RUNNING_ON_VALGRIND) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "WARNING: Running valgrind with the macOS native BLAS and LAPACK can fail. If it fails, try configuring with --download-fblaslapack or --download-f2cblaslapack"));
1098:   #endif
1099: #endif
1100:   /*
1101:       Set flag that we are completely initialized
1102:   */
1103:   PetscInitializeCalled = PETSC_TRUE;

1105:   PetscCall(PetscOptionsHasName(NULL, NULL, "-python", &flg));
1106:   if (flg) PetscCall(PetscPythonInitialize(NULL, NULL));

1108:   PetscCall(PetscOptionsHasName(NULL, NULL, "-mpi_linear_solver_server", &flg));
1109:   if (PetscDefined(USE_SINGLE_LIBRARY) && flg) PetscCall(PCMPIServerBegin());
1110:   else PetscCheck(!flg, PETSC_COMM_WORLD, PETSC_ERR_SUP, "PETSc configured using -with-single-library=0; -mpi_linear_solver_server not supported in that case");
1111:   PetscFunctionReturn(PETSC_SUCCESS);
1112: }

1114: // "Unknown section 'Environmental Variables'"
1115: // PetscClangLinter pragma disable: -fdoc-section-header-unknown
1116: /*@C
1117:   PetscInitialize - Initializes the PETSc database and MPI.
1118:   `PetscInitialize()` calls MPI_Init() if that has yet to be called,
1119:   so this routine should always be called near the beginning of
1120:   your program -- usually the very first line!

1122:   Collective on `MPI_COMM_WORLD` or `PETSC_COMM_WORLD` if it has been set

1124:   Input Parameters:
1125: + argc - count of number of command line arguments
1126: . args - the command line arguments
1127: . file - [optional] PETSc database file, append ":yaml" to filename to specify YAML options format.
1128:           Use NULL or empty string to not check for code specific file.
1129:           Also checks ~/.petscrc, .petscrc and petscrc.
1130:           Use -skip_petscrc in the code specific file (or command line) to skip ~/.petscrc, .petscrc and petscrc files.
1131: - help - [optional] Help message to print, use NULL for no message

1133:    If you wish PETSc code to run ONLY on a subcommunicator of `MPI_COMM_WORLD`, create that
1134:    communicator first and assign it to `PETSC_COMM_WORLD` BEFORE calling `PetscInitialize()`. Thus if you are running a
1135:    four process job and two processes will run PETSc and have `PetscInitialize()` and PetscFinalize() and two process will not,
1136:    then do this. If ALL processes in the job are using `PetscInitialize()` and `PetscFinalize()` then you don't need to do this, even
1137:    if different subcommunicators of the job are doing different things with PETSc.

1139:   Options Database Keys:
1140: + -help [intro]                                       - prints help method for each option; if intro is given the program stops after printing the introductory help message
1141: . -start_in_debugger [noxterm,dbx,xdb,gdb,...]        - Starts program in debugger
1142: . -on_error_attach_debugger [noxterm,dbx,xdb,gdb,...] - Starts debugger when error detected
1143: . -on_error_emacs <machinename>                       - causes emacsclient to jump to error file
1144: . -on_error_abort                                     - calls `abort()` when error detected (no traceback)
1145: . -on_error_mpiabort                                  - calls `MPI_abort()` when error detected
1146: . -error_output_stdout                                - prints PETSc error messages to stdout instead of the default stderr
1147: . -error_output_none                                  - does not print the error messages (but handles errors in the same way as if this was not called)
1148: . -debugger_ranks [rank1,rank2,...]                   - Indicates ranks to start in debugger
1149: . -debugger_pause [sleeptime] (in seconds)            - Pauses debugger
1150: . -stop_for_debugger                                  - Print message on how to attach debugger manually to
1151:                         process and wait (-debugger_pause) seconds for attachment
1152: . -malloc_dump                                        - prints a list of all unfreed memory at the end of the run
1153: . -malloc_test                                        - like -malloc_dump -malloc_debug, but only active for debugging builds, ignored in optimized build. May want to set in PETSC_OPTIONS environmental variable
1154: . -malloc_view                                        - show a list of all allocated memory during `PetscFinalize()`
1155: . -malloc_view_threshold <t>                          - only list memory allocations of size greater than t with -malloc_view
1156: . -malloc_requested_size                              - malloc logging will record the requested size rather than size after alignment
1157: . -fp_trap                                            - Stops on floating point exceptions
1158: . -no_signal_handler                                  - Indicates not to trap error signals
1159: . -shared_tmp                                         - indicates /tmp directory is shared by all processors
1160: . -not_shared_tmp                                     - each processor has own /tmp
1161: . -tmp                                                - alternative name of /tmp directory
1162: . -get_total_flops                                    - returns total flops done by all processors
1163: - -memory_view                                        - Print memory usage at end of run

1165:   Options Database Keys for Option Database:
1166: + -skip_petscrc           - skip the default option files ~/.petscrc, .petscrc, petscrc
1167: . -options_monitor        - monitor all set options to standard output for the whole program run
1168: - -options_monitor_cancel - cancel options monitoring hard-wired using `PetscOptionsMonitorSet()`

1170:    Options -options_monitor_{all,cancel} are
1171:    position-independent and apply to all options set since the PETSc start.
1172:    They can be used also in option files.

1174:    See `PetscOptionsMonitorSet()` to do monitoring programmatically.

1176:   Options Database Keys for Profiling:
1177:    See Users-Manual: ch_profiling for details.
1178: + -info [filename][:[~]<list,of,classnames>[:[~]self]] - Prints verbose information. See `PetscInfo()`.
1179: . -log_sync                                            - Enable barrier synchronization for all events. This option is useful to debug imbalance within each event,
1180:         however it slows things down and gives a distorted view of the overall runtime.
1181: . -log_trace [filename]                                - Print traces of all PETSc calls to the screen (useful to determine where a program
1182:         hangs without running in the debugger).  See `PetscLogTraceBegin()`.
1183: . -log_view [:filename:format][,[:filename:format]...] - Prints summary of flop and timing information to screen or file, see `PetscLogView()` (up to 4 viewers)
1184: . -log_view_memory                                     - Includes in the summary from -log_view the memory used in each event, see `PetscLogView()`.
1185: . -log_view_gpu_time                                   - Includes in the summary from -log_view the time used in each GPU kernel, see `PetscLogView().
1186: . -log_exclude: <vec,mat,pc,ksp,snes>                  - excludes subset of object classes from logging
1187: . -log [filename]                                      - Logs profiling information in a dump file, see `PetscLogDump()`.
1188: . -log_all [filename]                                  - Same as `-log`.
1189: . -log_mpe [filename]                                  - Creates a logfile viewable by the utility Jumpshot (in MPICH distribution)
1190: . -log_perfstubs                                       - Starts a log handler with the perfstubs interface (which is used by TAU)
1191: . -log_nvtx                                            - Starts an nvtx log handler for use with Nsight
1192: . -viewfromoptions on,off                              - Enable or disable `XXXSetFromOptions()` calls, for applications with many small solves turn this off
1193: - -check_pointer_intensity 0,1,2                       - if pointers are checked for validity (debug version only), using 0 will result in faster code

1195:   Options Database Keys for SAWs:
1196: + -saws_port <portnumber>        - port number to publish SAWs data, default is 8080
1197: . -saws_port_auto_select         - have SAWs select a new unique port number where it publishes the data, the URL is printed to the screen
1198:                                    this is useful when you are running many jobs that utilize SAWs at the same time
1199: . -saws_log <filename>           - save a log of all SAWs communication
1200: . -saws_https <certificate file> - have SAWs use HTTPS instead of HTTP
1201: - -saws_root <directory>         - allow SAWs to have access to the given directory to search for requested resources and files

1203:   Environmental Variables:
1204: +   `PETSC_TMP` - alternative tmp directory
1205: .   `PETSC_SHARED_TMP` - tmp is shared by all processes
1206: .   `PETSC_NOT_SHARED_TMP` - each process has its own private tmp
1207: .   `PETSC_OPTIONS` - a string containing additional options for petsc in the form of command line "-key value" pairs
1208: .   `PETSC_OPTIONS_YAML` - (requires configuring PETSc to use libyaml) a string containing additional options for petsc in the form of a YAML document
1209: .   `PETSC_VIEWER_SOCKET_PORT` - socket number to use for socket viewer
1210: -   `PETSC_VIEWER_SOCKET_MACHINE` - machine to use for socket viewer to connect to

1212:   Level: beginner

1214:   Note:
1215:   If for some reason you must call `MPI_Init()` separately, call
1216:   it before `PetscInitialize()`.

1218:   Fortran Notes:
1219:   In Fortran this routine can be called with
1220: .vb
1221:        call PetscInitialize(ierr)
1222:        call PetscInitialize(file,ierr) or
1223:        call PetscInitialize(file,help,ierr)
1224: .ve

1226:   If your main program is C but you call Fortran code that also uses PETSc you need to call `PetscInitializeFortran()` soon after
1227:   calling `PetscInitialize()`.

1229:   Options Database Key for Developers:
1230: . -checkfunctionlist - automatically checks that function lists associated with objects are correctly cleaned up. Produces messages of the form:
1231:     "function name: MatInodeGetInodeSizes_C" if they are not cleaned up. This flag is always set for the test harness (in framework.py)

1233: .seealso: `PetscFinalize()`, `PetscInitializeFortran()`, `PetscGetArgs()`, `PetscInitializeNoArguments()`, `PetscLogGpuTime()`
1234: @*/
1235: PetscErrorCode PetscInitialize(int *argc, char ***args, const char file[], const char help[])
1236: {
1237:   PetscMPIInt flag;
1238:   const char *prog = "Unknown Name", *mpienv;

1240:   PetscFunctionBegin;
1241:   if (PetscInitializeCalled) PetscFunctionReturn(PETSC_SUCCESS);
1242:   PetscCallMPI(MPI_Initialized(&flag));
1243:   if (!flag) {
1244:     PetscCheck(PETSC_COMM_WORLD == MPI_COMM_NULL, PETSC_COMM_SELF, PETSC_ERR_SUP, "You cannot set PETSC_COMM_WORLD if you have not initialized MPI first");
1245:     PetscCall(PetscPreMPIInit_Private());
1246: #if defined(PETSC_HAVE_MPI_INIT_THREAD)
1247:     {
1248:       PetscMPIInt provided;
1249:       PetscCallMPI(MPI_Init_thread(argc, args, PETSC_MPI_THREAD_REQUIRED == PETSC_DECIDE ? MPI_THREAD_FUNNELED : PETSC_MPI_THREAD_REQUIRED, &provided));
1250:       PetscCheck(PETSC_MPI_THREAD_REQUIRED == PETSC_DECIDE || provided >= PETSC_MPI_THREAD_REQUIRED, PETSC_COMM_SELF, PETSC_ERR_MPI, "The MPI implementation's provided thread level is less than what you required");
1251:       if (PETSC_MPI_THREAD_REQUIRED == PETSC_DECIDE) PETSC_MPI_THREAD_REQUIRED = MPI_THREAD_FUNNELED; // assign it a valid value after check-up
1252:     }
1253: #else
1254:     PetscCallMPI(MPI_Init(argc, args));
1255: #endif
1256:     if (PetscDefined(HAVE_MPIUNI)) {
1257:       mpienv = getenv("PMI_SIZE");
1258:       if (!mpienv) mpienv = getenv("OMPI_COMM_WORLD_SIZE");
1259:       if (mpienv) {
1260:         PetscInt isize;
1261:         PetscCall(PetscOptionsStringToInt(mpienv, &isize));
1262:         if (isize != 1) printf("You are using an MPI-uni (sequential) install of PETSc but trying to launch parallel jobs; you need full MPI version of PETSc\n");
1263:         PetscCheck(isize == 1, MPI_COMM_SELF, PETSC_ERR_MPI, "You are using an MPI-uni (sequential) install of PETSc but trying to launch parallel jobs; you need full MPI version of PETSc");
1264:       }
1265:     }
1266:     PetscBeganMPI = PETSC_TRUE;
1267:   }

1269:   if (argc && *argc) prog = **args;
1270:   if (argc && args) {
1271:     PetscGlobalArgc = *argc;
1272:     PetscGlobalArgs = *args;
1273:   }
1274:   PetscCall(PetscInitialize_Common(prog, file, help, PETSC_FALSE, PETSC_FALSE, 0));
1275:   PetscFunctionReturn(PETSC_SUCCESS);
1276: }

1278: PETSC_INTERN PetscObject *PetscObjects;
1279: PETSC_INTERN PetscInt     PetscObjectsCounts;
1280: PETSC_INTERN PetscInt     PetscObjectsMaxCounts;
1281: PETSC_INTERN PetscBool    PetscObjectsLog;

1283: /*
1284:     Frees all the MPI types and operations that PETSc may have created
1285: */
1286: PetscErrorCode PetscFreeMPIResources(void)
1287: {
1288:   PetscFunctionBegin;
1289: #if defined(PETSC_HAVE_REAL___FLOAT128)
1290:   PetscCallMPI(MPI_Type_free(&MPIU___FLOAT128));
1291:   #if defined(PETSC_HAVE_COMPLEX)
1292:   PetscCallMPI(MPI_Type_free(&MPIU___COMPLEX128));
1293:   #endif
1294: #endif
1295: #if defined(PETSC_HAVE_REAL___FP16)
1296:   PetscCallMPI(MPI_Type_free(&MPIU___FP16));
1297: #endif

1299: #if defined(PETSC_USE_REAL___FLOAT128) || defined(PETSC_USE_REAL___FP16)
1300:   PetscCallMPI(MPI_Op_free(&MPIU_SUM));
1301:   PetscCallMPI(MPI_Op_free(&MPIU_MAX));
1302:   PetscCallMPI(MPI_Op_free(&MPIU_MIN));
1303: #elif defined(PETSC_HAVE_REAL___FLOAT128) || defined(PETSC_HAVE_REAL___FP16)
1304:   PetscCallMPI(MPI_Op_free(&MPIU_SUM___FP16___FLOAT128));
1305: #endif

1307:   PetscCallMPI(MPI_Type_free(&MPIU_2SCALAR));
1308:   PetscCallMPI(MPI_Type_free(&MPIU_REAL_INT));
1309:   PetscCallMPI(MPI_Type_free(&MPIU_SCALAR_INT));
1310: #if defined(PETSC_USE_64BIT_INDICES)
1311:   PetscCallMPI(MPI_Type_free(&MPIU_2INT));
1312: #endif
1313:   PetscCallMPI(MPI_Type_free(&MPI_4INT));
1314:   PetscCallMPI(MPI_Type_free(&MPIU_4INT));
1315:   PetscCallMPI(MPI_Op_free(&MPIU_MAXSUM_OP));
1316:   PetscCallMPI(MPI_Op_free(&Petsc_Garbage_SetIntersectOp));
1317:   PetscFunctionReturn(PETSC_SUCCESS);
1318: }

1320: PETSC_INTERN PetscErrorCode PetscLogFinalize(void);

1322: /*@C
1323:   PetscFinalize - Checks for options to be called at the conclusion
1324:   of the program. `MPI_Finalize()` is called only if the user had not
1325:   called `MPI_Init()` before calling `PetscInitialize()`.

1327:   Collective on `PETSC_COMM_WORLD`

1329:   Options Database Keys:
1330: + -options_view                    - Calls `PetscOptionsView()`
1331: . -options_left                    - Prints unused options that remain in the database
1332: . -objects_dump [all]              - Prints list of objects allocated by the user that have not been freed, the option all cause all outstanding objects to be listed
1333: . -mpidump                         - Calls PetscMPIDump()
1334: . -malloc_dump <optional filename> - Calls `PetscMallocDump()`, displays all memory allocated that has not been freed
1335: . -memory_view                     - Prints total memory usage
1336: - -malloc_view <optional filename> - Prints list of all memory allocated and in what functions

1338:   Level: beginner

1340:   Note:
1341:   See `PetscInitialize()` for other runtime options.

1343: .seealso: `PetscInitialize()`, `PetscOptionsView()`, `PetscMallocDump()`, `PetscMPIDump()`, `PetscEnd()`
1344: @*/
1345: PetscErrorCode PetscFinalize(void)
1346: {
1347:   PetscMPIInt rank;
1348:   PetscInt    nopt;
1349:   PetscBool   flg1 = PETSC_FALSE, flg2 = PETSC_FALSE, flg3 = PETSC_FALSE;
1350:   PetscBool   flg;
1351:   char        mname[PETSC_MAX_PATH_LEN];

1353:   PetscFunctionBegin;
1354:   PetscCheck(PetscInitializeCalled, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "PetscInitialize() must be called before PetscFinalize()");
1355:   PetscCall(PetscInfo(NULL, "PetscFinalize() called\n"));

1357:   PetscCall(PetscOptionsHasName(NULL, NULL, "-mpi_linear_solver_server", &flg));
1358:   if (PetscDefined(USE_SINGLE_LIBRARY) && flg) PetscCall(PCMPIServerEnd());

1360:   /* Clean up Garbage automatically on COMM_SELF and COMM_WORLD at finalize */
1361:   {
1362:     union
1363:     {
1364:       MPI_Comm comm;
1365:       void    *ptr;
1366:     } ucomm;
1367:     PetscMPIInt flg;
1368:     void       *tmp;

1370:     PetscCallMPI(MPI_Comm_get_attr(PETSC_COMM_SELF, Petsc_InnerComm_keyval, &ucomm, &flg));
1371:     if (flg) PetscCallMPI(MPI_Comm_get_attr(ucomm.comm, Petsc_Garbage_HMap_keyval, &tmp, &flg));
1372:     if (flg) PetscCall(PetscGarbageCleanup(PETSC_COMM_SELF));
1373:     PetscCallMPI(MPI_Comm_get_attr(PETSC_COMM_WORLD, Petsc_InnerComm_keyval, &ucomm, &flg));
1374:     if (flg) PetscCallMPI(MPI_Comm_get_attr(ucomm.comm, Petsc_Garbage_HMap_keyval, &tmp, &flg));
1375:     if (flg) PetscCall(PetscGarbageCleanup(PETSC_COMM_WORLD));
1376:   }

1378:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
1379: #if defined(PETSC_HAVE_ADIOS)
1380:   PetscCallExternal(adios_read_finalize_method, ADIOS_READ_METHOD_BP_AGGREGATE);
1381:   PetscCallExternal(adios_finalize, rank);
1382: #endif
1383:   PetscCall(PetscOptionsHasName(NULL, NULL, "-citations", &flg));
1384:   if (flg) {
1385:     char *cits, filename[PETSC_MAX_PATH_LEN];
1386:     FILE *fd = PETSC_STDOUT;

1388:     PetscCall(PetscOptionsGetString(NULL, NULL, "-citations", filename, sizeof(filename), NULL));
1389:     if (filename[0]) PetscCall(PetscFOpen(PETSC_COMM_WORLD, filename, "w", &fd));
1390:     PetscCall(PetscSegBufferGet(PetscCitationsList, 1, &cits));
1391:     cits[0] = 0;
1392:     PetscCall(PetscSegBufferExtractAlloc(PetscCitationsList, &cits));
1393:     PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "If you publish results based on this computation please cite the following:\n"));
1394:     PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "===========================================================================\n"));
1395:     PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "%s", cits));
1396:     PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "===========================================================================\n"));
1397:     PetscCall(PetscFClose(PETSC_COMM_WORLD, fd));
1398:     PetscCall(PetscFree(cits));
1399:   }
1400:   PetscCall(PetscSegBufferDestroy(&PetscCitationsList));

1402: #if defined(PETSC_SERIALIZE_FUNCTIONS)
1403:   PetscCall(PetscFPTDestroy());
1404: #endif

1406: #if defined(PETSC_HAVE_SAWS)
1407:   flg = PETSC_FALSE;
1408:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-saw_options", &flg, NULL));
1409:   if (flg) PetscCall(PetscOptionsSAWsDestroy());
1410: #endif

1412: #if defined(PETSC_HAVE_X)
1413:   flg1 = PETSC_FALSE;
1414:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-x_virtual", &flg1, NULL));
1415:   if (flg1) {
1416:     /*  this is a crude hack, but better than nothing */
1417:     PetscCall(PetscPOpen(PETSC_COMM_WORLD, NULL, "pkill -9 Xvfb", "r", NULL));
1418:   }
1419: #endif

1421: #if !defined(PETSC_HAVE_THREADSAFETY)
1422:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-memory_view", &flg2, NULL));
1423:   if (flg2) PetscCall(PetscMemoryView(PETSC_VIEWER_STDOUT_WORLD, "Summary of Memory Usage in PETSc\n"));
1424: #endif

1426:   if (PetscDefined(USE_LOG)) {
1427:     flg1 = PETSC_FALSE;
1428:     PetscCall(PetscOptionsGetBool(NULL, NULL, "-get_total_flops", &flg1, NULL));
1429:     if (flg1) {
1430:       PetscLogDouble flops = 0;
1431:       PetscCallMPI(MPI_Reduce(&petsc_TotalFlops, &flops, 1, MPI_DOUBLE, MPI_SUM, 0, PETSC_COMM_WORLD));
1432:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Total flops over all processors %g\n", flops));
1433:     }
1434:   }

1436:   if (PetscDefined(USE_LOG) && PetscDefined(HAVE_MPE)) {
1437:     mname[0] = 0;
1438:     PetscCall(PetscOptionsGetString(NULL, NULL, "-log_mpe", mname, sizeof(mname), &flg1));
1439:     if (flg1) PetscCall(PetscLogMPEDump(mname[0] ? mname : NULL));
1440:   }

1442:   // Free all objects registered with PetscObjectRegisterDestroy() such as PETSC_VIEWER_XXX_().
1443:   PetscCall(PetscObjectRegisterDestroyAll());

1445:   if (PetscDefined(USE_LOG)) {
1446:     PetscCall(PetscOptionsPushGetViewerOff(PETSC_FALSE));
1447:     PetscCall(PetscLogViewFromOptions());
1448:     PetscCall(PetscOptionsPopGetViewerOff());
1449:     //  It should be turned on with PetscLogGpuTime() and never turned off except in this place
1450:     PetscLogGpuTimeFlag = PETSC_FALSE;

1452:     // Free any objects created by the last block of code.
1453:     PetscCall(PetscObjectRegisterDestroyAll());

1455:     mname[0] = 0;
1456:     PetscCall(PetscOptionsGetString(NULL, NULL, "-log_all", mname, sizeof(mname), &flg1));
1457:     PetscCall(PetscOptionsGetString(NULL, NULL, "-log", mname, sizeof(mname), &flg2));
1458:     if (flg1 || flg2) PetscCall(PetscLogDump(mname));
1459:   }

1461:   flg1 = PETSC_FALSE;
1462:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-no_signal_handler", &flg1, NULL));
1463:   if (!flg1) PetscCall(PetscPopSignalHandler());
1464:   flg1 = PETSC_FALSE;
1465:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-mpidump", &flg1, NULL));
1466:   if (flg1) PetscCall(PetscMPIDump(stdout));
1467:   flg1 = PETSC_FALSE;
1468:   flg2 = PETSC_FALSE;
1469:   /* preemptive call to avoid listing this option in options table as unused */
1470:   PetscCall(PetscOptionsHasName(NULL, NULL, "-malloc_dump", &flg1));
1471:   PetscCall(PetscOptionsHasName(NULL, NULL, "-objects_dump", &flg1));
1472:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-options_view", &flg2, NULL));

1474:   if (flg2) {
1475:     PetscViewer viewer;
1476:     PetscCall(PetscViewerCreate(PETSC_COMM_WORLD, &viewer));
1477:     PetscCall(PetscViewerSetType(viewer, PETSCVIEWERASCII));
1478:     PetscCall(PetscOptionsView(NULL, viewer));
1479:     PetscCall(PetscViewerDestroy(&viewer));
1480:   }

1482:   /* to prevent PETSc -options_left from warning */
1483:   PetscCall(PetscOptionsHasName(NULL, NULL, "-nox", &flg1));
1484:   PetscCall(PetscOptionsHasName(NULL, NULL, "-nox_warning", &flg1));

1486:   flg3 = PETSC_FALSE; /* default value is required */
1487:   PetscCall(PetscOptionsGetBool(NULL, NULL, "-options_left", &flg3, &flg1));
1488:   if (!flg1) flg3 = PETSC_TRUE;
1489:   if (flg3) {
1490:     if (!flg2 && flg1) { /* have not yet printed the options */
1491:       PetscViewer viewer;
1492:       PetscCall(PetscViewerCreate(PETSC_COMM_WORLD, &viewer));
1493:       PetscCall(PetscViewerSetType(viewer, PETSCVIEWERASCII));
1494:       PetscCall(PetscOptionsView(NULL, viewer));
1495:       PetscCall(PetscViewerDestroy(&viewer));
1496:     }
1497:     PetscCall(PetscOptionsAllUsed(NULL, &nopt));
1498:     if (nopt) {
1499:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "WARNING! There are options you set that were not used!\n"));
1500:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "WARNING! could be spelling mistake, etc!\n"));
1501:       if (nopt == 1) {
1502:         PetscCall(PetscPrintf(PETSC_COMM_WORLD, "There is one unused database option. It is:\n"));
1503:       } else {
1504:         PetscCall(PetscPrintf(PETSC_COMM_WORLD, "There are %" PetscInt_FMT " unused database options. They are:\n", nopt));
1505:       }
1506:     } else if (flg3 && flg1) {
1507:       PetscCall(PetscPrintf(PETSC_COMM_WORLD, "There are no unused options.\n"));
1508:     }
1509:     PetscCall(PetscOptionsLeft(NULL));
1510:   }

1512: #if defined(PETSC_HAVE_SAWS)
1513:   if (!PetscGlobalRank) {
1514:     PetscCall(PetscStackSAWsViewOff());
1515:     PetscCallSAWs(SAWs_Finalize, ());
1516:   }
1517: #endif

1519:   /*
1520:        List all objects the user may have forgot to free
1521:   */
1522:   if (PetscDefined(USE_LOG) && PetscObjectsLog) {
1523:     PetscCall(PetscOptionsHasName(NULL, NULL, "-objects_dump", &flg1));
1524:     if (flg1) {
1525:       MPI_Comm local_comm;
1526:       char     string[64];

1528:       PetscCall(PetscOptionsGetString(NULL, NULL, "-objects_dump", string, sizeof(string), NULL));
1529:       PetscCallMPI(MPI_Comm_dup(PETSC_COMM_WORLD, &local_comm));
1530:       PetscCall(PetscSequentialPhaseBegin_Private(local_comm, 1));
1531:       PetscCall(PetscObjectsDump(stdout, (string[0] == 'a') ? PETSC_TRUE : PETSC_FALSE));
1532:       PetscCall(PetscSequentialPhaseEnd_Private(local_comm, 1));
1533:       PetscCallMPI(MPI_Comm_free(&local_comm));
1534:     }
1535:   }

1537:   PetscObjectsCounts    = 0;
1538:   PetscObjectsMaxCounts = 0;
1539:   PetscCall(PetscFree(PetscObjects));

1541:   /*
1542:      Destroy any packages that registered a finalize
1543:   */
1544:   PetscCall(PetscRegisterFinalizeAll());

1546:   PetscCall(PetscLogFinalize());

1548:   /*
1549:      Print PetscFunctionLists that have not been properly freed
1550:   */
1551:   if (PetscPrintFunctionList) PetscCall(PetscFunctionListPrintAll());

1553:   if (petsc_history) {
1554:     PetscCall(PetscCloseHistoryFile(&petsc_history));
1555:     petsc_history = NULL;
1556:   }
1557:   PetscCall(PetscOptionsHelpPrintedDestroy(&PetscOptionsHelpPrintedSingleton));
1558:   PetscCall(PetscInfoDestroy());

1560: #if !defined(PETSC_HAVE_THREADSAFETY)
1561:   if (!(PETSC_RUNNING_ON_VALGRIND)) {
1562:     char  fname[PETSC_MAX_PATH_LEN];
1563:     char  sname[PETSC_MAX_PATH_LEN];
1564:     FILE *fd;
1565:     int   err;

1567:     flg2 = PETSC_FALSE;
1568:     flg3 = PETSC_FALSE;
1569:     if (PetscDefined(USE_DEBUG)) PetscCall(PetscOptionsGetBool(NULL, NULL, "-malloc_test", &flg2, NULL));
1570:     PetscCall(PetscOptionsGetBool(NULL, NULL, "-malloc_debug", &flg3, NULL));
1571:     fname[0] = 0;
1572:     PetscCall(PetscOptionsGetString(NULL, NULL, "-malloc_dump", fname, sizeof(fname), &flg1));
1573:     if (flg1 && fname[0]) {
1574:       PetscCall(PetscSNPrintf(sname, sizeof(sname), "%s_%d", fname, rank));
1575:       fd = fopen(sname, "w");
1576:       PetscCheck(fd, PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open log file: %s", sname);
1577:       PetscCall(PetscMallocDump(fd));
1578:       err = fclose(fd);
1579:       PetscCheck(!err, PETSC_COMM_SELF, PETSC_ERR_SYS, "fclose() failed on file");
1580:     } else if (flg1 || flg2 || flg3) {
1581:       MPI_Comm local_comm;

1583:       PetscCallMPI(MPI_Comm_dup(PETSC_COMM_WORLD, &local_comm));
1584:       PetscCall(PetscSequentialPhaseBegin_Private(local_comm, 1));
1585:       PetscCall(PetscMallocDump(stdout));
1586:       PetscCall(PetscSequentialPhaseEnd_Private(local_comm, 1));
1587:       PetscCallMPI(MPI_Comm_free(&local_comm));
1588:     }
1589:     fname[0] = 0;
1590:     PetscCall(PetscOptionsGetString(NULL, NULL, "-malloc_view", fname, sizeof(fname), &flg1));
1591:     if (flg1 && fname[0]) {
1592:       PetscCall(PetscSNPrintf(sname, sizeof(sname), "%s_%d", fname, rank));
1593:       fd = fopen(sname, "w");
1594:       PetscCheck(fd, PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open log file: %s", sname);
1595:       PetscCall(PetscMallocView(fd));
1596:       err = fclose(fd);
1597:       PetscCheck(!err, PETSC_COMM_SELF, PETSC_ERR_SYS, "fclose() failed on file");
1598:     } else if (flg1) {
1599:       MPI_Comm local_comm;

1601:       PetscCallMPI(MPI_Comm_dup(PETSC_COMM_WORLD, &local_comm));
1602:       PetscCall(PetscSequentialPhaseBegin_Private(local_comm, 1));
1603:       PetscCall(PetscMallocView(stdout));
1604:       PetscCall(PetscSequentialPhaseEnd_Private(local_comm, 1));
1605:       PetscCallMPI(MPI_Comm_free(&local_comm));
1606:     }
1607:   }
1608: #endif

1610:   /*
1611:      Close any open dynamic libraries
1612:   */
1613:   PetscCall(PetscFinalize_DynamicLibraries());

1615:   /* Can be destroyed only after all the options are used */
1616:   PetscCall(PetscOptionsDestroyDefault());

1618:   PetscGlobalArgc = 0;
1619:   PetscGlobalArgs = NULL;

1621: #if defined(PETSC_HAVE_KOKKOS)
1622:   if (PetscBeganKokkos) {
1623:     PetscCall(PetscKokkosFinalize_Private());
1624:     PetscBeganKokkos       = PETSC_FALSE;
1625:     PetscKokkosInitialized = PETSC_FALSE;
1626:   }
1627: #endif

1629: #if defined(PETSC_HAVE_NVSHMEM)
1630:   if (PetscBeganNvshmem) {
1631:     PetscCall(PetscNvshmemFinalize());
1632:     PetscBeganNvshmem = PETSC_FALSE;
1633:   }
1634: #endif

1636:   PetscCall(PetscFreeMPIResources());

1638:   /*
1639:      Destroy any known inner MPI_Comm's and attributes pointing to them
1640:      Note this will not destroy any new communicators the user has created.

1642:      If all PETSc objects were not destroyed those left over objects will have hanging references to
1643:      the MPI_Comms that were freed; but that is ok because those PETSc objects will never be used again
1644:  */
1645:   {
1646:     PetscCommCounter *counter;
1647:     PetscMPIInt       flg;
1648:     MPI_Comm          icomm;
1649:     union
1650:     {
1651:       MPI_Comm comm;
1652:       void    *ptr;
1653:     } ucomm;
1654:     PetscCallMPI(MPI_Comm_get_attr(PETSC_COMM_SELF, Petsc_InnerComm_keyval, &ucomm, &flg));
1655:     if (flg) {
1656:       icomm = ucomm.comm;
1657:       PetscCallMPI(MPI_Comm_get_attr(icomm, Petsc_Counter_keyval, &counter, &flg));
1658:       PetscCheck(flg, PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "Inner MPI_Comm does not have expected tag/name counter, problem with corrupted memory");

1660:       PetscCallMPI(MPI_Comm_delete_attr(PETSC_COMM_SELF, Petsc_InnerComm_keyval));
1661:       PetscCallMPI(MPI_Comm_delete_attr(icomm, Petsc_Counter_keyval));
1662:       PetscCallMPI(MPI_Comm_free(&icomm));
1663:     }
1664:     PetscCallMPI(MPI_Comm_get_attr(PETSC_COMM_WORLD, Petsc_InnerComm_keyval, &ucomm, &flg));
1665:     if (flg) {
1666:       icomm = ucomm.comm;
1667:       PetscCallMPI(MPI_Comm_get_attr(icomm, Petsc_Counter_keyval, &counter, &flg));
1668:       PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_ARG_CORRUPT, "Inner MPI_Comm does not have expected tag/name counter, problem with corrupted memory");

1670:       PetscCallMPI(MPI_Comm_delete_attr(PETSC_COMM_WORLD, Petsc_InnerComm_keyval));
1671:       PetscCallMPI(MPI_Comm_delete_attr(icomm, Petsc_Counter_keyval));
1672:       PetscCallMPI(MPI_Comm_free(&icomm));
1673:     }
1674:   }

1676:   PetscCallMPI(MPI_Comm_free_keyval(&Petsc_Counter_keyval));
1677:   PetscCallMPI(MPI_Comm_free_keyval(&Petsc_InnerComm_keyval));
1678:   PetscCallMPI(MPI_Comm_free_keyval(&Petsc_OuterComm_keyval));
1679:   PetscCallMPI(MPI_Comm_free_keyval(&Petsc_ShmComm_keyval));
1680:   PetscCallMPI(MPI_Comm_free_keyval(&Petsc_CreationIdx_keyval));
1681:   PetscCallMPI(MPI_Comm_free_keyval(&Petsc_Garbage_HMap_keyval));

1683:   // Free keyvals which may be silently created by some routines
1684:   if (Petsc_SharedWD_keyval != MPI_KEYVAL_INVALID) PetscCallMPI(MPI_Comm_free_keyval(&Petsc_SharedWD_keyval));
1685:   if (Petsc_SharedTmp_keyval != MPI_KEYVAL_INVALID) PetscCallMPI(MPI_Comm_free_keyval(&Petsc_SharedTmp_keyval));

1687:   PetscCall(PetscSpinlockDestroy(&PetscViewerASCIISpinLockOpen));
1688:   PetscCall(PetscSpinlockDestroy(&PetscViewerASCIISpinLockStdout));
1689:   PetscCall(PetscSpinlockDestroy(&PetscViewerASCIISpinLockStderr));
1690:   PetscCall(PetscSpinlockDestroy(&PetscCommSpinLock));

1692:   if (PetscBeganMPI) {
1693:     PetscMPIInt flag;
1694:     PetscCallMPI(MPI_Finalized(&flag));
1695:     PetscCheck(!flag, PETSC_COMM_SELF, PETSC_ERR_LIB, "MPI_Finalize() has already been called, even though MPI_Init() was called by PetscInitialize()");
1696:     /* wait until the very last moment to disable error handling */
1697:     PetscErrorHandlingInitialized = PETSC_FALSE;
1698:     PetscCallMPI(MPI_Finalize());
1699:   } else PetscErrorHandlingInitialized = PETSC_FALSE;

1701:   /*

1703:      Note: In certain cases PETSC_COMM_WORLD is never MPI_Comm_free()ed because
1704:    the communicator has some outstanding requests on it. Specifically if the
1705:    flag PETSC_HAVE_BROKEN_REQUEST_FREE is set (for IBM MPI implementation). See
1706:    src/vec/utils/vpscat.c. Due to this the memory allocated in PetscCommDuplicate()
1707:    is never freed as it should be. Thus one may obtain messages of the form
1708:    [ 1] 8 bytes PetscCommDuplicate() line 645 in src/sys/mpiu.c indicating the
1709:    memory was not freed.

1711: */
1712:   PetscCall(PetscMallocClear());
1713:   PetscCall(PetscStackReset());

1715:   PetscInitializeCalled = PETSC_FALSE;
1716:   PetscFinalizeCalled   = PETSC_TRUE;
1717: #if defined(PETSC_USE_COVERAGE)
1718:   /*
1719:      flush gcov, otherwise during CI the flushing continues into the next pipeline resulting in git not being able to delete directories since the
1720:      gcov files are still being added to the directories as git tries to remove the directories.
1721:    */
1722:   __gcov_flush();
1723: #endif
1724:   /* To match PetscFunctionBegin() at the beginning of this function */
1725:   PetscStackClearTop;
1726:   return PETSC_SUCCESS;
1727: }

1729: #if defined(PETSC_MISSING_LAPACK_lsame_)
1730: PETSC_EXTERN int lsame_(char *a, char *b)
1731: {
1732:   if (*a == *b) return 1;
1733:   if (*a + 32 == *b) return 1;
1734:   if (*a - 32 == *b) return 1;
1735:   return 0;
1736: }
1737: #endif

1739: #if defined(PETSC_MISSING_LAPACK_lsame)
1740: PETSC_EXTERN int lsame(char *a, char *b)
1741: {
1742:   if (*a == *b) return 1;
1743:   if (*a + 32 == *b) return 1;
1744:   if (*a - 32 == *b) return 1;
1745:   return 0;
1746: }
1747: #endif