Actual source code: pinit.c

  1: #define PETSC_DLL
  2: /*
  3:    This file defines the initialization of PETSc, including PetscInitialize()
  4: */

 6:  #include petsc.h
 7:  #include petscsys.h

  9: #if defined(PETSC_USE_LOG)
 10: EXTERN PetscErrorCode PetscLogBegin_Private(void);
 11: #endif

 14: /* -----------------------------------------------------------------------------------------*/


 18: EXTERN PetscErrorCode PetscInitialize_DynamicLibraries(void);
 19: EXTERN PetscErrorCode PetscFinalize_DynamicLibraries(void);
 20: EXTERN PetscErrorCode PetscFListDestroyAll(void);
 21: EXTERN PetscErrorCode PetscSequentialPhaseBegin_Private(MPI_Comm,int);
 22: EXTERN PetscErrorCode PetscSequentialPhaseEnd_Private(MPI_Comm,int);
 23: EXTERN PetscErrorCode PetscLogCloseHistoryFile(FILE **);
 24: EXTERN PetscErrorCode PetscOptionsHelpDestroyList(void);

 26: /* this is used by the _, __, and ___ macros (see include/petscerror.h) */
 27: PetscErrorCode __g0;

 29: /* user may set this BEFORE calling PetscInitialize() */
 30: MPI_Comm PETSC_COMM_WORLD = 0;

 32: /*
 33:      Declare and set all the string names of the PETSc enums
 34: */
 35: const char *PetscTruths[]    = {"FALSE","TRUE","PetscTruth","PETSC_",0};
 36: const char *PetscDataTypes[] = {"INT", "DOUBLE", "COMPLEX",
 37:                                 "LONG","SHORT",  "FLOAT",
 38:                                 "CHAR","LOGICAL","ENUM","TRUTH","LONGDOUBLE","PetscDataType","PETSC_",0};

 40: PetscTruth PetscPreLoadingUsed = PETSC_FALSE;
 41: PetscTruth PetscPreLoadingOn   = PETSC_FALSE;

 43: /*
 44:        Checks the options database for initializations related to the 
 45:     PETSc components
 46: */
 49: PetscErrorCode  PetscOptionsCheckInitial_Components(void)
 50: {
 51:   PetscTruth flg1;

 55:   PetscOptionsHasName(PETSC_NULL,"-help",&flg1);
 56:   if (flg1) {
 57: #if defined (PETSC_USE_LOG)
 58:     MPI_Comm   comm = PETSC_COMM_WORLD;
 59:     (*PetscHelpPrintf)(comm,"------Additional PETSc component options--------\n");
 60:     (*PetscHelpPrintf)(comm," -log_summary_exclude: <vec,mat,pc.ksp,snes>\n");
 61:     (*PetscHelpPrintf)(comm," -info_exclude: <null,vec,mat,pc,ksp,snes,ts>\n");
 62:     (*PetscHelpPrintf)(comm,"-----------------------------------------------\n");
 63: #endif
 64:   }
 65:   return(0);
 66: }

 70: /*@C
 71:       PetscInitializeNoArguments - Calls PetscInitialize() from C/C++ without
 72:         the command line arguments.

 74:    Collective
 75:   
 76:    Level: advanced

 78: .seealso: PetscInitialize(), PetscInitializeFortran()
 79: @*/
 80: PetscErrorCode  PetscInitializeNoArguments(void)
 81: {
 83:   int            argc = 0;
 84:   char           **args = 0;

 87:   PetscInitialize(&argc,&args,PETSC_NULL,PETSC_NULL);
 88:   PetscFunctionReturn(ierr);
 89: }

 93: /*@
 94:       PetscInitialized - Determine whether PETSc is initialized.
 95:   
 96: 7   Level: beginner

 98: .seealso: PetscInitialize(), PetscInitializeNoArguments(), PetscInitializeFortran()
 99: @*/
100: PetscErrorCode  PetscInitialized(PetscTruth *isInitialized)
101: {
104:   *isInitialized = PetscInitializeCalled;
105:   return(0);
106: }

110: /*@
111:       PetscFinalized - Determine whether PetscFinalize() has been called yet
112:   
113:    Level: developer

115: .seealso: PetscInitialize(), PetscInitializeNoArguments(), PetscInitializeFortran()
116: @*/
117: PetscErrorCode  PetscFinalized(PetscTruth *isFinalized)
118: {
121:   *isFinalized = PetscFinalizeCalled;
122:   return(0);
123: }

125: EXTERN PetscErrorCode        PetscOptionsCheckInitial_Private(void);

128: /*
129:        This function is the MPI reduction operation used to compute the sum of the 
130:    first half of the datatype and the max of the second half.
131: */
132: MPI_Op PetscMaxSum_Op = 0;

137: void  MPIAPI PetscMaxSum_Local(void *in,void *out,int *cnt,MPI_Datatype *datatype)
138: {
139:   PetscInt *xin = (PetscInt*)in,*xout = (PetscInt*)out,i,count = *cnt;

142:   if (*datatype != MPIU_2INT) {
143:     (*PetscErrorPrintf)("Can only handle MPIU_2INT data types");
144:     MPI_Abort(MPI_COMM_WORLD,1);
145:   }

147:   for (i=0; i<count; i++) {
148:     xout[2*i]    = PetscMax(xout[2*i],xin[2*i]);
149:     xout[2*i+1] += xin[2*i+1];
150:   }
151:   PetscStackPop;
152:   return;
153: }

156: /*
157:     Returns the max of the first entry owned by this processor and the
158: sum of the second entry.

160:     The reason nprocs[2*i] contains lengths nprocs[2*i+1] contains flag of 1 if length is nonzero
161: is so that the PetscMaxSum_Op() can set TWO values, if we passed in only nprocs[i] with lengths
162: there would be no place to store the both needed results.
163: */
166: PetscErrorCode  PetscMaxSum(MPI_Comm comm,const PetscInt nprocs[],PetscInt *max,PetscInt *sum)
167: {
168:   PetscMPIInt    size,rank;
169:   PetscInt       *work;

173:   MPI_Comm_size(comm,&size);
174:   MPI_Comm_rank(comm,&rank);
175:   PetscMalloc(2*size*sizeof(PetscInt),&work);
176:   MPI_Allreduce((void*)nprocs,work,size,MPIU_2INT,PetscMaxSum_Op,comm);
177:   *max   = work[2*rank];
178:   *sum   = work[2*rank+1];
179:   PetscFree(work);
180:   return(0);
181: }

183: /* ----------------------------------------------------------------------------*/
184: MPI_Op  PetscADMax_Op = 0;

189: void  MPIAPI PetscADMax_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
190: {
191:   PetscScalar *xin = (PetscScalar *)in,*xout = (PetscScalar*)out;
192:   PetscInt    i,count = *cnt;

195:   if (*datatype != MPIU_2SCALAR) {
196:     (*PetscErrorPrintf)("Can only handle MPIU_2SCALAR data (i.e. double or complex) types");
197:     MPI_Abort(MPI_COMM_WORLD,1);
198:   }

200:   for (i=0; i<count; i++) {
201:     if (PetscRealPart(xout[2*i]) < PetscRealPart(xin[2*i])) {
202:       xout[2*i]   = xin[2*i];
203:       xout[2*i+1] = xin[2*i+1];
204:     }
205:   }

207:   PetscStackPop;
208:   return;
209: }

212: MPI_Op  PetscADMin_Op = 0;

217: void  MPIAPI PetscADMin_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
218: {
219:   PetscScalar *xin = (PetscScalar *)in,*xout = (PetscScalar*)out;
220:   PetscInt    i,count = *cnt;

223:   if (*datatype != MPIU_2SCALAR) {
224:     (*PetscErrorPrintf)("Can only handle MPIU_2SCALAR data (i.e. double or complex) types");
225:     MPI_Abort(MPI_COMM_WORLD,1);
226:   }

228:   for (i=0; i<count; i++) {
229:     if (PetscRealPart(xout[2*i]) > PetscRealPart(xin[2*i])) {
230:       xout[2*i]   = xin[2*i];
231:       xout[2*i+1] = xin[2*i+1];
232:     }
233:   }

235:   PetscStackPop;
236:   return;
237: }
239: /* ---------------------------------------------------------------------------------------*/

241: #if defined(PETSC_USE_COMPLEX)
242: MPI_Op PetscSum_Op = 0;

247: void  PetscSum_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
248: {
249:   PetscScalar *xin = (PetscScalar *)in,*xout = (PetscScalar*)out;
250:   PetscInt    i,count = *cnt;

253:   if (*datatype != MPIU_SCALAR) {
254:     (*PetscErrorPrintf)("Can only handle MPIU_SCALAR data (i.e. double or complex) types");
255:     MPI_Abort(MPI_COMM_WORLD,1);
256:   }

258:   for (i=0; i<count; i++) {
259:     xout[i] += xin[i];
260:   }

262:   PetscStackPop;
263:   return;
264: }
266: #endif

268: #if defined(PETSC_USE_PETSC_MPI_EXTERNAL32)
269: #if !defined(PETSC_WORDS_BIGENDIAN)
275: #endif
276: #endif

278: int  PetscGlobalArgc   = 0;
279: char **PetscGlobalArgs = 0;

283: /*@C
284:    PetscGetArgs - Allows you to access the raw command line arguments anywhere
285:      after PetscInitialize() is called but before PetscFinalize().

287:    Not Collective

289:    Output Parameters:
290: +  argc - count of number of command line arguments
291: -  args - the command line arguments

293:    Level: intermediate

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

299:    Concepts: command line arguments
300:    
301: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArguments()

303: @*/
304: PetscErrorCode  PetscGetArgs(int *argc,char ***args)
305: {
307:   if (!PetscInitializeCalled && PetscFinalizeCalled) {
308:     SETERRQ(PETSC_ERR_ORDER,"You must call after PetscInitialize() but before PetscFinalize()");
309:   }
310:   *argc = PetscGlobalArgc;
311:   *args = PetscGlobalArgs;
312:   return(0);
313: }

317: /*@C
318:    PetscGetArguments - Allows you to access the  command line arguments anywhere
319:      after PetscInitialize() is called but before PetscFinalize().

321:    Not Collective

323:    Output Parameters:
324: .  args - the command line arguments

326:    Level: intermediate

328:    Notes:
329:       This does NOT start with the program name and IS null terminated (final arg is void)

331:    Concepts: command line arguments
332:    
333: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscFreeArguments()

335: @*/
336: PetscErrorCode  PetscGetArguments(char ***args)
337: {
338:   PetscInt       i,argc = PetscGlobalArgc;

342:   if (!PetscInitializeCalled && PetscFinalizeCalled) {
343:     SETERRQ(PETSC_ERR_ORDER,"You must call after PetscInitialize() but before PetscFinalize()");
344:   }
345:   if (!argc) {*args = 0; return(0);}
346:   PetscMalloc(argc*sizeof(char*),args);
347:   for (i=0; i<argc-1; i++) {
348:     PetscStrallocpy(PetscGlobalArgs[i+1],&(*args)[i]);
349:   }
350:   (*args)[argc-1] = 0;
351:   return(0);
352: }

356: /*@C
357:    PetscFreeArguments - Frees the memory obtained with PetscGetArguments()

359:    Not Collective

361:    Output Parameters:
362: .  args - the command line arguments 

364:    Level: intermediate

366:    Concepts: command line arguments
367:    
368: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscGetArguments()

370: @*/
371: PetscErrorCode  PetscFreeArguments(char **args)
372: {
373:   PetscInt       i = 0;

377:   if (!args) {return(0);}
378:   while (args[i]) {
379:     PetscFree(args[i]);
380:     i++;
381:   }
382:   PetscFree(args);
383:   return(0);
384: }

388: /*@C
389:    PetscInitialize - Initializes the PETSc database and MPI. 
390:    PetscInitialize() calls MPI_Init() if that has yet to be called,
391:    so this routine should always be called near the beginning of 
392:    your program -- usually the very first line! 

394:    Collective on MPI_COMM_WORLD or PETSC_COMM_WORLD if it has been set

396:    Input Parameters:
397: +  argc - count of number of command line arguments
398: .  args - the command line arguments
399: .  file - [optional] PETSc database file, also checks ~username/.petscrc and .petscrc use PETSC_NULL to not check for
400:           code specific file. Use -skip_petscrc in the code specific file to skip the .petscrc files
401: -  help - [optional] Help message to print, use PETSC_NULL for no message

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

409:    Options Database Keys:
410: +  -start_in_debugger [noxterm,dbx,xdb,gdb,...] - Starts program in debugger
411: .  -on_error_attach_debugger [noxterm,dbx,xdb,gdb,...] - Starts debugger when error detected
412: .  -on_error_emacs <machinename> causes emacsclient to jump to error file
413: .  -on_error_abort calls abort() when error detected (no traceback)
414: .  -on_error_mpiabort calls MPI_abort() when error detected
415: .  -error_output_stderr prints error messages to stderr instead of the default stdout
416: .  -error_output_none does not print the error messages (but handles errors in the same way as if this was not called)
417: .  -debugger_nodes [node1,node2,...] - Indicates nodes to start in debugger
418: .  -debugger_pause [sleeptime] (in seconds) - Pauses debugger
419: .  -stop_for_debugger - Print message on how to attach debugger manually to 
420:                         process and wait (-debugger_pause) seconds for attachment
421: .  -malloc - Indicates use of PETSc error-checking malloc (on by default for debug version of libraries)
422: .  -malloc no - Indicates not to use error-checking malloc
423: .  -malloc_debug - check for memory corruption at EVERY malloc or free
424: .  -fp_trap - Stops on floating point exceptions (Note that on the
425:               IBM RS6000 this slows code by at least a factor of 10.)
426: .  -no_signal_handler - Indicates not to trap error signals
427: .  -shared_tmp - indicates /tmp directory is shared by all processors
428: .  -not_shared_tmp - each processor has own /tmp
429: .  -tmp - alternative name of /tmp directory
430: .  -get_total_flops - returns total flops done by all processors
431: -  -memory_info - Print memory usage at end of run

433:    Options Database Keys for Profiling:
434:    See the Profiling chapter of the users manual for details.
435: +  -log_trace [filename] - Print traces of all PETSc calls
436:         to the screen (useful to determine where a program
437:         hangs without running in the debugger).  See PetscLogTraceBegin().
438: .  -info <optional filename> - Prints verbose information to the screen
439: -  -info_exclude <null,vec,mat,pc,ksp,snes,ts> - Excludes some of the verbose messages

441:    Environmental Variables:
442: +   PETSC_TMP - alternative tmp directory
443: .   PETSC_SHARED_TMP - tmp is shared by all processes
444: .   PETSC_NOT_SHARED_TMP - each process has its own private tmp
445: .   PETSC_VIEWER_SOCKET_PORT - socket number to use for socket viewer
446: -   PETSC_VIEWER_SOCKET_MACHINE - machine to use for socket viewer to connect to


449:    Level: beginner

451:    Notes:
452:    If for some reason you must call MPI_Init() separately, call
453:    it before PetscInitialize().

455:    Fortran Version:
456:    In Fortran this routine has the format
457: $       call PetscInitialize(file,ierr)

459: +   ierr - error return code
460: -  file - [optional] PETSc database file, also checks ~username/.petscrc and .petscrc use PETSC_NULL_CHARACTER to not check for
461:           code specific file. Use -skip_petscrc in the code specific file to skip the .petscrc files
462:            
463:    Important Fortran Note:
464:    In Fortran, you MUST use PETSC_NULL_CHARACTER to indicate a
465:    null character string; you CANNOT just use PETSC_NULL as 
466:    in the C version.  See the users manual for details.


469:    Concepts: initializing PETSc
470:    
471: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs()

473: @*/
474: PetscErrorCode  PetscInitialize(int *argc,char ***args,const char file[],const char help[])
475: {
477:   PetscMPIInt    flag, size;
478:   PetscInt       nodesize;
479:   PetscTruth     flg;
480:   char           hostname[256];

483:   if (PetscInitializeCalled) return(0);
484:   if (PetscFinalizeCalled)
485:     SETERRQ(1,"You can not initialize PETSc a second time")

487:   /* these must be initialized in a routine, not as a constant declaration*/
488:   PETSC_STDOUT = stdout;
489:   PETSC_STDERR = stderr;

491:   PetscOptionsCreate();

493:   /*
494:      We initialize the program name here (before MPI_Init()) because MPICH has a bug in 
495:      it that it sets args[0] on all processors to be args[0] on the first processor.
496:   */
497:   if (argc && *argc) {
498:     PetscSetProgramName(**args);
499:   } else {
500:     PetscSetProgramName("Unknown Name");
501:   }

503:   MPI_Initialized(&flag);
504:   if (!flag) {
505:     if (PETSC_COMM_WORLD) SETERRQ(PETSC_ERR_SUP,"You cannot set PETSC_COMM_WORLD if you have not initialized MPI first");
506:     MPI_Init(argc,args);
507:     PetscBeganMPI = PETSC_TRUE;
508:   }
509:   if (argc && args) {
510:     PetscGlobalArgc = *argc;
511:     PetscGlobalArgs = *args;
512:   }
513:   PetscInitializeCalled = PETSC_TRUE;
514:   PetscFinalizeCalled   = PETSC_FALSE;

516:   if (!PETSC_COMM_WORLD) {
517:     PETSC_COMM_WORLD = MPI_COMM_WORLD;
518:   }

520:   /* Done after init due to a bug in MPICH-GM? */
521:   PetscErrorPrintfInitialize();

523:   MPI_Comm_rank(MPI_COMM_WORLD,&PetscGlobalRank);
524:   MPI_Comm_size(MPI_COMM_WORLD,&PetscGlobalSize);

526: #if defined(PETSC_USE_COMPLEX)
527:   /* 
528:      Initialized the global complex variable; this is because with 
529:      shared libraries the constructors for global variables
530:      are not called; at least on IRIX.
531:   */
532:   {
533: #if defined(PETSC_CLANGUAGE_CXX)
534:     PetscScalar ic(0.0,1.0);
535:     PETSC_i = ic;
536: #else
537:     PETSC_i = I;
538: #endif
539:   }

541:   MPI_Type_contiguous(2,MPIU_REAL,&MPIU_COMPLEX);
542:   MPI_Type_commit(&MPIU_COMPLEX);
543:   MPI_Op_create(PetscSum_Local,1,&PetscSum_Op);
544: #endif

546:   /*
547:      Create the PETSc MPI reduction operator that sums of the first
548:      half of the entries and maxes the second half.
549:   */
550:   MPI_Op_create(PetscMaxSum_Local,1,&PetscMaxSum_Op);

552:   MPI_Type_contiguous(2,MPIU_SCALAR,&MPIU_2SCALAR);
553:   MPI_Type_commit(&MPIU_2SCALAR);
554:   MPI_Op_create(PetscADMax_Local,1,&PetscADMax_Op);
555:   MPI_Op_create(PetscADMin_Local,1,&PetscADMin_Op);

557:   MPI_Type_contiguous(2,MPIU_INT,&MPIU_2INT);
558:   MPI_Type_commit(&MPIU_2INT);

560:   /*
561:      Build the options database
562:   */
563:   PetscOptionsInsert(argc,args,file);

565: 
566:   /*
567:      Print main application help message
568:   */
569:   PetscOptionsHasName(PETSC_NULL,"-help",&flg);
570:   if (help && flg) {
571:     PetscPrintf(PETSC_COMM_WORLD,help);
572:   }
573:   PetscOptionsCheckInitial_Private();
574: 
575:   /* SHOULD PUT IN GUARDS: Make sure logging is initialized, even if we do not print it out */
576: #if defined(PETSC_USE_LOG)
577:   PetscLogBegin_Private();
578: #endif

580:   /*
581:      Load the dynamic libraries (on machines that support them), this registers all
582:      the solvers etc. (On non-dynamic machines this initializes the PetscDraw and PetscViewer classes)
583:   */
584:   PetscInitialize_DynamicLibraries();

586:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
587:   PetscInfo1(0,"PETSc successfully started: number of processors = %d\n",size);
588:   PetscGetHostName(hostname,256);
589:   PetscInfo1(0,"Running on machine: %s\n",hostname);

591:   PetscOptionsCheckInitial_Components();
592:   /* Check the options database for options related to the options database itself */
593:   PetscOptionsSetFromOptions();

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

599:       Currently not used because it is not supported by MPICH.
600:   */
601: #if !defined(PETSC_WORDS_BIGENDIAN)
602:   MPI_Register_datarep((char *)"petsc",PetscDataRep_read_conv_fn,PetscDataRep_write_conv_fn,PetscDataRep_extent_fn,PETSC_NULL);
603: #endif  
604: #endif

606:   PetscOptionsGetInt(PETSC_NULL,"-openmp_spawn_size",&nodesize,&flg);
607:   if (flg) {
608: #if defined(PETSC_HAVE_MPI_COMM_SPAWN)
609:     PetscOpenMPSpawn((PetscMPIInt) nodesize);
610: #else
611:     SETERRQ(PETSC_ERR_SUP,"PETSc built without MPI 2 (MPI_Comm_spawn) support, use -openmp_merge_size instead");
612: #endif
613:   } else {
614:     PetscOptionsGetInt(PETSC_NULL,"-openmp_merge_size",&nodesize,&flg);
615:     if (flg) {
616:       PetscOpenMPMerge((PetscMPIInt) nodesize);
617:     }
618:   }
619:   PetscOptionsHasName(PETSC_NULL,"-python",&flg);
620:   if (flg) {PetscPythonInitialize(PETSC_NULL,PETSC_NULL);}

622:   PetscFunctionReturn(ierr);
623: }


628: /*@C 
629:    PetscFinalize - Checks for options to be called at the conclusion
630:    of the program. MPI_Finalize() is called only if the user had not
631:    called MPI_Init() before calling PetscInitialize().

633:    Collective on PETSC_COMM_WORLD

635:    Options Database Keys:
636: +  -options_table - Calls PetscOptionsPrint()
637: .  -options_left - Prints unused options that remain in the database
638: .  -options_left no - Does not print unused options that remain in the database
639: .  -mpidump - Calls PetscMPIDump()
640: .  -malloc_dump - Calls PetscMallocDump()
641: .  -malloc_info - Prints total memory usage
642: -  -malloc_log - Prints summary of memory usage

644:    Options Database Keys for Profiling:
645:    See the Profiling chapter of the users manual for details.
646: +  -log_summary [filename] - Prints summary of flop and timing
647:         information to screen. If the filename is specified the
648:         summary is written to the file. (for code compiled with 
649:         PETSC_USE_LOG).  See PetscLogPrintSummary().
650: .  -log_all [filename] - Logs extensive profiling information
651:         (for code compiled with PETSC_USE_LOG). See PetscLogDump(). 
652: .  -log [filename] - Logs basic profiline information (for
653:         code compiled with PETSC_USE_LOG).  See PetscLogDump().
654: .  -log_sync - Log the synchronization in scatters, inner products
655:         and norms
656: -  -log_mpe [filename] - Creates a logfile viewable by the 
657:       utility Upshot/Nupshot (in MPICH distribution)

659:    Level: beginner

661:    Note:
662:    See PetscInitialize() for more general runtime options.

664: .seealso: PetscInitialize(), PetscOptionsPrint(), PetscMallocDump(), PetscMPIDump(), PetscEnd()
665: @*/
666: PetscErrorCode  PetscFinalize(void)
667: {
669:   PetscMPIInt    rank;
670:   int            nopt;
671:   PetscTruth     flg1,flg2,flg3;
673: 

676:   if (!PetscInitializeCalled) {
677:     (*PetscErrorPrintf)("PetscInitialize() must be called before PetscFinalize()\n");
678:     return(0);
679:   }
680:   PetscOpenMPFinalize();

682:   PetscOptionsHasName(PETSC_NULL,"-python",&flg1);
683:   if (flg1) {PetscPythonFinalize();}

685:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
686:   PetscOptionsHasName(PETSC_NULL,"-malloc_info",&flg2);
687:   if (!flg2) {
688:     PetscOptionsHasName(PETSC_NULL,"-memory_info",&flg2);
689:   }
690:   if (flg2) {
691:     PetscMemoryShowUsage(PETSC_VIEWER_STDOUT_WORLD,"Summary of Memory Usage in PETSc\n");
692:   }

694:   /*
695:      Free all objects registered with PetscObjectRegisterDestroy() such as PETSC_VIEWER_XXX_().
696:   */
697:   PetscObjectRegisterDestroyAll();

699:   /*
700:        Free all the registered create functions, such as KSPList, VecList, SNESList, etc
701:   */
702:   PetscFListDestroyAll();

704:   /* 
705:      Destroy any packages that registered a finalize 
706:   */
707:   PetscRegisterFinalizeAll();

709:   /*
710:      Destroy all the function registration lists created
711:   */
712:   PetscFinalize_DynamicLibraries();

714: #if defined(PETSC_USE_LOG)
715:   PetscOptionsHasName(PETSC_NULL,"-get_total_flops",&flg1);
716:   if (flg1) {
717:     PetscLogDouble flops = 0;
718:     MPI_Reduce(&_TotalFlops,&flops,1,MPI_DOUBLE,MPI_SUM,0,PETSC_COMM_WORLD);
719:     PetscPrintf(PETSC_COMM_WORLD,"Total flops over all processors %g\n",flops);
720:   }
721: #endif

723:   PetscOptionsHelpDestroyList();

725: #if defined(PETSC_USE_DEBUG)
726:   if (PetscStackActive) {
727:     PetscStackDestroy();
728:   }
729: #endif

731: #if defined(PETSC_USE_LOG)
732:   {
733:     char mname[PETSC_MAX_PATH_LEN];
734: #if defined(PETSC_HAVE_MPE)
735:     mname[0] = 0;
736:     PetscOptionsGetString(PETSC_NULL,"-log_mpe",mname,PETSC_MAX_PATH_LEN,&flg1);
737:     if (flg1){
738:       if (mname[0]) {PetscLogMPEDump(mname);}
739:       else          {PetscLogMPEDump(0);}
740:     }
741: #endif
742:     mname[0] = 0;
743:     PetscOptionsGetString(PETSC_NULL,"-log_summary",mname,PETSC_MAX_PATH_LEN,&flg1);
744:     if (flg1) {
745:       if (mname[0])  {PetscLogPrintSummary(PETSC_COMM_WORLD,mname);}
746:       else           {PetscLogPrintSummary(PETSC_COMM_WORLD,0);}
747:     }

749:     PetscOptionsGetString(PETSC_NULL,"-log_detailed",mname,PETSC_MAX_PATH_LEN,&flg1);
750:     if (flg1) {
751:       if (mname[0])  {PetscLogPrintDetailed(PETSC_COMM_WORLD,mname);}
752:       else           {PetscLogPrintDetailed(PETSC_COMM_WORLD,0);}
753:     }

755:     mname[0] = 0;
756:     PetscOptionsGetString(PETSC_NULL,"-log_all",mname,PETSC_MAX_PATH_LEN,&flg1);
757:     PetscOptionsGetString(PETSC_NULL,"-log",mname,PETSC_MAX_PATH_LEN,&flg2);
758:     if (flg1 || flg2){
759:       if (mname[0]) PetscLogDump(mname);
760:       else          PetscLogDump(0);
761:     }
762:     PetscLogDestroy();
763:   }
764: #endif
765:   PetscOptionsHasName(PETSC_NULL,"-no_signal_handler",&flg1);
766:   if (!flg1) { PetscPopSignalHandler();}
767:   PetscOptionsHasName(PETSC_NULL,"-mpidump",&flg1);
768:   if (flg1) {
769:     PetscMPIDump(stdout);
770:   }
771:   PetscOptionsHasName(PETSC_NULL,"-malloc_dump",&flg1);
772:   PetscOptionsHasName(PETSC_NULL,"-options_table",&flg2);
773:   if (flg2) {
774:     if (!rank) {PetscOptionsPrint(stdout);}
775:   }

777:   /* to prevent PETSc -options_left from warning */
778:   PetscOptionsHasName(PETSC_NULL,"-nox",&flg1);CHKERRQ(ierr)
779:   PetscOptionsHasName(PETSC_NULL,"-nox_warning",&flg1);CHKERRQ(ierr)

781:   if (!PetscOpenMPWorker) { /* worker processes skip this because they do not usually process options */
782:     flg3 = PETSC_FALSE; /* default value is required */
783:     PetscOptionsGetTruth(PETSC_NULL,"-options_left",&flg3,&flg1);
784:     PetscOptionsAllUsed(&nopt);
785:     if (flg3) {
786:       if (!flg2) { /* have not yet printed the options */
787:         PetscOptionsPrint(stdout);
788:       }
789:       if (!nopt) {
790:         PetscPrintf(PETSC_COMM_WORLD,"There are no unused options.\n");
791:       } else if (nopt == 1) {
792:         PetscPrintf(PETSC_COMM_WORLD,"There is one unused database option. It is:\n");
793:       } else {
794:         PetscPrintf(PETSC_COMM_WORLD,"There are %d unused database options. They are:\n",nopt);
795:       }
796:     }
797: #if defined(PETSC_USE_DEBUG)
798:     if (nopt && !flg3 && !flg1) {
799:       PetscPrintf(PETSC_COMM_WORLD,"WARNING! There are options you set that were not used!\n");
800:       PetscPrintf(PETSC_COMM_WORLD,"WARNING! could be spelling mistake, etc!\n");
801:       PetscOptionsLeft();
802:     } else if (nopt && flg3) {
803: #else 
804:     if (nopt && flg3) {
805: #endif
806:       PetscOptionsLeft();
807:     }
808:   }

810:   PetscOptionsHasName(PETSC_NULL,"-log_history",&flg1);
811:   if (flg1) {
812:     PetscLogCloseHistoryFile(&petsc_history);
813:     petsc_history = 0;
814:   }

816:   PetscInfoAllow(PETSC_FALSE,PETSC_NULL);

818:   PetscOptionsHasName(PETSC_NULL,"-malloc_dump",&flg1);
819:   PetscOptionsHasName(PETSC_NULL,"-malloc_log",&flg3);
820:   if (flg1) {
821:     char fname[PETSC_MAX_PATH_LEN];
822:     FILE *fd;
823:     int  err;

825:     fname[0] = 0;
826:     PetscOptionsGetString(PETSC_NULL,"-malloc_dump",fname,250,&flg1);
827:     if (flg1 && fname[0]) {
828:       char sname[PETSC_MAX_PATH_LEN];

830:       sprintf(sname,"%s_%d",fname,rank);
831:       fd   = fopen(sname,"w"); if (!fd) SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open log file: %s",sname);
832:       PetscMallocDump(fd);
833:       err = fclose(fd);
834:       if (err) SETERRQ(PETSC_ERR_SYS,"fclose() failed on file");
835:     } else {
836:       MPI_Comm local_comm;

838:       MPI_Comm_dup(MPI_COMM_WORLD,&local_comm);
839:       PetscSequentialPhaseBegin_Private(local_comm,1);
840:         PetscMallocDump(stdout);
841:       PetscSequentialPhaseEnd_Private(local_comm,1);
842:       MPI_Comm_free(&local_comm);
843:     }
844:   }
845:   if (flg3) {
846:     char fname[PETSC_MAX_PATH_LEN];
847:     FILE *fd;
848: 
849:     fname[0] = 0;
850:     PetscOptionsGetString(PETSC_NULL,"-malloc_log",fname,250,&flg1);
851:     if (flg1 && fname[0]) {
852:       char sname[PETSC_MAX_PATH_LEN];
853:       int  err;

855:       sprintf(sname,"%s_%d",fname,rank);
856:       fd   = fopen(sname,"w"); if (!fd) SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open log file: %s",sname);
857:       PetscMallocDumpLog(fd);
858:       err = fclose(fd);
859:       if (err) SETERRQ(PETSC_ERR_SYS,"fclose() failed on file");
860:     } else {
861:       PetscMallocDumpLog(stdout);
862:     }
863:   }
864:   /* Can be destroyed only after all the options are used */
865:   PetscOptionsDestroy();

867:   PetscGlobalArgc = 0;
868:   PetscGlobalArgs = 0;

870: #if defined(PETSC_USE_COMPLEX)
871:   MPI_Op_free(&PetscSum_Op);
872:   MPI_Type_free(&MPIU_COMPLEX);
873: #endif
874:   MPI_Type_free(&MPIU_2SCALAR);
875:   MPI_Type_free(&MPIU_2INT);
876:   MPI_Op_free(&PetscMaxSum_Op);
877:   MPI_Op_free(&PetscADMax_Op);
878:   MPI_Op_free(&PetscADMin_Op);

880:   PetscInfo(0,"PETSc successfully ended!\n");
881:   if (PetscBeganMPI) {
882: #if defined(PETSC_HAVE_MPI_FINALIZED)
883:     PetscMPIInt flag;
884:     MPI_Finalized(&flag);
885:     if (flag) SETERRQ(PETSC_ERR_LIB,"MPI_Finalize() has already been called, even though MPI_Init() was called by PetscInitialize()");
886: #endif
887:     MPI_Finalize();
888:   }

890:   if(PETSC_ZOPEFD != NULL){
891:     if (PETSC_ZOPEFD != PETSC_STDOUT) fprintf(PETSC_ZOPEFD, "<<<end>>>");
892:     else fprintf(PETSC_STDOUT, "<<<end>>>");}

894: /*

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

904: */
905:   PetscMallocClear();
906:   PetscInitializeCalled = PETSC_FALSE;
907:   PetscFinalizeCalled   = PETSC_TRUE;
908:   PetscFunctionReturn(ierr);
909: }

911: /*
912:      These may be used in code that ADIC is to be used on
913: */

917: /*@
918:       PetscGlobalMax - Computes the maximum value over several processors

920:      Collective on MPI_Comm

922:    Input Parameters:
923: +   local - the local value
924: -   comm - the processors that find the maximum

926:    Output Parameter:
927: .   result - the maximum value
928:   
929:    Level: intermediate

931:    Notes:
932:      These functions are to be used inside user functions that are to be processed with 
933:    ADIC. PETSc will automatically provide differentiated versions of these functions

935: .seealso: PetscGlobalMin(), PetscGlobalSum()
936: @*/
937: PetscErrorCode  PetscGlobalMax(PetscReal* local,PetscReal* result,MPI_Comm comm)
938: {
939:   return MPI_Allreduce(local,result,1,MPIU_REAL,MPI_MAX,comm);
940: }

944: /*@
945:       PetscGlobalMin - Computes the minimum value over several processors

947:      Collective on MPI_Comm

949:    Input Parameters:
950: +   local - the local value
951: -   comm - the processors that find the minimum

953:    Output Parameter:
954: .   result - the minimum value
955:   
956:    Level: intermediate

958:    Notes:
959:      These functions are to be used inside user functions that are to be processed with 
960:    ADIC. PETSc will automatically provide differentiated versions of these functions

962: .seealso: PetscGlobalMax(), PetscGlobalSum()
963: @*/
964: PetscErrorCode  PetscGlobalMin(PetscReal* local,PetscReal* result,MPI_Comm comm)
965: {
966:   return MPI_Allreduce(local,result,1,MPIU_REAL,MPI_MIN,comm);
967: }

971: /*@
972:       PetscGlobalSum - Computes the sum over sever processors

974:      Collective on MPI_Comm

976:    Input Parameters:
977: +   local - the local value
978: -   comm - the processors that find the sum

980:    Output Parameter:
981: .   result - the sum
982:   
983:    Level: intermediate

985:    Notes:
986:      These functions are to be used inside user functions that are to be processed with 
987:    ADIC. PETSc will automatically provide differentiated versions of these functions

989: .seealso: PetscGlobalMin(), PetscGlobalMax()
990: @*/
991: PetscErrorCode  PetscGlobalSum(PetscScalar* local,PetscScalar* result,MPI_Comm comm)
992: {
993:   return MPI_Allreduce(local,result,1,MPIU_SCALAR,PetscSum_Op,comm);
994: }