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

 13: /* -----------------------------------------------------------------------------------------*/


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

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

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

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

 39: PetscCookie PETSC_LARGEST_COOKIE = PETSC_SMALLEST_COOKIE;
 40: PetscCookie PETSC_OBJECT_COOKIE = 0;

 42: PetscTruth PetscPreLoadingUsed = PETSC_FALSE;
 43: PetscTruth PetscPreLoadingOn   = PETSC_FALSE;

 45: PetscErrorCode  PetscCookieRegister(PetscCookie *cookie)
 46: {
 47:   if (*cookie == PETSC_DECIDE || *cookie == PETSC_NULL) {
 48:     *cookie = ++PETSC_LARGEST_COOKIE;
 49:   } else if (*cookie > 0) {
 50:     /* Need to check here for montonicity and insert if necessary */
 51:     return 0;
 52:   } else {
 53:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Invalid suggested cookie %d", (int)*cookie);
 54:   }
 55:   return 0;
 56: }

 58: /*
 59:        Checks the options database for initializations related to the 
 60:     PETSc components
 61: */
 64: PetscErrorCode  PetscOptionsCheckInitial_Components(void)
 65: {
 66:   PetscTruth flg1;

 70:   PetscOptionsHasName(PETSC_NULL,"-help",&flg1);
 71:   if (flg1) {
 72: #if defined (PETSC_USE_LOG)
 73:     MPI_Comm   comm = PETSC_COMM_WORLD;
 74:     (*PetscHelpPrintf)(comm,"------Additional PETSc component options--------\n");
 75:     (*PetscHelpPrintf)(comm," -log_summary_exclude: <vec,mat,pc.ksp,snes>\n");
 76:     (*PetscHelpPrintf)(comm," -info_exclude: <null,vec,mat,pc,ksp,snes,ts>\n");
 77:     (*PetscHelpPrintf)(comm,"-----------------------------------------------\n");
 78: #endif
 79:   }
 80:   return(0);
 81: }

 85: /*@C
 86:       PetscInitializeNoArguments - Calls PetscInitialize() from C/C++ without
 87:         the command line arguments.

 89:    Collective
 90:   
 91:    Level: advanced

 93: .seealso: PetscInitialize(), PetscInitializeFortran()
 94: @*/
 95: PetscErrorCode  PetscInitializeNoArguments(void)
 96: {
 98:   int            argc = 0;
 99:   char           **args = 0;

102:   PetscInitialize(&argc,&args,PETSC_NULL,PETSC_NULL);
103:   PetscFunctionReturn(ierr);
104: }

108: /*@
109:       PetscInitialized - Determine whether PETSc is initialized.
110:   
111:    Level: beginner

113: .seealso: PetscInitialize(), PetscInitializeNoArguments(), PetscInitializeFortran()
114: @*/
115: PetscErrorCode  PetscInitialized(PetscTruth *isInitialized)
116: {
119:   *isInitialized = PetscInitializeCalled;
120:   return(0);
121: }

125: /*@
126:       PetscFinalized - Determine whether PetscFinalize() has been called yet
127:   
128:    Level: developer

130: .seealso: PetscInitialize(), PetscInitializeNoArguments(), PetscInitializeFortran()
131: @*/
132: PetscErrorCode  PetscFinalized(PetscTruth *isFinalized)
133: {
136:   *isFinalized = PetscFinalizeCalled;
137:   return(0);
138: }

140: EXTERN PetscErrorCode        PetscOptionsCheckInitial_Private(void);

143: /*
144:        This function is the MPI reduction operation used to compute the sum of the 
145:    first half of the datatype and the max of the second half.
146: */
147: MPI_Op PetscMaxSum_Op = 0;

152: void  PetscMaxSum_Local(void *in,void *out,int *cnt,MPI_Datatype *datatype)
153: {
154:   PetscInt *xin = (PetscInt*)in,*xout = (PetscInt*)out,i,count = *cnt;

157:   if (*datatype != MPIU_2INT) {
158:     (*PetscErrorPrintf)("Can only handle MPIU_2INT data types");
159:     MPI_Abort(MPI_COMM_WORLD,1);
160:   }

162:   for (i=0; i<count; i++) {
163:     xout[2*i]    = PetscMax(xout[2*i],xin[2*i]);
164:     xout[2*i+1] += xin[2*i+1];
165:   }
166:   PetscStackPop;
167:   return;
168: }

171: /*
172:     Returns the max of the first entry owned by this processor and the
173: sum of the second entry.
174: */
177: PetscErrorCode  PetscMaxSum(MPI_Comm comm,const PetscInt nprocs[],PetscInt *max,PetscInt *sum)
178: {
179:   PetscMPIInt    size,rank;
180:   PetscInt       *work;

184:   MPI_Comm_size(comm,&size);
185:   MPI_Comm_rank(comm,&rank);
186:   PetscMalloc(2*size*sizeof(PetscInt),&work);
187:   MPI_Allreduce((void*)nprocs,work,size,MPIU_2INT,PetscMaxSum_Op,comm);
188:   *max   = work[2*rank];
189:   *sum   = work[2*rank+1];
190:   PetscFree(work);
191:   return(0);
192: }

194: /* ----------------------------------------------------------------------------*/
195: MPI_Op  PetscADMax_Op = 0;

200: void  PetscADMax_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
201: {
202:   PetscScalar *xin = (PetscScalar *)in,*xout = (PetscScalar*)out;
203:   PetscInt    i,count = *cnt;

206:   if (*datatype != MPIU_2SCALAR) {
207:     (*PetscErrorPrintf)("Can only handle MPIU_2SCALAR data (i.e. double or complex) types");
208:     MPI_Abort(MPI_COMM_WORLD,1);
209:   }

211:   for (i=0; i<count; i++) {
212:     if (PetscRealPart(xout[2*i]) < PetscRealPart(xin[2*i])) {
213:       xout[2*i]   = xin[2*i];
214:       xout[2*i+1] = xin[2*i+1];
215:     }
216:   }

218:   PetscStackPop;
219:   return;
220: }

223: MPI_Op  PetscADMin_Op = 0;

228: void  PetscADMin_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
229: {
230:   PetscScalar *xin = (PetscScalar *)in,*xout = (PetscScalar*)out;
231:   PetscInt    i,count = *cnt;

234:   if (*datatype != MPIU_2SCALAR) {
235:     (*PetscErrorPrintf)("Can only handle MPIU_2SCALAR data (i.e. double or complex) types");
236:     MPI_Abort(MPI_COMM_WORLD,1);
237:   }

239:   for (i=0; i<count; i++) {
240:     if (PetscRealPart(xout[2*i]) > PetscRealPart(xin[2*i])) {
241:       xout[2*i]   = xin[2*i];
242:       xout[2*i+1] = xin[2*i+1];
243:     }
244:   }

246:   PetscStackPop;
247:   return;
248: }
250: /* ---------------------------------------------------------------------------------------*/

252: #if defined(PETSC_USE_COMPLEX)
253: MPI_Op PetscSum_Op = 0;

258: void  PetscSum_Local(void *in,void *out,PetscMPIInt *cnt,MPI_Datatype *datatype)
259: {
260:   PetscScalar *xin = (PetscScalar *)in,*xout = (PetscScalar*)out;
261:   PetscInt         i,count = *cnt;

264:   if (*datatype != MPIU_SCALAR) {
265:     (*PetscErrorPrintf)("Can only handle MPIU_SCALAR data (i.e. double or complex) types");
266:     MPI_Abort(MPI_COMM_WORLD,1);
267:   }

269:   for (i=0; i<count; i++) {
270:     xout[i] += xin[i];
271:   }

273:   PetscStackPop;
274:   return;
275: }
277: #endif

279: static int  PetscGlobalArgc   = 0;
280: static char **PetscGlobalArgs = 0;

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

288:    Not Collective

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

294:    Level: intermediate

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

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

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

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

322:    Not Collective

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

327:    Level: intermediate

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

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

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

343:   if (!PetscGlobalArgs) {
344:     SETERRQ(PETSC_ERR_ORDER,"You must call after PetscInitialize() but before PetscFinalize()");
345:   }
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:   while (args[i]) {
378:     PetscFree(args[i]);
379:     i++;
380:   }
381:   PetscFree(args);
382:   return(0);
383: }

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

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

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

402:    If you wish PETSc to run on a subcommunicator of MPI_COMM_WORLD, create that
403:    communicator first and assign it to PETSC_COMM_WORLD BEFORE calling PetscInitialize()

405:    Options Database Keys:
406: +  -start_in_debugger [noxterm,dbx,xdb,gdb,...] - Starts program in debugger
407: .  -on_error_attach_debugger [noxterm,dbx,xdb,gdb,...] - Starts debugger when error detected
408: .  -on_error_emacs <machinename> causes emacsclient to jump to error file
409: .  -debugger_nodes [node1,node2,...] - Indicates nodes to start in debugger
410: .  -debugger_pause [sleeptime] (in seconds) - Pauses debugger
411: .  -stop_for_debugger - Print message on how to attach debugger manually to 
412:                         process and wait (-debugger_pause) seconds for attachment
413: .  -malloc - Indicates use of PETSc error-checking malloc (on by default for debug version of libraries)
414: .  -malloc no - Indicates not to use error-checking malloc
415: .  -malloc_debug - check for memory corruption at EVERY malloc or free
416: .  -fp_trap - Stops on floating point exceptions (Note that on the
417:               IBM RS6000 this slows code by at least a factor of 10.)
418: .  -no_signal_handler - Indicates not to trap error signals
419: .  -shared_tmp - indicates /tmp directory is shared by all processors
420: .  -not_shared_tmp - each processor has own /tmp
421: .  -tmp - alternative name of /tmp directory
422: .  -get_total_flops - returns total flops done by all processors
423: -  -memory_info - Print memory usage at end of run

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

433:    Environmental Variables:
434: +   PETSC_TMP - alternative tmp directory
435: .   PETSC_SHARED_TMP - tmp is shared by all processes
436: .   PETSC_NOT_SHARED_TMP - each process has its own private tmp
437: .   PETSC_VIEWER_SOCKET_PORT - socket number to use for socket viewer
438: -   PETSC_VIEWER_SOCKET_MACHINE - machine to use for socket viewer to connect to


441:    Level: beginner

443:    Notes:
444:    If for some reason you must call MPI_Init() separately, call
445:    it before PetscInitialize().

447:    Fortran Version:
448:    In Fortran this routine has the format
449: $       call PetscInitialize(file,ierr)

451: +   ierr - error return code
452: -  file - [optional] PETSc database file, also checks ~username/.petscrc and .petscrc use PETSC_CHARACTER_NULL to not check for
453:           code specific file. Use -skip_petscrc in the code specific file to skip the .petscrc files
454:            
455:    Important Fortran Note:
456:    In Fortran, you MUST use PETSC_NULL_CHARACTER to indicate a
457:    null character string; you CANNOT just use PETSC_NULL as 
458:    in the C version.  See the users manual for details.


461:    Concepts: initializing PETSc
462:    
463: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs()

465: @*/
466: PetscErrorCode  PetscInitialize(int *argc,char ***args,const char file[],const char help[])
467: {
469:   PetscMPIInt    flag, size,nodesize;
470:   PetscTruth     flg;
471:   char           hostname[256];

474:   if (PetscInitializeCalled) return(0);

476:   /* this must be initialized in a routine, not as a constant declaration*/
477:   PETSC_STDOUT = stdout;

479:   PetscOptionsCreate();

481:   /*
482:      We initialize the program name here (before MPI_Init()) because MPICH has a bug in 
483:      it that it sets args[0] on all processors to be args[0] on the first processor.
484:   */
485:   if (argc && *argc) {
486:     PetscSetProgramName(**args);
487:   } else {
488:     PetscSetProgramName("Unknown Name");
489:   }


492:   MPI_Initialized(&flag);
493:   if (!flag) {
494:     if (PETSC_COMM_WORLD) SETERRQ(PETSC_ERR_SUP,"You cannot set PETSC_COMM_WORLD if you have not initialized MPI first");
495:     MPI_Init(argc,args);
496:     PetscBeganMPI = PETSC_TRUE;
497:   }
498:   if (argc && args) {
499:     PetscGlobalArgc = *argc;
500:     PetscGlobalArgs = *args;
501:   }
502:   PetscInitializeCalled = PETSC_TRUE;
503:   PetscFinalizeCalled   = PETSC_FALSE;

505:   if (!PETSC_COMM_WORLD) {
506:     PETSC_COMM_WORLD = MPI_COMM_WORLD;
507:   }

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

512:   MPI_Comm_rank(MPI_COMM_WORLD,&PetscGlobalRank);
513:   MPI_Comm_size(MPI_COMM_WORLD,&PetscGlobalSize);

515: #if defined(PETSC_USE_COMPLEX)
516:   /* 
517:      Initialized the global complex variable; this is because with 
518:      shared libraries the constructors for global variables
519:      are not called; at least on IRIX.
520:   */
521:   {
522: #if defined(PETSC_CLANGUAGE_CXX)
523:     PetscScalar ic(0.0,1.0);
524:     PETSC_i = ic;
525: #else
526:     PetscScalar ic;
527:     ic = 1.I;
528:     PETSC_i = ic;
529: #endif
530:   }

532:   MPI_Type_contiguous(2,MPIU_REAL,&MPIU_COMPLEX);
533:   MPI_Type_commit(&MPIU_COMPLEX);
534:   MPI_Op_create(PetscSum_Local,1,&PetscSum_Op);
535: #endif

537:   /*
538:      Create the PETSc MPI reduction operator that sums of the first
539:      half of the entries and maxes the second half.
540:   */
541:   MPI_Op_create(PetscMaxSum_Local,1,&PetscMaxSum_Op);

543:   MPI_Type_contiguous(2,MPIU_SCALAR,&MPIU_2SCALAR);
544:   MPI_Type_commit(&MPIU_2SCALAR);
545:   MPI_Op_create(PetscADMax_Local,1,&PetscADMax_Op);
546:   MPI_Op_create(PetscADMin_Local,1,&PetscADMin_Op);

548:   MPI_Type_contiguous(2,MPIU_INT,&MPIU_2INT);
549:   MPI_Type_commit(&MPIU_2INT);

551:   /*
552:      Build the options database and check for user setup requests
553:   */
554:   PetscOptionsInsert(argc,args,file);

556:   /*
557:      Print main application help message
558:   */
559:   PetscOptionsHasName(PETSC_NULL,"-help",&flg);
560:   if (help && flg) {
561:     PetscPrintf(PETSC_COMM_WORLD,help);
562:   }
563:   PetscOptionsCheckInitial_Private();

565:   /* SHOULD PUT IN GUARDS: Make sure logging is initialized, even if we do not print it out */
566: #if defined(PETSC_USE_LOG)
567:   PetscLogBegin_Private();
568: #endif

570:   /*
571:      Load the dynamic libraries (on machines that support them), this registers all
572:      the solvers etc. (On non-dynamic machines this initializes the PetscDraw and PetscViewer classes)
573:   */
574:   PetscInitialize_DynamicLibraries();

576:   MPI_Comm_size(PETSC_COMM_WORLD,&size);
577:   PetscInfo1(0,"PETSc successfully started: number of processors = %d\n",size);
578:   PetscGetHostName(hostname,256);
579:   PetscInfo1(0,"Running on machine: %s\n",hostname);

581:   PetscOptionsCheckInitial_Components();
582:   /* Check the options database for options related to the options database itself */
583:   PetscOptionsSetFromOptions();

585:   PetscOptionsGetInt(PETSC_NULL,"-openmp_spawn_size",&nodesize,&flg);
586:   if (flg) {
587: #if defined(PETSC_HAVE_MPI_COMM_SPAWN)
588:     PetscOpenMPSpawn(nodesize);
589: #else
590:     SETERRQ(PETSC_ERR_SUP,"PETSc built without MPI 2 (MPI_Comm_spawn) support, use -openmp_node_size instead");
591: #endif
592:   } else {
593:     PetscOptionsGetInt(PETSC_NULL,"-openmp_merge_size",&nodesize,&flg);
594:     if (flg) {
595:       PetscOpenMPMerge(nodesize);
596:     }
597:   }

599:   PetscFunctionReturn(ierr);
600: }


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

610:    Collective on PETSC_COMM_WORLD

612:    Options Database Keys:
613: +  -options_table - Calls PetscOptionsPrint()
614: .  -options_left - Prints unused options that remain in the database
615: .  -options_left no - Does not print unused options that remain in the database
616: .  -mpidump - Calls PetscMPIDump()
617: .  -malloc_dump - Calls PetscMallocDump()
618: .  -malloc_info - Prints total memory usage
619: -  -malloc_log - Prints summary of memory usage

621:    Options Database Keys for Profiling:
622:    See the Profiling chapter of the users manual for details.
623: +  -log_summary [filename] - Prints summary of flop and timing
624:         information to screen. If the filename is specified the
625:         summary is written to the file. (for code compiled with 
626:         PETSC_USE_LOG).  See PetscLogPrintSummary().
627: .  -log_all [filename] - Logs extensive profiling information
628:         (for code compiled with PETSC_USE_LOG). See PetscLogDump(). 
629: .  -log [filename] - Logs basic profiline information (for
630:         code compiled with PETSC_USE_LOG).  See PetscLogDump().
631: .  -log_sync - Log the synchronization in scatters, inner products
632:         and norms
633: -  -log_mpe [filename] - Creates a logfile viewable by the 
634:       utility Upshot/Nupshot (in MPICH distribution)

636:    Level: beginner

638:    Note:
639:    See PetscInitialize() for more general runtime options.

641: .seealso: PetscInitialize(), PetscOptionsPrint(), PetscMallocDump(), PetscMPIDump(), PetscEnd()
642: @*/
643: PetscErrorCode  PetscFinalize(void)
644: {
646:   PetscMPIInt    rank;
647:   int            nopt;
648:   PetscTruth     flg1,flg2,flg3;
649: 

652:   if (!PetscInitializeCalled) {
653:     (*PetscErrorPrintf)("PetscInitialize() must be called before PetscFinalize()\n");
654:     return(0);
655:   }

657:   PetscOpenMPFinalize();

659:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
660:   PetscOptionsHasName(PETSC_NULL,"-malloc_info",&flg2);
661:   if (!flg2) {
662:     PetscOptionsHasName(PETSC_NULL,"-memory_info",&flg2);
663:   }
664:   if (flg2) {
665:     PetscMemoryShowUsage(PETSC_VIEWER_STDOUT_WORLD,"Summary of Memory Usage in PETSc\n");
666:   }

668:   /* Destroy auxiliary packages */
669: #if defined(PETSC_HAVE_MATHEMATICA)
670:   PetscViewerMathematicaFinalizePackage();
671: #endif
672: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
673: #if defined(PETSC_HAVE_SIEVE)
674:   {
675:     PetscErrorCode (*func)(void);
676:     char lib[PETSC_MAX_PATH_LEN];

678:     PetscStrcpy(lib,"${PETSC_LIB_DIR}");
679:     PetscStrcat(lib,"/libpetscdm");
680:     PetscDLLibrarySym(PETSC_COMM_WORLD,&DLLibrariesLoaded,lib,"DMFinalizePackage",(void **) &func);
681:     if (func) {
682:       (*func)();
683:     }
684:   }
685: #endif
686: #endif

688:   /*
689:      Destroy all the function registration lists created
690:   */
691:   PetscFinalize_DynamicLibraries();

693: #if defined(PETSC_USE_LOG)
694:   PetscOptionsHasName(PETSC_NULL,"-get_total_flops",&flg1);
695:   if (flg1) {
696:     PetscLogDouble flops = 0;
697:     MPI_Reduce(&_TotalFlops,&flops,1,MPI_DOUBLE,MPI_SUM,0,PETSC_COMM_WORLD);
698:     PetscPrintf(PETSC_COMM_WORLD,"Total flops over all processors %g\n",flops);
699:   }
700: #endif

702:   /*
703:      Free all objects registered with PetscObjectRegisterDestroy() such ast
704:     PETSC_VIEWER_XXX_().
705:   */
706:   PetscObjectRegisterDestroyAll();
707:   PetscOptionsHelpDestroyList();

709: #if defined(PETSC_USE_DEBUG)
710:   if (PetscStackActive) {
711:     PetscStackDestroy();
712:   }
713: #endif

715: #if defined(PETSC_USE_LOG)
716:   {
717:     char mname[PETSC_MAX_PATH_LEN];
718: #if defined(PETSC_HAVE_MPE)
719:     mname[0] = 0;
720:     PetscOptionsGetString(PETSC_NULL,"-log_mpe",mname,PETSC_MAX_PATH_LEN,&flg1);
721:     if (flg1){
722:       if (mname[0]) {PetscLogMPEDump(mname);}
723:       else          {PetscLogMPEDump(0);}
724:     }
725: #endif
726:     mname[0] = 0;
727:     PetscOptionsGetString(PETSC_NULL,"-log_summary",mname,PETSC_MAX_PATH_LEN,&flg1);
728:     if (flg1) {
729:       if (mname[0])  {PetscLogPrintSummary(PETSC_COMM_WORLD,mname);}
730:       else           {PetscLogPrintSummary(PETSC_COMM_WORLD,0);}
731:     }

733:     mname[0] = 0;
734:     PetscOptionsGetString(PETSC_NULL,"-log_all",mname,PETSC_MAX_PATH_LEN,&flg1);
735:     PetscOptionsGetString(PETSC_NULL,"-log",mname,PETSC_MAX_PATH_LEN,&flg2);
736:     if (flg1 || flg2){
737:       if (mname[0]) PetscLogDump(mname);
738:       else          PetscLogDump(0);
739:     }
740:     PetscLogDestroy();
741:   }
742: #endif
743:   PetscOptionsHasName(PETSC_NULL,"-no_signal_handler",&flg1);
744:   if (!flg1) { PetscPopSignalHandler();}
745:   PetscOptionsHasName(PETSC_NULL,"-mpidump",&flg1);
746:   if (flg1) {
747:     PetscMPIDump(stdout);
748:   }
749:   PetscOptionsHasName(PETSC_NULL,"-malloc_dump",&flg1);
750:   PetscOptionsHasName(PETSC_NULL,"-options_table",&flg2);
751:   if (flg2) {
752:     if (!rank) {PetscOptionsPrint(stdout);}
753:   }

755:   /* to prevent PETSc -options_left from warning */
756:   PetscOptionsHasName(PETSC_NULL,"-nox_warning",&flg1);CHKERRQ(ierr)
757:   PetscOptionsHasName(PETSC_NULL,"-error_output_stderr",&flg1);

759:   flg3 = PETSC_FALSE; /* default value is required */
760:   PetscOptionsGetTruth(PETSC_NULL,"-options_left",&flg3,&flg1);
761:   PetscOptionsAllUsed(&nopt);
762:   if (flg3) {
763:     if (!flg2) { /* have not yet printed the options */
764:       PetscOptionsPrint(stdout);
765:     }
766:     if (!nopt) {
767:       PetscPrintf(PETSC_COMM_WORLD,"There are no unused options.\n");
768:     } else if (nopt == 1) {
769:       PetscPrintf(PETSC_COMM_WORLD,"There is one unused database option. It is:\n");
770:     } else {
771:       PetscPrintf(PETSC_COMM_WORLD,"There are %d unused database options. They are:\n",nopt);
772:     }
773:   }
774: #if defined(PETSC_USE_DEBUG)
775:   if (nopt && !flg3 && !flg1) {
776:     PetscPrintf(PETSC_COMM_WORLD,"WARNING! There are options you set that were not used!\n");
777:     PetscPrintf(PETSC_COMM_WORLD,"WARNING! could be spelling mistake, etc!\n");
778:     PetscOptionsLeft();
779:   } else if (nopt && flg3) {
780: #else 
781:   if (nopt && flg3) {
782: #endif
783:     PetscOptionsLeft();
784:   }

786:   PetscOptionsHasName(PETSC_NULL,"-log_history",&flg1);
787:   if (flg1) {
788:     PetscLogCloseHistoryFile(&petsc_history);
789:     petsc_history = 0;
790:   }

792:   PetscInfoAllow(PETSC_FALSE,PETSC_NULL);

794:   /*
795:        Free all the registered create functions, such as KSPList, VecList, SNESList, etc
796:   */
797:   PetscFListDestroyAll();

799:   PetscOptionsHasName(PETSC_NULL,"-malloc_dump",&flg1);
800:   PetscOptionsHasName(PETSC_NULL,"-malloc_log",&flg3);
801:   if (flg1) {
802:     char fname[PETSC_MAX_PATH_LEN];
803:     FILE *fd;
804: 
805:     fname[0] = 0;
806:     PetscOptionsGetString(PETSC_NULL,"-malloc_dump",fname,250,&flg1);
807:     if (flg1 && fname[0]) {
808:       char sname[PETSC_MAX_PATH_LEN];

810:       sprintf(sname,"%s_%d",fname,rank);
811:       fd   = fopen(sname,"w"); if (!fd) SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open log file: %s",sname);
812:       PetscMallocDump(fd);
813:       fclose(fd);
814:     } else {
815:       MPI_Comm local_comm;

817:       MPI_Comm_dup(MPI_COMM_WORLD,&local_comm);
818:       PetscSequentialPhaseBegin_Private(local_comm,1);
819:         PetscMallocDump(stdout);
820:       PetscSequentialPhaseEnd_Private(local_comm,1);
821:       MPI_Comm_free(&local_comm);
822:     }
823:   }
824:   if (flg3) {
825:     char fname[PETSC_MAX_PATH_LEN];
826:     FILE *fd;
827: 
828:     fname[0] = 0;
829:     PetscOptionsGetString(PETSC_NULL,"-malloc_log",fname,250,&flg1);
830:     if (flg1 && fname[0]) {
831:       char sname[PETSC_MAX_PATH_LEN];

833:       sprintf(sname,"%s_%d",fname,rank);
834:       fd   = fopen(sname,"w"); if (!fd) SETERRQ1(PETSC_ERR_FILE_OPEN,"Cannot open log file: %s",sname);
835:       PetscMallocDumpLog(fd);
836:       fclose(fd);
837:     } else {
838:       PetscMallocDumpLog(stdout);
839:     }
840:   }
841:   /* Can be destroyed only after all the options are used */
842:   PetscOptionsDestroy();

844:   PetscGlobalArgc = 0;
845:   PetscGlobalArgs = 0;

847: #if defined(PETSC_USE_COMPLEX)
848:   MPI_Op_free(&PetscSum_Op);
849:   MPI_Type_free(&MPIU_COMPLEX);
850: #endif
851:   MPI_Type_free(&MPIU_2SCALAR);
852:   MPI_Type_free(&MPIU_2INT);
853:   MPI_Op_free(&PetscMaxSum_Op);
854:   MPI_Op_free(&PetscADMax_Op);
855:   MPI_Op_free(&PetscADMin_Op);

857:   PetscInfo(0,"PETSc successfully ended!\n");
858:   if (PetscBeganMPI) {
859:     MPI_Finalize();
860:   }

862: /*

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

872: */
873:   PetscClearMalloc();
874:   PetscInitializeCalled = PETSC_FALSE;
875:   PetscFinalizeCalled   = PETSC_TRUE;
876:   PetscFunctionReturn(ierr);
877: }

879: /*
880:      These may be used in code that ADIC is to be used on
881: */

885: /*@C
886:       PetscGlobalMax - Computes the maximum value over several processors

888:      Collective on MPI_Comm

890:    Input Parameters:
891: +   local - the local value
892: -   comm - the processors that find the maximum

894:    Output Parameter:
895: .   result - the maximum value
896:   
897:    Level: intermediate

899:    Notes:
900:      These functions are to be used inside user functions that are to be processed with 
901:    ADIC. PETSc will automatically provide differentiated versions of these functions

903: .seealso: PetscGlobalMin(), PetscGlobalSum()
904: @*/
905: PetscErrorCode  PetscGlobalMax(PetscReal* local,PetscReal* result,MPI_Comm comm)
906: {
907:   return MPI_Allreduce(local,result,1,MPIU_REAL,MPI_MAX,comm);
908: }

912: /*@C
913:       PetscGlobalMin - Computes the minimum value over several processors

915:      Collective on MPI_Comm

917:    Input Parameters:
918: +   local - the local value
919: -   comm - the processors that find the minimum

921:    Output Parameter:
922: .   result - the minimum value
923:   
924:    Level: intermediate

926:    Notes:
927:      These functions are to be used inside user functions that are to be processed with 
928:    ADIC. PETSc will automatically provide differentiated versions of these functions

930: .seealso: PetscGlobalMax(), PetscGlobalSum()
931: @*/
932: PetscErrorCode  PetscGlobalMin(PetscReal* local,PetscReal* result,MPI_Comm comm)
933: {
934:   return MPI_Allreduce(local,result,1,MPIU_REAL,MPI_MIN,comm);
935: }

939: /*@C
940:       PetscGlobalSum - Computes the sum over sever processors

942:      Collective on MPI_Comm

944:    Input Parameters:
945: +   local - the local value
946: -   comm - the processors that find the sum

948:    Output Parameter:
949: .   result - the sum
950:   
951:    Level: intermediate

953:    Notes:
954:      These functions are to be used inside user functions that are to be processed with 
955:    ADIC. PETSc will automatically provide differentiated versions of these functions

957: .seealso: PetscGlobalMin(), PetscGlobalMax()
958: @*/
959: PetscErrorCode  PetscGlobalSum(PetscScalar* local,PetscScalar* result,MPI_Comm comm)
960: {
961:   return MPI_Allreduce(local,result,1,MPIU_SCALAR,PetscSum_Op,comm);
962: }