Actual source code: pinit.c
petsc-3.4.0 2013-05-13
2: /*
3: This file defines the initialization of PETSc, including PetscInitialize()
4: */
5: #define PETSC_DESIRE_COMPLEX
6: #include <petsc-private/petscimpl.h> /*I "petscsys.h" I*/
7: #include <petscviewer.h>
9: #if defined(PETSC_HAVE_CUDA)
10: #include <cublas.h>
11: #endif
13: #include <petscthreadcomm.h>
15: #if defined(PETSC_USE_LOG)
16: extern PetscErrorCode PetscLogBegin_Private(void);
17: #endif
18: extern PetscBool PetscHMPIWorker;
21: #if defined(PETSC_SERIALIZE_FUNCTIONS)
22: PetscFPT PetscFPTData = 0;
23: #endif
25: /* -----------------------------------------------------------------------------------------*/
27: extern FILE *petsc_history;
29: extern PetscErrorCode PetscInitialize_DynamicLibraries(void);
30: extern PetscErrorCode PetscFinalize_DynamicLibraries(void);
31: extern PetscErrorCode PetscFunctionListPrintAll(void);
32: extern PetscErrorCode PetscSequentialPhaseBegin_Private(MPI_Comm,int);
33: extern PetscErrorCode PetscSequentialPhaseEnd_Private(MPI_Comm,int);
34: extern PetscErrorCode PetscCloseHistoryFile(FILE**);
36: /* user may set this BEFORE calling PetscInitialize() */
37: MPI_Comm PETSC_COMM_WORLD = MPI_COMM_NULL;
39: PetscMPIInt Petsc_Counter_keyval = MPI_KEYVAL_INVALID;
40: PetscMPIInt Petsc_InnerComm_keyval = MPI_KEYVAL_INVALID;
41: PetscMPIInt Petsc_OuterComm_keyval = MPI_KEYVAL_INVALID;
43: /*
44: Declare and set all the string names of the PETSc enums
45: */
46: const char *const PetscBools[] = {"FALSE","TRUE","PetscBool","PETSC_",0};
47: const char *const PetscCopyModes[] = {"COPY_VALUES","OWN_POINTER","USE_POINTER","PetscCopyMode","PETSC_",0};
48: const char *const PetscDataTypes[] = {"INT","DOUBLE","COMPLEX","LONG","SHORT","FLOAT",
49: "CHAR","LOGICAL","ENUM","BOOL","LONGDOUBLE","OBJECT","FUNCTION","PetscDataType","PETSC_",0};
51: PetscBool PetscPreLoadingUsed = PETSC_FALSE;
52: PetscBool PetscPreLoadingOn = PETSC_FALSE;
54: /* pthread_key for PetscStack */
55: #if defined(PETSC_HAVE_PTHREADCLASSES) && !defined(PETSC_PTHREAD_LOCAL)
56: pthread_key_t petscstack;
57: #endif
59: /*
60: Checks the options database for initializations related to the
61: PETSc components
62: */
65: PetscErrorCode PetscOptionsCheckInitial_Components(void)
66: {
67: PetscBool flg1;
71: PetscOptionsHasName(NULL,"-help",&flg1);
72: if (flg1) {
73: #if defined(PETSC_USE_LOG)
74: MPI_Comm comm = PETSC_COMM_WORLD;
75: (*PetscHelpPrintf)(comm,"------Additional PETSc component options--------\n");
76: (*PetscHelpPrintf)(comm," -log_summary_exclude: <vec,mat,pc.ksp,snes>\n");
77: (*PetscHelpPrintf)(comm," -info_exclude: <null,vec,mat,pc,ksp,snes,ts>\n");
78: (*PetscHelpPrintf)(comm,"-----------------------------------------------\n");
79: #endif
80: }
81: return(0);
82: }
86: /*
87: PetscInitializeNoPointers - Calls PetscInitialize() from C/C++ without the pointers to argc and args
89: Collective
91: Level: advanced
93: Notes: this is called only by the PETSc MATLAB and Julia interface. Even though it might start MPI it sets the flag to
94: indicate that it did NOT start MPI so that the PetscFinalize() does not end MPI, thus allowing PetscInitialize() to
95: be called multiple times from MATLAB and Julia without the problem of trying to initialize MPI more than once.
97: Turns off PETSc signal handling because that can interact with MATLAB's signal handling causing random crashes.
99: .seealso: PetscInitialize(), PetscInitializeFortran(), PetscInitializeNoArguments()
100: */
101: PetscErrorCode PetscInitializeNoPointers(int argc,char **args,const char *filename,const char *help)
102: {
104: int myargc = argc;
105: char **myargs = args;
108: PetscInitialize(&myargc,&myargs,filename,help);
109: PetscPopSignalHandler();
110: PetscBeganMPI = PETSC_FALSE;
111: PetscFunctionReturn(ierr);
112: }
116: /*
117: Used by MATLAB and Julia interface to get communicator
118: */
119: PetscErrorCode PetscGetPETSC_COMM_SELF(MPI_Comm *comm)
120: {
122: *comm = PETSC_COMM_SELF;
123: return(0);
124: }
128: /*@C
129: PetscInitializeNoArguments - Calls PetscInitialize() from C/C++ without
130: the command line arguments.
132: Collective
134: Level: advanced
136: .seealso: PetscInitialize(), PetscInitializeFortran()
137: @*/
138: PetscErrorCode PetscInitializeNoArguments(void)
139: {
141: int argc = 0;
142: char **args = 0;
145: PetscInitialize(&argc,&args,NULL,NULL);
146: PetscFunctionReturn(ierr);
147: }
151: /*@
152: PetscInitialized - Determine whether PETSc is initialized.
154: 7 Level: beginner
156: .seealso: PetscInitialize(), PetscInitializeNoArguments(), PetscInitializeFortran()
157: @*/
158: PetscErrorCode PetscInitialized(PetscBool *isInitialized)
159: {
162: *isInitialized = PetscInitializeCalled;
163: return(0);
164: }
168: /*@
169: PetscFinalized - Determine whether PetscFinalize() has been called yet
171: Level: developer
173: .seealso: PetscInitialize(), PetscInitializeNoArguments(), PetscInitializeFortran()
174: @*/
175: PetscErrorCode PetscFinalized(PetscBool *isFinalized)
176: {
179: *isFinalized = PetscFinalizeCalled;
180: return(0);
181: }
183: extern PetscErrorCode PetscOptionsCheckInitial_Private(void);
185: /*
186: This function is the MPI reduction operation used to compute the sum of the
187: first half of the datatype and the max of the second half.
188: */
189: MPI_Op PetscMaxSum_Op = 0;
193: PETSC_EXTERN void MPIAPI PetscMaxSum_Local(void *in,void *out,int *cnt,MPI_Datatype *datatype)
194: {
195: PetscInt *xin = (PetscInt*)in,*xout = (PetscInt*)out,i,count = *cnt;
198: if (*datatype != MPIU_2INT) {
199: (*PetscErrorPrintf)("Can only handle MPIU_2INT data types");
200: MPI_Abort(MPI_COMM_WORLD,1);
201: }
203: for (i=0; i<count; i++) {
204: xout[2*i] = PetscMax(xout[2*i],xin[2*i]);
205: xout[2*i+1] += xin[2*i+1];
206: }
207: PetscFunctionReturnVoid();
208: }
210: /*
211: Returns the max of the first entry owned by this processor and the
212: sum of the second entry.
214: The reason nprocs[2*i] contains lengths nprocs[2*i+1] contains flag of 1 if length is nonzero
215: is so that the PetscMaxSum_Op() can set TWO values, if we passed in only nprocs[i] with lengths
216: there would be no place to store the both needed results.
217: */
220: PetscErrorCode PetscMaxSum(MPI_Comm comm,const PetscInt nprocs[],PetscInt *max,PetscInt *sum)
221: {
222: PetscMPIInt size,rank;
223: struct {PetscInt max,sum;} *work;
227: MPI_Comm_size(comm,&size);
228: MPI_Comm_rank(comm,&rank);
229: PetscMalloc(size*sizeof(*work),&work);
230: MPI_Allreduce((void*)nprocs,work,size,MPIU_2INT,PetscMaxSum_Op,comm);
231: *max = work[rank].max;
232: *sum = work[rank].sum;
233: PetscFree(work);
234: return(0);
235: }
237: /* ----------------------------------------------------------------------------*/
238: MPI_Op PetscADMax_Op = 0;
242: PETSC_EXTERN void MPIAPI PetscADMax_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
243: {
244: PetscScalar *xin = (PetscScalar*)in,*xout = (PetscScalar*)out;
245: PetscInt i,count = *cnt;
248: if (*datatype != MPIU_2SCALAR) {
249: (*PetscErrorPrintf)("Can only handle MPIU_2SCALAR data (i.e. double or complex) types");
250: MPI_Abort(MPI_COMM_WORLD,1);
251: }
253: for (i=0; i<count; i++) {
254: if (PetscRealPart(xout[2*i]) < PetscRealPart(xin[2*i])) {
255: xout[2*i] = xin[2*i];
256: xout[2*i+1] = xin[2*i+1];
257: }
258: }
259: PetscFunctionReturnVoid();
260: }
262: MPI_Op PetscADMin_Op = 0;
266: PETSC_EXTERN void MPIAPI PetscADMin_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
267: {
268: PetscScalar *xin = (PetscScalar*)in,*xout = (PetscScalar*)out;
269: PetscInt i,count = *cnt;
272: if (*datatype != MPIU_2SCALAR) {
273: (*PetscErrorPrintf)("Can only handle MPIU_2SCALAR data (i.e. double or complex) types");
274: MPI_Abort(MPI_COMM_WORLD,1);
275: }
277: for (i=0; i<count; i++) {
278: if (PetscRealPart(xout[2*i]) > PetscRealPart(xin[2*i])) {
279: xout[2*i] = xin[2*i];
280: xout[2*i+1] = xin[2*i+1];
281: }
282: }
283: PetscFunctionReturnVoid();
284: }
285: /* ---------------------------------------------------------------------------------------*/
287: #if (defined(PETSC_HAVE_COMPLEX) && !defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)) || defined(PETSC_USE_REAL___FLOAT128)
288: MPI_Op MPIU_SUM = 0;
292: PETSC_EXTERN void PetscSum_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
293: {
294: PetscInt i,count = *cnt;
297: if (*datatype == MPIU_REAL) {
298: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
299: for (i=0; i<count; i++) xout[i] += xin[i];
300: }
301: #if defined(PETSC_HAVE_COMPLEX)
302: else if (*datatype == MPIU_COMPLEX) {
303: PetscComplex *xin = (PetscComplex*)in,*xout = (PetscComplex*)out;
304: for (i=0; i<count; i++) xout[i] += xin[i];
305: }
306: #endif
307: else {
308: (*PetscErrorPrintf)("Can only handle MPIU_REAL or MPIU_COMPLEX data types");
309: MPI_Abort(MPI_COMM_WORLD,1);
310: }
311: PetscFunctionReturnVoid();
312: }
313: #endif
315: #if defined(PETSC_USE_REAL___FLOAT128)
316: MPI_Op MPIU_MAX = 0;
317: MPI_Op MPIU_MIN = 0;
321: PETSC_EXTERN void PetscMax_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
322: {
323: PetscInt i,count = *cnt;
326: if (*datatype == MPIU_REAL) {
327: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
328: for (i=0; i<count; i++) xout[i] = PetscMax(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++) {
334: xout[i] = PetscRealPartComplex(xout[i])<PetscRealPartComplex(xin[i]) ? xin[i] : xout[i];
335: }
336: }
337: #endif
338: else {
339: (*PetscErrorPrintf)("Can only handle MPIU_REAL or MPIU_COMPLEX data types");
340: MPI_Abort(MPI_COMM_WORLD,1);
341: }
342: PetscFunctionReturnVoid();
343: }
347: PETSC_EXTERN void PetscMin_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
348: {
349: PetscInt i,count = *cnt;
352: if (*datatype == MPIU_REAL) {
353: PetscReal *xin = (PetscReal*)in,*xout = (PetscReal*)out;
354: for (i=0; i<count; i++) xout[i] = PetscMin(xout[i],xin[i]);
355: }
356: #if defined(PETSC_HAVE_COMPLEX)
357: else if (*datatype == MPIU_COMPLEX) {
358: PetscComplex *xin = (PetscComplex*)in,*xout = (PetscComplex*)out;
359: for (i=0; i<count; i++) {
360: xout[i] = PetscRealPartComplex(xout[i])>PetscRealPartComplex(xin[i]) ? xin[i] : xout[i];
361: }
362: }
363: #endif
364: else {
365: (*PetscErrorPrintf)("Can only handle MPIU_REAL or MPIU_SCALAR data (i.e. double or complex) types");
366: MPI_Abort(MPI_COMM_WORLD,1);
367: }
368: PetscFunctionReturnVoid();
369: }
370: #endif
374: /*
375: Private routine to delete internal tag/name counter storage when a communicator is freed.
377: 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.
379: Note: this is declared extern "C" because it is passed to MPI_Keyval_create()
381: */
382: PETSC_EXTERN PetscMPIInt MPIAPI Petsc_DelCounter(MPI_Comm comm,PetscMPIInt keyval,void *count_val,void *extra_state)
383: {
387: PetscInfo1(0,"Deleting counter data in an MPI_Comm %ld\n",(long)comm);if (ierr) PetscFunctionReturn((PetscMPIInt)ierr);
388: PetscFree(count_val);if (ierr) PetscFunctionReturn((PetscMPIInt)ierr);
389: PetscFunctionReturn(MPI_SUCCESS);
390: }
394: /*
395: This does not actually free anything, it simply marks when a reference count to an internal or external MPI_Comm reaches zero and the
396: the external MPI_Comm drops its reference to the internal or external MPI_Comm
398: This is called by MPI, not by users. This is called when MPI_Comm_free() is called on the communicator.
400: Note: this is declared extern "C" because it is passed to MPI_Keyval_create()
402: */
403: PETSC_EXTERN PetscMPIInt MPIAPI Petsc_DelComm(MPI_Comm comm,PetscMPIInt keyval,void *attr_val,void *extra_state)
404: {
406: PetscMPIInt flg;
407: MPI_Comm icomm;
408: void *ptr;
411: MPI_Attr_get(comm,Petsc_InnerComm_keyval,&ptr,&flg);
412: if (flg) {
413: /* Use PetscMemcpy() because casting from pointer to integer of different size is not allowed with some compilers */
414: PetscMemcpy(&icomm,&ptr,sizeof(MPI_Comm));
415: MPI_Attr_get(icomm,Petsc_OuterComm_keyval,&ptr,&flg);
416: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"Inner MPI_Comm does not have expected reference to outer comm");
417: MPI_Attr_delete(icomm,Petsc_OuterComm_keyval);
418: PetscInfo1(0,"User MPI_Comm m %ld is being freed, removing reference from inner PETSc comm to this outer comm\n",(long)comm);if (ierr) PetscFunctionReturn((PetscMPIInt)ierr);
419: } else {
420: PetscInfo1(0,"Removing reference to PETSc communicator imbedded in a user MPI_Comm m %ld\n",(long)comm);if (ierr) PetscFunctionReturn((PetscMPIInt)ierr);
421: }
422: PetscFunctionReturn(MPI_SUCCESS);
423: }
425: #if defined(PETSC_USE_PETSC_MPI_EXTERNAL32)
426: #if !defined(PETSC_WORDS_BIGENDIAN)
427: PETSC_EXTERN PetscMPIInt PetscDataRep_extent_fn(MPI_Datatype,MPI_Aint*,void*);
428: PETSC_EXTERN PetscMPIInt PetscDataRep_read_conv_fn(void*, MPI_Datatype,PetscMPIInt,void*,MPI_Offset,void*);
429: PETSC_EXTERN PetscMPIInt PetscDataRep_write_conv_fn(void*, MPI_Datatype,PetscMPIInt,void*,MPI_Offset,void*);
430: #endif
431: #endif
433: int PetscGlobalArgc = 0;
434: char **PetscGlobalArgs = 0;
438: /*@C
439: PetscGetArgs - Allows you to access the raw command line arguments anywhere
440: after PetscInitialize() is called but before PetscFinalize().
442: Not Collective
444: Output Parameters:
445: + argc - count of number of command line arguments
446: - args - the command line arguments
448: Level: intermediate
450: Notes:
451: This is usually used to pass the command line arguments into other libraries
452: that are called internally deep in PETSc or the application.
454: The first argument contains the program name as is normal for C arguments.
456: Concepts: command line arguments
458: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArguments()
460: @*/
461: PetscErrorCode PetscGetArgs(int *argc,char ***args)
462: {
464: if (!PetscInitializeCalled && PetscFinalizeCalled) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"You must call after PetscInitialize() but before PetscFinalize()");
465: *argc = PetscGlobalArgc;
466: *args = PetscGlobalArgs;
467: return(0);
468: }
472: /*@C
473: PetscGetArguments - Allows you to access the command line arguments anywhere
474: after PetscInitialize() is called but before PetscFinalize().
476: Not Collective
478: Output Parameters:
479: . args - the command line arguments
481: Level: intermediate
483: Notes:
484: This does NOT start with the program name and IS null terminated (final arg is void)
486: Concepts: command line arguments
488: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscFreeArguments()
490: @*/
491: PetscErrorCode PetscGetArguments(char ***args)
492: {
493: PetscInt i,argc = PetscGlobalArgc;
497: if (!PetscInitializeCalled && PetscFinalizeCalled) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ORDER,"You must call after PetscInitialize() but before PetscFinalize()");
498: if (!argc) {*args = 0; return(0);}
499: PetscMalloc(argc*sizeof(char*),args);
500: for (i=0; i<argc-1; i++) {
501: PetscStrallocpy(PetscGlobalArgs[i+1],&(*args)[i]);
502: }
503: (*args)[argc-1] = 0;
504: return(0);
505: }
509: /*@C
510: PetscFreeArguments - Frees the memory obtained with PetscGetArguments()
512: Not Collective
514: Output Parameters:
515: . args - the command line arguments
517: Level: intermediate
519: Concepts: command line arguments
521: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscGetArguments()
523: @*/
524: PetscErrorCode PetscFreeArguments(char **args)
525: {
526: PetscInt i = 0;
530: if (!args) return(0);
531: while (args[i]) {
532: PetscFree(args[i]);
533: i++;
534: }
535: PetscFree(args);
536: return(0);
537: }
541: /*@C
542: PetscInitialize - Initializes the PETSc database and MPI.
543: PetscInitialize() calls MPI_Init() if that has yet to be called,
544: so this routine should always be called near the beginning of
545: your program -- usually the very first line!
547: Collective on MPI_COMM_WORLD or PETSC_COMM_WORLD if it has been set
549: Input Parameters:
550: + argc - count of number of command line arguments
551: . args - the command line arguments
552: . file - [optional] PETSc database file, also checks ~username/.petscrc and .petscrc use NULL to not check for
553: code specific file. Use -skip_petscrc in the code specific file to skip the .petscrc files
554: - help - [optional] Help message to print, use NULL for no message
556: If you wish PETSc code to run ONLY on a subcommunicator of MPI_COMM_WORLD, create that
557: communicator first and assign it to PETSC_COMM_WORLD BEFORE calling PetscInitialize(). Thus if you are running a
558: four process job and two processes will run PETSc and have PetscInitialize() and PetscFinalize() and two process will not,
559: then do this. If ALL processes in the job are using PetscInitialize() and PetscFinalize() then you don't need to do this, even
560: if different subcommunicators of the job are doing different things with PETSc.
562: Options Database Keys:
563: + -start_in_debugger [noxterm,dbx,xdb,gdb,...] - Starts program in debugger
564: . -on_error_attach_debugger [noxterm,dbx,xdb,gdb,...] - Starts debugger when error detected
565: . -on_error_emacs <machinename> causes emacsclient to jump to error file
566: . -on_error_abort calls abort() when error detected (no traceback)
567: . -on_error_mpiabort calls MPI_abort() when error detected
568: . -error_output_stderr prints error messages to stderr instead of the default stdout
569: . -error_output_none does not print the error messages (but handles errors in the same way as if this was not called)
570: . -debugger_nodes [node1,node2,...] - Indicates nodes to start in debugger
571: . -debugger_pause [sleeptime] (in seconds) - Pauses debugger
572: . -stop_for_debugger - Print message on how to attach debugger manually to
573: process and wait (-debugger_pause) seconds for attachment
574: . -malloc - Indicates use of PETSc error-checking malloc (on by default for debug version of libraries)
575: . -malloc no - Indicates not to use error-checking malloc
576: . -malloc_debug - check for memory corruption at EVERY malloc or free
577: . -malloc_test - like -malloc_dump -malloc_debug, but only active for debugging builds
578: . -fp_trap - Stops on floating point exceptions (Note that on the
579: IBM RS6000 this slows code by at least a factor of 10.)
580: . -no_signal_handler - Indicates not to trap error signals
581: . -shared_tmp - indicates /tmp directory is shared by all processors
582: . -not_shared_tmp - each processor has own /tmp
583: . -tmp - alternative name of /tmp directory
584: . -get_total_flops - returns total flops done by all processors
585: . -memory_info - Print memory usage at end of run
586: - -server <port> - start PETSc webserver (default port is 8080)
588: Options Database Keys for Profiling:
589: See the <a href="../../docs/manual.pdf#nameddest=Chapter 10 Profiling">profiling chapter of the users manual</a> for details.
590: + -info <optional filename> - Prints verbose information to the screen
591: . -info_exclude <null,vec,mat,pc,ksp,snes,ts> - Excludes some of the verbose messages
592: . -log_sync - Log the synchronization in scatters, inner products and norms
593: . -log_trace [filename] - Print traces of all PETSc calls to the screen (useful to determine where a program
594: hangs without running in the debugger). See PetscLogTraceBegin().
595: . -log_summary [filename] - Prints summary of flop and timing information to screen. If the filename is specified the
596: summary is written to the file. See PetscLogView().
597: . -log_summary_python [filename] - Prints data on of flop and timing usage to a file or screen. See PetscLogPrintSViewPython().
598: . -log_all [filename] - Logs extensive profiling information See PetscLogDump().
599: . -log [filename] - Logs basic profiline information See PetscLogDump().
600: - -log_mpe [filename] - Creates a logfile viewable by the utility Jumpshot (in MPICH distribution)
602: Only one of -log_trace, -log_summary, -log_all, -log, or -log_mpe may be used at a time
604: Environmental Variables:
605: + PETSC_TMP - alternative tmp directory
606: . PETSC_SHARED_TMP - tmp is shared by all processes
607: . PETSC_NOT_SHARED_TMP - each process has its own private tmp
608: . PETSC_VIEWER_SOCKET_PORT - socket number to use for socket viewer
609: - PETSC_VIEWER_SOCKET_MACHINE - machine to use for socket viewer to connect to
612: Level: beginner
614: Notes:
615: If for some reason you must call MPI_Init() separately, call
616: it before PetscInitialize().
618: Fortran Version:
619: In Fortran this routine has the format
620: $ call PetscInitialize(file,ierr)
622: + ierr - error return code
623: - file - [optional] PETSc database file, also checks ~username/.petscrc and .petscrc use NULL_CHARACTER to not check for
624: code specific file. Use -skip_petscrc in the code specific file to skip the .petscrc files
626: Important Fortran Note:
627: In Fortran, you MUST use NULL_CHARACTER to indicate a
628: null character string; you CANNOT just use NULL as
629: in the C version. See the <a href="../../docs/manual.pdf">users manual</a> for details.
631: If your main program is C but you call Fortran code that also uses PETSc you need to call PetscInitializeFortran() soon after
632: calling PetscInitialize().
634: Concepts: initializing PETSc
636: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscInitializeNoArguments()
638: @*/
639: PetscErrorCode PetscInitialize(int *argc,char ***args,const char file[],const char help[])
640: {
642: PetscMPIInt flag, size;
643: PetscInt nodesize;
644: PetscBool flg;
645: char hostname[256];
648: if (PetscInitializeCalled) return(0);
650: /* these must be initialized in a routine, not as a constant declaration*/
651: PETSC_STDOUT = stdout;
652: PETSC_STDERR = stderr;
654: PetscOptionsCreate();
656: /*
657: We initialize the program name here (before MPI_Init()) because MPICH has a bug in
658: it that it sets args[0] on all processors to be args[0] on the first processor.
659: */
660: if (argc && *argc) {
661: PetscSetProgramName(**args);
662: } else {
663: PetscSetProgramName("Unknown Name");
664: }
666: MPI_Initialized(&flag);
667: if (!flag) {
668: if (PETSC_COMM_WORLD != MPI_COMM_NULL) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"You cannot set PETSC_COMM_WORLD if you have not initialized MPI first");
669: #if defined(PETSC_HAVE_MPI_INIT_THREAD)
670: {
671: PetscMPIInt provided;
672: MPI_Init_thread(argc,args,MPI_THREAD_FUNNELED,&provided);
673: }
674: #else
675: MPI_Init(argc,args);
676: #endif
677: PetscBeganMPI = PETSC_TRUE;
678: }
679: if (argc && args) {
680: PetscGlobalArgc = *argc;
681: PetscGlobalArgs = *args;
682: }
683: PetscFinalizeCalled = PETSC_FALSE;
685: if (PETSC_COMM_WORLD == MPI_COMM_NULL) PETSC_COMM_WORLD = MPI_COMM_WORLD;
686: MPI_Comm_set_errhandler(PETSC_COMM_WORLD,MPI_ERRORS_RETURN);
688: /* Done after init due to a bug in MPICH-GM? */
689: PetscErrorPrintfInitialize();
691: MPI_Comm_rank(MPI_COMM_WORLD,&PetscGlobalRank);
692: MPI_Comm_size(MPI_COMM_WORLD,&PetscGlobalSize);
694: MPIU_BOOL = MPI_INT;
695: MPIU_ENUM = MPI_INT;
697: /*
698: Initialized the global complex variable; this is because with
699: shared libraries the constructors for global variables
700: are not called; at least on IRIX.
701: */
702: #if defined(PETSC_HAVE_COMPLEX)
703: {
704: #if defined(PETSC_CLANGUAGE_CXX)
705: PetscComplex ic(0.0,1.0);
706: PETSC_i = ic;
707: #elif defined(PETSC_CLANGUAGE_C)
708: PETSC_i = _Complex_I;
709: #endif
710: }
712: #if !defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)
713: MPI_Type_contiguous(2,MPI_DOUBLE,&MPIU_C_DOUBLE_COMPLEX);
714: MPI_Type_commit(&MPIU_C_DOUBLE_COMPLEX);
715: MPI_Type_contiguous(2,MPI_FLOAT,&MPIU_C_COMPLEX);
716: MPI_Type_commit(&MPIU_C_COMPLEX);
717: #endif
718: #endif /* PETSC_HAVE_COMPLEX */
720: /*
721: Create the PETSc MPI reduction operator that sums of the first
722: half of the entries and maxes the second half.
723: */
724: MPI_Op_create(PetscMaxSum_Local,1,&PetscMaxSum_Op);
726: #if defined(PETSC_USE_REAL___FLOAT128)
727: MPI_Type_contiguous(2,MPI_DOUBLE,&MPIU___FLOAT128);
728: MPI_Type_commit(&MPIU___FLOAT128);
729: #if defined(PETSC_HAVE_COMPLEX)
730: MPI_Type_contiguous(4,MPI_DOUBLE,&MPIU___COMPLEX128);
731: MPI_Type_commit(&MPIU___COMPLEX128);
732: #endif
733: MPI_Op_create(PetscMax_Local,1,&MPIU_MAX);
734: MPI_Op_create(PetscMin_Local,1,&MPIU_MIN);
735: #endif
737: #if (defined(PETSC_HAVE_COMPLEX) && !defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)) || defined(PETSC_USE_REAL___FLOAT128)
738: MPI_Op_create(PetscSum_Local,1,&MPIU_SUM);
739: #endif
741: MPI_Type_contiguous(2,MPIU_SCALAR,&MPIU_2SCALAR);
742: MPI_Type_commit(&MPIU_2SCALAR);
743: MPI_Op_create(PetscADMax_Local,1,&PetscADMax_Op);
744: MPI_Op_create(PetscADMin_Local,1,&PetscADMin_Op);
746: #if defined(PETSC_USE_64BIT_INDICES) || !defined(MPI_2INT)
747: MPI_Type_contiguous(2,MPIU_INT,&MPIU_2INT);
748: MPI_Type_commit(&MPIU_2INT);
749: #endif
751: /*
752: Attributes to be set on PETSc communicators
753: */
754: MPI_Keyval_create(MPI_NULL_COPY_FN,Petsc_DelCounter,&Petsc_Counter_keyval,(void*)0);
755: MPI_Keyval_create(MPI_NULL_COPY_FN,Petsc_DelComm,&Petsc_InnerComm_keyval,(void*)0);
756: MPI_Keyval_create(MPI_NULL_COPY_FN,Petsc_DelComm,&Petsc_OuterComm_keyval,(void*)0);
758: /*
759: Build the options database
760: */
761: PetscOptionsInsert(argc,args,file);
764: /*
765: Print main application help message
766: */
767: PetscOptionsHasName(NULL,"-help",&flg);
768: if (help && flg) {
769: PetscPrintf(PETSC_COMM_WORLD,help);
770: }
771: PetscOptionsCheckInitial_Private();
773: /* SHOULD PUT IN GUARDS: Make sure logging is initialized, even if we do not print it out */
774: #if defined(PETSC_USE_LOG)
775: PetscLogBegin_Private();
776: #endif
778: /*
779: Load the dynamic libraries (on machines that support them), this registers all
780: the solvers etc. (On non-dynamic machines this initializes the PetscDraw and PetscViewer classes)
781: */
782: PetscInitialize_DynamicLibraries();
784: MPI_Comm_size(PETSC_COMM_WORLD,&size);
785: PetscInfo1(0,"PETSc successfully started: number of processors = %d\n",size);
786: PetscGetHostName(hostname,256);
787: PetscInfo1(0,"Running on machine: %s\n",hostname);
789: PetscOptionsCheckInitial_Components();
790: /* Check the options database for options related to the options database itself */
791: PetscOptionsSetFromOptions();
793: #if defined(PETSC_USE_PETSC_MPI_EXTERNAL32)
794: /*
795: Tell MPI about our own data representation converter, this would/should be used if extern32 is not supported by the MPI
797: Currently not used because it is not supported by MPICH.
798: */
799: #if !defined(PETSC_WORDS_BIGENDIAN)
800: MPI_Register_datarep((char*)"petsc",PetscDataRep_read_conv_fn,PetscDataRep_write_conv_fn,PetscDataRep_extent_fn,NULL);
801: #endif
802: #endif
804: PetscOptionsGetInt(NULL,"-hmpi_spawn_size",&nodesize,&flg);
805: if (flg) {
806: #if defined(PETSC_HAVE_MPI_COMM_SPAWN)
807: PetscHMPISpawn((PetscMPIInt) nodesize); /* worker nodes never return from here; they go directly to PetscEnd() */
808: #else
809: SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"PETSc built without MPI 2 (MPI_Comm_spawn) support, use -hmpi_merge_size instead");
810: #endif
811: } else {
812: PetscOptionsGetInt(NULL,"-hmpi_merge_size",&nodesize,&flg);
813: if (flg) {
814: PetscHMPIMerge((PetscMPIInt) nodesize,NULL,NULL);
815: if (PetscHMPIWorker) { /* if worker then never enter user code */
816: PetscInitializeCalled = PETSC_TRUE;
817: PetscEnd();
818: }
819: }
820: }
822: #if defined(PETSC_HAVE_CUDA)
823: {
824: PetscMPIInt p;
825: for (p = 0; p < PetscGlobalSize; ++p) {
826: if (p == PetscGlobalRank) cublasInit();
827: MPI_Barrier(PETSC_COMM_WORLD);
828: }
829: }
830: #endif
832: PetscOptionsHasName(NULL,"-python",&flg);
833: if (flg) {
834: PetscInitializeCalled = PETSC_TRUE;
835: PetscPythonInitialize(NULL,NULL);
836: }
838: PetscThreadCommInitializePackage();
840: /*
841: Setup building of stack frames for all function calls
842: */
843: #if defined(PETSC_USE_DEBUG)
844: PetscThreadLocalRegister((PetscThreadKey*)&petscstack); /* Creates petscstack_key if needed */
845: PetscStackCreate();
846: #endif
848: #if defined(PETSC_SERIALIZE_FUNCTIONS)
849: PetscFPTCreate(10000);
850: #endif
852: /*
853: Once we are completedly initialized then we can set this variables
854: */
855: PetscInitializeCalled = PETSC_TRUE;
856: return(0);
857: }
859: extern PetscObject *PetscObjects;
860: extern PetscInt PetscObjectsCounts, PetscObjectsMaxCounts;
864: /*@C
865: PetscFinalize - Checks for options to be called at the conclusion
866: of the program. MPI_Finalize() is called only if the user had not
867: called MPI_Init() before calling PetscInitialize().
869: Collective on PETSC_COMM_WORLD
871: Options Database Keys:
872: + -options_table - Calls PetscOptionsView()
873: . -options_left - Prints unused options that remain in the database
874: . -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
875: . -mpidump - Calls PetscMPIDump()
876: . -malloc_dump - Calls PetscMallocDump()
877: . -malloc_info - Prints total memory usage
878: - -malloc_log - Prints summary of memory usage
880: Level: beginner
882: Note:
883: See PetscInitialize() for more general runtime options.
885: .seealso: PetscInitialize(), PetscOptionsView(), PetscMallocDump(), PetscMPIDump(), PetscEnd()
886: @*/
887: PetscErrorCode PetscFinalize(void)
888: {
890: PetscMPIInt rank;
891: PetscInt nopt;
892: PetscBool flg1 = PETSC_FALSE,flg2 = PETSC_FALSE,flg3 = PETSC_FALSE;
893: #if defined(PETSC_HAVE_AMS)
894: PetscBool flg = PETSC_FALSE;
895: #endif
896: #if defined(PETSC_USE_LOG)
897: char mname[PETSC_MAX_PATH_LEN];
898: #endif
901: if (!PetscInitializeCalled) {
902: printf("PetscInitialize() must be called before PetscFinalize()\n");
903: PetscFunctionReturn(PETSC_ERR_ARG_WRONGSTATE);
904: }
905: PetscInfo(NULL,"PetscFinalize() called\n");
907: #if defined(PETSC_SERIALIZE_FUNCTIONS)
908: PetscFPTDestroy();
909: #endif
912: #if defined(PETSC_HAVE_AMS)
913: PetscOptionsGetBool(NULL,"-options_gui",&flg,NULL);
914: if (flg) {
915: PetscOptionsAMSDestroy();
916: }
917: #endif
919: #if defined(PETSC_HAVE_SERVER)
920: flg1 = PETSC_FALSE;
921: PetscOptionsGetBool(NULL,"-server",&flg1,NULL);
922: if (flg1) {
923: /* this is a crude hack, but better than nothing */
924: PetscPOpen(PETSC_COMM_WORLD,NULL,"pkill -9 petscwebserver","r",NULL);
925: }
926: #endif
928: PetscHMPIFinalize();
930: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
931: PetscOptionsGetBool(NULL,"-malloc_info",&flg2,NULL);
932: if (!flg2) {
933: flg2 = PETSC_FALSE;
934: PetscOptionsGetBool(NULL,"-memory_info",&flg2,NULL);
935: }
936: if (flg2) {
937: PetscMemoryShowUsage(PETSC_VIEWER_STDOUT_WORLD,"Summary of Memory Usage in PETSc\n");
938: }
940: #if defined(PETSC_USE_LOG)
941: flg1 = PETSC_FALSE;
942: PetscOptionsGetBool(NULL,"-get_total_flops",&flg1,NULL);
943: if (flg1) {
944: PetscLogDouble flops = 0;
945: MPI_Reduce(&petsc_TotalFlops,&flops,1,MPI_DOUBLE,MPI_SUM,0,PETSC_COMM_WORLD);
946: PetscPrintf(PETSC_COMM_WORLD,"Total flops over all processors %g\n",flops);
947: }
948: #endif
951: #if defined(PETSC_USE_LOG)
952: #if defined(PETSC_HAVE_MPE)
953: mname[0] = 0;
955: PetscOptionsGetString(NULL,"-log_mpe",mname,PETSC_MAX_PATH_LEN,&flg1);
956: if (flg1) {
957: if (mname[0]) {PetscLogMPEDump(mname);}
958: else {PetscLogMPEDump(0);}
959: }
960: #endif
961: mname[0] = 0;
963: PetscOptionsGetString(NULL,"-log_summary",mname,PETSC_MAX_PATH_LEN,&flg1);
964: if (flg1) {
965: PetscViewer viewer;
966: if (mname[0]) {
967: PetscViewerASCIIOpen(PETSC_COMM_WORLD,mname,&viewer);
968: PetscLogView(viewer);
969: PetscViewerDestroy(&viewer);
970: } else {
971: viewer = PETSC_VIEWER_STDOUT_WORLD;
972: PetscLogView(viewer);
973: }
974: }
976: mname[0] = 0;
978: PetscOptionsGetString(NULL,"-log_summary_python",mname,PETSC_MAX_PATH_LEN,&flg1);
979: if (flg1) {
980: PetscViewer viewer;
981: if (mname[0]) {
982: PetscViewerASCIIOpen(PETSC_COMM_WORLD,mname,&viewer);
983: PetscLogViewPython(viewer);
984: PetscViewerDestroy(&viewer);
985: } else {
986: viewer = PETSC_VIEWER_STDOUT_WORLD;
987: PetscLogViewPython(viewer);
988: }
989: }
991: PetscOptionsGetString(NULL,"-log_detailed",mname,PETSC_MAX_PATH_LEN,&flg1);
992: if (flg1) {
993: if (mname[0]) {PetscLogPrintDetailed(PETSC_COMM_WORLD,mname);}
994: else {PetscLogPrintDetailed(PETSC_COMM_WORLD,0);}
995: }
997: mname[0] = 0;
999: PetscOptionsGetString(NULL,"-log_all",mname,PETSC_MAX_PATH_LEN,&flg1);
1000: PetscOptionsGetString(NULL,"-log",mname,PETSC_MAX_PATH_LEN,&flg2);
1001: if (flg1 || flg2) {
1002: if (mname[0]) PetscLogDump(mname);
1003: else PetscLogDump(0);
1004: }
1005: #endif
1007: /*
1008: Free all objects registered with PetscObjectRegisterDestroy() such as PETSC_VIEWER_XXX_().
1009: */
1010: PetscObjectRegisterDestroyAll();
1012: PetscStackDestroy();
1014: flg1 = PETSC_FALSE;
1015: PetscOptionsGetBool(NULL,"-no_signal_handler",&flg1,NULL);
1016: if (!flg1) { PetscPopSignalHandler();}
1017: flg1 = PETSC_FALSE;
1018: PetscOptionsGetBool(NULL,"-mpidump",&flg1,NULL);
1019: if (flg1) {
1020: PetscMPIDump(stdout);
1021: }
1022: flg1 = PETSC_FALSE;
1023: flg2 = PETSC_FALSE;
1024: /* preemptive call to avoid listing this option in options table as unused */
1025: PetscOptionsHasName(NULL,"-malloc_dump",&flg1);
1026: PetscOptionsHasName(NULL,"-objects_dump",&flg1);
1027: PetscOptionsGetBool(NULL,"-options_table",&flg2,NULL);
1029: if (flg2) {
1030: PetscViewer viewer;
1031: PetscViewerASCIIGetStdout(PETSC_COMM_WORLD,&viewer);
1032: PetscOptionsView(viewer);
1033: PetscViewerDestroy(&viewer);
1034: }
1036: /* to prevent PETSc -options_left from warning */
1037: PetscOptionsHasName(NULL,"-nox",&flg1);
1038: PetscOptionsHasName(NULL,"-nox_warning",&flg1);
1040: if (!PetscHMPIWorker) { /* worker processes skip this because they do not usually process options */
1041: flg3 = PETSC_FALSE; /* default value is required */
1042: PetscOptionsGetBool(NULL,"-options_left",&flg3,&flg1);
1043: PetscOptionsAllUsed(&nopt);
1044: if (flg3) {
1045: if (!flg2) { /* have not yet printed the options */
1046: PetscViewer viewer;
1047: PetscViewerASCIIGetStdout(PETSC_COMM_WORLD,&viewer);
1048: PetscOptionsView(viewer);
1049: PetscViewerDestroy(&viewer);
1050: }
1051: if (!nopt) {
1052: PetscPrintf(PETSC_COMM_WORLD,"There are no unused options.\n");
1053: } else if (nopt == 1) {
1054: PetscPrintf(PETSC_COMM_WORLD,"There is one unused database option. It is:\n");
1055: } else {
1056: PetscPrintf(PETSC_COMM_WORLD,"There are %D unused database options. They are:\n",nopt);
1057: }
1058: }
1059: #if defined(PETSC_USE_DEBUG)
1060: if (nopt && !flg3 && !flg1) {
1061: PetscPrintf(PETSC_COMM_WORLD,"WARNING! There are options you set that were not used!\n");
1062: PetscPrintf(PETSC_COMM_WORLD,"WARNING! could be spelling mistake, etc!\n");
1063: PetscOptionsLeft();
1064: } else if (nopt && flg3) {
1065: #else
1066: if (nopt && flg3) {
1067: #endif
1068: PetscOptionsLeft();
1069: }
1070: }
1072: {
1073: PetscThreadComm tcomm_world;
1074: PetscGetThreadCommWorld(&tcomm_world);
1075: /* Free global thread communicator */
1076: PetscThreadCommDestroy(&tcomm_world);
1077: }
1079: /*
1080: List all objects the user may have forgot to free
1081: */
1082: PetscOptionsHasName(NULL,"-objects_dump",&flg1);
1083: if (flg1) {
1084: MPI_Comm local_comm;
1085: char string[64];
1087: PetscOptionsGetString(NULL,"-objects_dump",string,64,NULL);
1088: MPI_Comm_dup(MPI_COMM_WORLD,&local_comm);
1089: PetscSequentialPhaseBegin_Private(local_comm,1);
1090: PetscObjectsDump(stdout,(string[0] == 'a') ? PETSC_TRUE : PETSC_FALSE);
1091: PetscSequentialPhaseEnd_Private(local_comm,1);
1092: MPI_Comm_free(&local_comm);
1093: }
1094: PetscObjectsCounts = 0;
1095: PetscObjectsMaxCounts = 0;
1097: PetscFree(PetscObjects);
1099: #if defined(PETSC_USE_LOG)
1100: PetscLogDestroy();
1101: #endif
1103: /*
1104: Destroy any packages that registered a finalize
1105: */
1106: PetscRegisterFinalizeAll();
1108: /*
1109: Destroy all the function registration lists created
1110: */
1111: PetscFinalize_DynamicLibraries();
1113: /*
1114: Print PetscFunctionLists that have not been properly freed
1116: PetscFunctionListPrintAll();
1117: */
1119: if (petsc_history) {
1120: PetscCloseHistoryFile(&petsc_history);
1121: petsc_history = 0;
1122: }
1124: PetscInfoAllow(PETSC_FALSE,NULL);
1126: {
1127: char fname[PETSC_MAX_PATH_LEN];
1128: FILE *fd;
1129: int err;
1131: fname[0] = 0;
1133: PetscOptionsGetString(NULL,"-malloc_dump",fname,250,&flg1);
1134: flg2 = PETSC_FALSE;
1135: PetscOptionsGetBool(NULL,"-malloc_test",&flg2,NULL);
1136: #if defined(PETSC_USE_DEBUG)
1137: if (PETSC_RUNNING_ON_VALGRIND) flg2 = PETSC_FALSE;
1138: #else
1139: flg2 = PETSC_FALSE; /* Skip reporting for optimized builds regardless of -malloc_test */
1140: #endif
1141: if (flg1 && fname[0]) {
1142: char sname[PETSC_MAX_PATH_LEN];
1144: sprintf(sname,"%s_%d",fname,rank);
1145: fd = fopen(sname,"w"); if (!fd) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Cannot open log file: %s",sname);
1146: PetscMallocDump(fd);
1147: err = fclose(fd);
1148: if (err) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fclose() failed on file");
1149: } else if (flg1 || flg2) {
1150: MPI_Comm local_comm;
1152: MPI_Comm_dup(MPI_COMM_WORLD,&local_comm);
1153: PetscSequentialPhaseBegin_Private(local_comm,1);
1154: PetscMallocDump(stdout);
1155: PetscSequentialPhaseEnd_Private(local_comm,1);
1156: MPI_Comm_free(&local_comm);
1157: }
1158: }
1160: {
1161: char fname[PETSC_MAX_PATH_LEN];
1162: FILE *fd = NULL;
1164: fname[0] = 0;
1166: PetscOptionsGetString(NULL,"-malloc_log",fname,250,&flg1);
1167: PetscOptionsHasName(NULL,"-malloc_log_threshold",&flg2);
1168: if (flg1 && fname[0]) {
1169: int err;
1171: if (!rank) {
1172: fd = fopen(fname,"w");
1173: if (!fd) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN,"Cannot open log file: %s",fname);
1174: }
1175: PetscMallocDumpLog(fd);
1176: if (fd) {
1177: err = fclose(fd);
1178: if (err) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fclose() failed on file");
1179: }
1180: } else if (flg1 || flg2) {
1181: PetscMallocDumpLog(stdout);
1182: }
1183: }
1184: /* Can be destroyed only after all the options are used */
1185: PetscOptionsDestroy();
1187: PetscGlobalArgc = 0;
1188: PetscGlobalArgs = 0;
1190: #if defined(PETSC_USE_REAL___FLOAT128)
1191: MPI_Type_free(&MPIU___FLOAT128);
1192: #if defined(PETSC_HAVE_COMPLEX)
1193: MPI_Type_free(&MPIU___COMPLEX128);
1194: #endif
1195: MPI_Op_free(&MPIU_MAX);
1196: MPI_Op_free(&MPIU_MIN);
1197: #endif
1199: #if defined(PETSC_HAVE_COMPLEX)
1200: #if !defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)
1201: MPI_Type_free(&MPIU_C_DOUBLE_COMPLEX);
1202: MPI_Type_free(&MPIU_C_COMPLEX);
1203: #endif
1204: #endif
1206: #if (defined(PETSC_HAVE_COMPLEX) && !defined(PETSC_HAVE_MPI_C_DOUBLE_COMPLEX)) || defined(PETSC_USE_REAL___FLOAT128)
1207: MPI_Op_free(&MPIU_SUM);
1208: #endif
1210: MPI_Type_free(&MPIU_2SCALAR);
1211: #if defined(PETSC_USE_64BIT_INDICES) || !defined(MPI_2INT)
1212: MPI_Type_free(&MPIU_2INT);
1213: #endif
1214: MPI_Op_free(&PetscMaxSum_Op);
1215: MPI_Op_free(&PetscADMax_Op);
1216: MPI_Op_free(&PetscADMin_Op);
1218: /*
1219: Destroy any known inner MPI_Comm's and attributes pointing to them
1220: Note this will not destroy any new communicators the user has created.
1222: If all PETSc objects were not destroyed those left over objects will have hanging references to
1223: the MPI_Comms that were freed; but that is ok because those PETSc objects will never be used again
1224: */
1225: {
1226: PetscCommCounter *counter;
1227: PetscMPIInt flg;
1228: MPI_Comm icomm;
1229: void *ptr;
1230: MPI_Attr_get(PETSC_COMM_SELF,Petsc_InnerComm_keyval,&ptr,&flg);
1231: if (flg) {
1232: /* Use PetscMemcpy() because casting from pointer to integer of different size is not allowed with some compilers */
1233: PetscMemcpy(&icomm,&ptr,sizeof(MPI_Comm));
1234: MPI_Attr_get(icomm,Petsc_Counter_keyval,&counter,&flg);
1235: if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"Inner MPI_Comm does not have expected tag/name counter, problem with corrupted memory");
1237: MPI_Attr_delete(PETSC_COMM_SELF,Petsc_InnerComm_keyval);
1238: MPI_Attr_delete(icomm,Petsc_Counter_keyval);
1239: MPI_Comm_free(&icomm);
1240: }
1241: MPI_Attr_get(PETSC_COMM_WORLD,Petsc_InnerComm_keyval,&ptr,&flg);
1242: if (flg) {
1243: /* Use PetscMemcpy() because casting from pointer to integer of different size is not allowed with some compilers */
1244: PetscMemcpy(&icomm,&ptr,sizeof(MPI_Comm));
1245: MPI_Attr_get(icomm,Petsc_Counter_keyval,&counter,&flg);
1246: if (!flg) SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_ARG_CORRUPT,"Inner MPI_Comm does not have expected tag/name counter, problem with corrupted memory");
1248: MPI_Attr_delete(PETSC_COMM_WORLD,Petsc_InnerComm_keyval);
1249: MPI_Attr_delete(icomm,Petsc_Counter_keyval);
1250: MPI_Comm_free(&icomm);
1251: }
1252: }
1254: MPI_Keyval_free(&Petsc_Counter_keyval);
1255: MPI_Keyval_free(&Petsc_InnerComm_keyval);
1256: MPI_Keyval_free(&Petsc_OuterComm_keyval);
1258: #if defined(PETSC_HAVE_CUDA)
1259: {
1260: PetscInt p;
1261: for (p = 0; p < PetscGlobalSize; ++p) {
1262: if (p == PetscGlobalRank) cublasShutdown();
1263: MPI_Barrier(PETSC_COMM_WORLD);
1264: }
1265: }
1266: #endif
1268: if (PetscBeganMPI) {
1269: #if defined(PETSC_HAVE_MPI_FINALIZED)
1270: PetscMPIInt flag;
1271: MPI_Finalized(&flag);
1272: if (flag) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"MPI_Finalize() has already been called, even though MPI_Init() was called by PetscInitialize()");
1273: #endif
1274: MPI_Finalize();
1275: }
1276: /*
1278: Note: In certain cases PETSC_COMM_WORLD is never MPI_Comm_free()ed because
1279: the communicator has some outstanding requests on it. Specifically if the
1280: flag PETSC_HAVE_BROKEN_REQUEST_FREE is set (for IBM MPI implementation). See
1281: src/vec/utils/vpscat.c. Due to this the memory allocated in PetscCommDuplicate()
1282: is never freed as it should be. Thus one may obtain messages of the form
1283: [ 1] 8 bytes PetscCommDuplicate() line 645 in src/sys/mpiu.c indicating the
1284: memory was not freed.
1286: */
1287: PetscMallocClear();
1289: PetscInitializeCalled = PETSC_FALSE;
1290: PetscFinalizeCalled = PETSC_TRUE;
1291: PetscFunctionReturn(ierr);
1292: }
1294: #if defined(PETSC_MISSING_LAPACK_lsame_)
1295: PETSC_EXTERN int lsame_(char *a,char *b)
1296: {
1297: if (*a == *b) return 1;
1298: if (*a + 32 == *b) return 1;
1299: if (*a - 32 == *b) return 1;
1300: return 0;
1301: }
1302: #endif
1304: #if defined(PETSC_MISSING_LAPACK_lsame)
1305: PETSC_EXTERN int lsame(char *a,char *b)
1306: {
1307: if (*a == *b) return 1;
1308: if (*a + 32 == *b) return 1;
1309: if (*a - 32 == *b) return 1;
1310: return 0;
1311: }
1312: #endif