Actual source code: plog.c


  2: /*
  3:       PETSc code to log object creation and destruction and PETSc events.

  5:       This provides the public API used by the rest of PETSc and by users.

  7:       These routines use a private API that is not used elsewhere in PETSc and is not
  8:       accessible to users. The private API is defined in logimpl.h and the utils directory.

 10: */
 11: #include <petsc/private/logimpl.h>
 12: #include <petsctime.h>
 13: #include <petscviewer.h>

 15: PetscErrorCode PetscLogObjectParent(PetscObject p,PetscObject c)
 16: {
 17:   if (!c || !p) return 0;
 18:   c->parent   = p;
 19:   c->parentid = p->id;
 20:   return 0;
 21: }

 23: /*@C
 24:    PetscLogObjectMemory - Adds to an object a count of additional amount of memory that is used by the object.

 26:    Not collective.

 28:    Input Parameters:
 29: +  obj  - the PETSc object
 30: -  mem  - the amount of memory that is being added to the object

 32:    Level: developer

 34:    Developer Notes:
 35:     Currently we do not always do a good job of associating all memory allocations with an object.

 37: .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscInitializeNoArguments()

 39: @*/
 40: PetscErrorCode PetscLogObjectMemory(PetscObject p,PetscLogDouble m)
 41: {
 42:   if (!p) return 0;
 43:   p->mem += m;
 44:   return 0;
 45: }

 47: PetscLogEvent PETSC_LARGEST_EVENT = PETSC_EVENT;

 49: #if defined(PETSC_USE_LOG)
 50: #include <petscmachineinfo.h>
 51: #include <petscconfiginfo.h>

 53: /* used in the MPI_XXX() count macros in petsclog.h */

 55: /* Action and object logging variables */
 56: Action    *petsc_actions            = NULL;
 57: Object    *petsc_objects            = NULL;
 58: PetscBool petsc_logActions          = PETSC_FALSE;
 59: PetscBool petsc_logObjects          = PETSC_FALSE;
 60: int       petsc_numActions          = 0, petsc_maxActions = 100;
 61: int       petsc_numObjects          = 0, petsc_maxObjects = 100;
 62: int       petsc_numObjectsDestroyed = 0;

 64: /* Global counters */
 65: PetscLogDouble petsc_BaseTime        = 0.0;
 66: PetscLogDouble petsc_TotalFlops      = 0.0;  /* The number of flops */
 67: PetscLogDouble petsc_tmp_flops       = 0.0;  /* The incremental number of flops */
 68: PetscLogDouble petsc_send_ct         = 0.0;  /* The number of sends */
 69: PetscLogDouble petsc_recv_ct         = 0.0;  /* The number of receives */
 70: PetscLogDouble petsc_send_len        = 0.0;  /* The total length of all sent messages */
 71: PetscLogDouble petsc_recv_len        = 0.0;  /* The total length of all received messages */
 72: PetscLogDouble petsc_isend_ct        = 0.0;  /* The number of immediate sends */
 73: PetscLogDouble petsc_irecv_ct        = 0.0;  /* The number of immediate receives */
 74: PetscLogDouble petsc_isend_len       = 0.0;  /* The total length of all immediate send messages */
 75: PetscLogDouble petsc_irecv_len       = 0.0;  /* The total length of all immediate receive messages */
 76: PetscLogDouble petsc_wait_ct         = 0.0;  /* The number of waits */
 77: PetscLogDouble petsc_wait_any_ct     = 0.0;  /* The number of anywaits */
 78: PetscLogDouble petsc_wait_all_ct     = 0.0;  /* The number of waitalls */
 79: PetscLogDouble petsc_sum_of_waits_ct = 0.0;  /* The total number of waits */
 80: PetscLogDouble petsc_allreduce_ct    = 0.0;  /* The number of reductions */
 81: PetscLogDouble petsc_gather_ct       = 0.0;  /* The number of gathers and gathervs */
 82: PetscLogDouble petsc_scatter_ct      = 0.0;  /* The number of scatters and scattervs */
 83: #if defined(PETSC_HAVE_DEVICE)
 84: PetscLogDouble petsc_ctog_ct         = 0.0;  /* The total number of CPU to GPU copies */
 85: PetscLogDouble petsc_gtoc_ct         = 0.0;  /* The total number of GPU to CPU copies */
 86: PetscLogDouble petsc_ctog_sz         = 0.0;  /* The total size of CPU to GPU copies */
 87: PetscLogDouble petsc_gtoc_sz         = 0.0;  /* The total size of GPU to CPU copies */
 88: PetscLogDouble petsc_ctog_ct_scalar  = 0.0;  /* The total number of CPU to GPU copies */
 89: PetscLogDouble petsc_gtoc_ct_scalar  = 0.0;  /* The total number of GPU to CPU copies */
 90: PetscLogDouble petsc_ctog_sz_scalar  = 0.0;  /* The total size of CPU to GPU copies */
 91: PetscLogDouble petsc_gtoc_sz_scalar  = 0.0;  /* The total size of GPU to CPU copies */
 92: PetscLogDouble petsc_gflops          = 0.0;  /* The flops done on a GPU */
 93: PetscLogDouble petsc_gtime           = 0.0;  /* The time spent on a GPU */

 95: #if defined(PETSC_USE_DEBUG)
 96: PetscBool petsc_gtime_inuse = PETSC_FALSE;
 97: #endif
 98: #endif

100: /* Logging functions */
101: PetscErrorCode (*PetscLogPHC)(PetscObject) = NULL;
102: PetscErrorCode (*PetscLogPHD)(PetscObject) = NULL;
103: PetscErrorCode (*PetscLogPLB)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL;
104: PetscErrorCode (*PetscLogPLE)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL;

106: /* Tracing event logging variables */
107: FILE             *petsc_tracefile            = NULL;
108: int              petsc_tracelevel            = 0;
109: const char       *petsc_traceblanks          = "                                                                                                    ";
110: char             petsc_tracespace[128]       = " ";
111: PetscLogDouble   petsc_tracetime             = 0.0;
112: static PetscBool PetscLogInitializeCalled = PETSC_FALSE;

114: PETSC_INTERN PetscErrorCode PetscLogInitialize(void)
115: {
116:   int            stage;
117:   PetscBool      opt;

121:   if (PetscLogInitializeCalled) return(0);
122:   PetscLogInitializeCalled = PETSC_TRUE;

124:   PetscOptionsHasName(NULL,NULL, "-log_exclude_actions", &opt);
125:   if (opt) petsc_logActions = PETSC_FALSE;
126:   PetscOptionsHasName(NULL,NULL, "-log_exclude_objects", &opt);
127:   if (opt) petsc_logObjects = PETSC_FALSE;
128:   if (petsc_logActions) {
129:     PetscMalloc1(petsc_maxActions, &petsc_actions);
130:   }
131:   if (petsc_logObjects) {
132:     PetscMalloc1(petsc_maxObjects, &petsc_objects);
133:   }
134:   PetscLogPHC = PetscLogObjCreateDefault;
135:   PetscLogPHD = PetscLogObjDestroyDefault;
136:   /* Setup default logging structures */
137:   PetscStageLogCreate(&petsc_stageLog);
138:   PetscStageLogRegister(petsc_stageLog, "Main Stage", &stage);

140:   /* All processors sync here for more consistent logging */
141:   MPI_Barrier(PETSC_COMM_WORLD);
142:   PetscTime(&petsc_BaseTime);
143:   PetscLogStagePush(stage);
144:   return(0);
145: }

147: PETSC_INTERN PetscErrorCode PetscLogFinalize(void)
148: {
149:   PetscStageLog  stageLog;

153:   PetscFree(petsc_actions);
154:   PetscFree(petsc_objects);
155:   PetscLogNestedEnd();
156:   PetscLogSet(NULL, NULL);

158:   /* Resetting phase */
159:   PetscLogGetStageLog(&stageLog);
160:   PetscStageLogDestroy(stageLog);

162:   petsc_TotalFlops            = 0.0;
163:   petsc_numActions            = 0;
164:   petsc_numObjects            = 0;
165:   petsc_numObjectsDestroyed   = 0;
166:   petsc_maxActions            = 100;
167:   petsc_maxObjects            = 100;
168:   petsc_actions               = NULL;
169:   petsc_objects               = NULL;
170:   petsc_logActions            = PETSC_FALSE;
171:   petsc_logObjects            = PETSC_FALSE;
172:   petsc_BaseTime              = 0.0;
173:   petsc_TotalFlops            = 0.0;
174:   petsc_tmp_flops             = 0.0;
175:   petsc_send_ct               = 0.0;
176:   petsc_recv_ct               = 0.0;
177:   petsc_send_len              = 0.0;
178:   petsc_recv_len              = 0.0;
179:   petsc_isend_ct              = 0.0;
180:   petsc_irecv_ct              = 0.0;
181:   petsc_isend_len             = 0.0;
182:   petsc_irecv_len             = 0.0;
183:   petsc_wait_ct               = 0.0;
184:   petsc_wait_any_ct           = 0.0;
185:   petsc_wait_all_ct           = 0.0;
186:   petsc_sum_of_waits_ct       = 0.0;
187:   petsc_allreduce_ct          = 0.0;
188:   petsc_gather_ct             = 0.0;
189:   petsc_scatter_ct            = 0.0;
190:   #if defined(PETSC_HAVE_DEVICE)
191:   petsc_ctog_ct               = 0.0;
192:   petsc_gtoc_ct               = 0.0;
193:   petsc_ctog_sz               = 0.0;
194:   petsc_gtoc_sz               = 0.0;
195:   petsc_gflops                = 0.0;
196:   petsc_gtime                 = 0.0;
197:   #endif
198:   PETSC_LARGEST_EVENT         = PETSC_EVENT;
199:   PetscLogPHC                 = NULL;
200:   PetscLogPHD                 = NULL;
201:   petsc_tracefile             = NULL;
202:   petsc_tracelevel            = 0;
203:   petsc_traceblanks           = "                                                                                                    ";
204:   petsc_tracespace[0]         = ' '; petsc_tracespace[1] = 0;
205:   petsc_tracetime             = 0.0;
206:   PETSC_LARGEST_CLASSID       = PETSC_SMALLEST_CLASSID;
207:   PETSC_OBJECT_CLASSID        = 0;
208:   petsc_stageLog              = NULL;
209:   PetscLogInitializeCalled    = PETSC_FALSE;
210:   return(0);
211: }

213: /*@C
214:   PetscLogSet - Sets the logging functions called at the beginning and ending of every event.

216:   Not Collective

218:   Input Parameters:
219: + b - The function called at beginning of event
220: - e - The function called at end of event

222:   Level: developer

224: .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogTraceBegin()
225: @*/
226: PetscErrorCode  PetscLogSet(PetscErrorCode (*b)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject),
227:                             PetscErrorCode (*e)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject))
228: {
230:   PetscLogPLB = b;
231:   PetscLogPLE = e;
232:   return(0);
233: }

235: /*@C
236:   PetscLogDefaultBegin - Turns on logging of objects and events. This logs flop
237:   rates and object creation and should not slow programs down too much.
238:   This routine may be called more than once.

240:   Logically Collective over PETSC_COMM_WORLD

242:   Options Database Keys:
243: . -log_view [viewertype:filename:viewerformat] - Prints summary of flop and timing information to the
244:                   screen (for code configured with --with-log=1 (which is the default))

246:   Usage:
247: .vb
248:       PetscInitialize(...);
249:       PetscLogDefaultBegin();
250:        ... code ...
251:       PetscLogView(viewer); or PetscLogDump();
252:       PetscFinalize();
253: .ve

255:   Notes:
256:   PetscLogView(viewer) or PetscLogDump() actually cause the printing of
257:   the logging information.

259:   Level: advanced

261: .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogTraceBegin()
262: @*/
263: PetscErrorCode  PetscLogDefaultBegin(void)
264: {

268:   PetscLogSet(PetscLogEventBeginDefault, PetscLogEventEndDefault);
269:   return(0);
270: }

272: /*@C
273:   PetscLogAllBegin - Turns on extensive logging of objects and events. Logs
274:   all events. This creates large log files and slows the program down.

276:   Logically Collective on PETSC_COMM_WORLD

278:   Options Database Keys:
279: . -log_all - Prints extensive log information

281:   Usage:
282: .vb
283:      PetscInitialize(...);
284:      PetscLogAllBegin();
285:      ... code ...
286:      PetscLogDump(filename);
287:      PetscFinalize();
288: .ve

290:   Notes:
291:   A related routine is PetscLogDefaultBegin() (with the options key -log), which is
292:   intended for production runs since it logs only flop rates and object
293:   creation (and shouldn't significantly slow the programs).

295:   Level: advanced

297: .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogTraceBegin()
298: @*/
299: PetscErrorCode  PetscLogAllBegin(void)
300: {

304:   PetscLogSet(PetscLogEventBeginComplete, PetscLogEventEndComplete);
305:   return(0);
306: }

308: /*@C
309:   PetscLogTraceBegin - Activates trace logging.  Every time a PETSc event
310:   begins or ends, the event name is printed.

312:   Logically Collective on PETSC_COMM_WORLD

314:   Input Parameter:
315: . file - The file to print trace in (e.g. stdout)

317:   Options Database Key:
318: . -log_trace [filename] - Activates PetscLogTraceBegin()

320:   Notes:
321:   PetscLogTraceBegin() prints the processor number, the execution time (sec),
322:   then "Event begin:" or "Event end:" followed by the event name.

324:   PetscLogTraceBegin() allows tracing of all PETSc calls, which is useful
325:   to determine where a program is hanging without running in the
326:   debugger.  Can be used in conjunction with the -info option.

328:   Level: intermediate

330: .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogDefaultBegin()
331: @*/
332: PetscErrorCode  PetscLogTraceBegin(FILE *file)
333: {

337:   petsc_tracefile = file;

339:   PetscLogSet(PetscLogEventBeginTrace, PetscLogEventEndTrace);
340:   return(0);
341: }

343: /*@
344:   PetscLogActions - Determines whether actions are logged for the graphical viewer.

346:   Not Collective

348:   Input Parameter:
349: . flag - PETSC_TRUE if actions are to be logged

351:   Level: intermediate

353:   Note: Logging of actions continues to consume more memory as the program
354:   runs. Long running programs should consider turning this feature off.

356:   Options Database Keys:
357: . -log_exclude_actions - Turns off actions logging

359: .seealso: PetscLogStagePush(), PetscLogStagePop()
360: @*/
361: PetscErrorCode  PetscLogActions(PetscBool flag)
362: {
364:   petsc_logActions = flag;
365:   return(0);
366: }

368: /*@
369:   PetscLogObjects - Determines whether objects are logged for the graphical viewer.

371:   Not Collective

373:   Input Parameter:
374: . flag - PETSC_TRUE if objects are to be logged

376:   Level: intermediate

378:   Note: Logging of objects continues to consume more memory as the program
379:   runs. Long running programs should consider turning this feature off.

381:   Options Database Keys:
382: . -log_exclude_objects - Turns off objects logging

384: .seealso: PetscLogStagePush(), PetscLogStagePop()
385: @*/
386: PetscErrorCode  PetscLogObjects(PetscBool flag)
387: {
389:   petsc_logObjects = flag;
390:   return(0);
391: }

393: /*------------------------------------------------ Stage Functions --------------------------------------------------*/
394: /*@C
395:   PetscLogStageRegister - Attaches a character string name to a logging stage.

397:   Not Collective

399:   Input Parameter:
400: . sname - The name to associate with that stage

402:   Output Parameter:
403: . stage - The stage number

405:   Level: intermediate

407: .seealso: PetscLogStagePush(), PetscLogStagePop()
408: @*/
409: PetscErrorCode  PetscLogStageRegister(const char sname[],PetscLogStage *stage)
410: {
411:   PetscStageLog  stageLog;
412:   PetscLogEvent  event;

416:   PetscLogGetStageLog(&stageLog);
417:   PetscStageLogRegister(stageLog, sname, stage);
418:   /* Copy events already changed in the main stage, this sucks */
419:   PetscEventPerfLogEnsureSize(stageLog->stageInfo[*stage].eventLog, stageLog->eventLog->numEvents);
420:   for (event = 0; event < stageLog->eventLog->numEvents; event++) {
421:     PetscEventPerfInfoCopy(&stageLog->stageInfo[0].eventLog->eventInfo[event],&stageLog->stageInfo[*stage].eventLog->eventInfo[event]);
422:   }
423:   PetscClassPerfLogEnsureSize(stageLog->stageInfo[*stage].classLog, stageLog->classLog->numClasses);
424:   return(0);
425: }

427: /*@C
428:   PetscLogStagePush - This function pushes a stage on the stack.

430:   Not Collective

432:   Input Parameter:
433: . stage - The stage on which to log

435:   Usage:
436:   If the option -log_sumary is used to run the program containing the
437:   following code, then 2 sets of summary data will be printed during
438:   PetscFinalize().
439: .vb
440:       PetscInitialize(int *argc,char ***args,0,0);
441:       [stage 0 of code]
442:       PetscLogStagePush(1);
443:       [stage 1 of code]
444:       PetscLogStagePop();
445:       PetscBarrier(...);
446:       [more stage 0 of code]
447:       PetscFinalize();
448: .ve

450:   Notes:
451:   Use PetscLogStageRegister() to register a stage.

453:   Level: intermediate

455: .seealso: PetscLogStagePop(), PetscLogStageRegister(), PetscBarrier()
456: @*/
457: PetscErrorCode  PetscLogStagePush(PetscLogStage stage)
458: {
459:   PetscStageLog  stageLog;

463:   PetscLogGetStageLog(&stageLog);
464:   PetscStageLogPush(stageLog, stage);
465:   return(0);
466: }

468: /*@C
469:   PetscLogStagePop - This function pops a stage from the stack.

471:   Not Collective

473:   Usage:
474:   If the option -log_sumary is used to run the program containing the
475:   following code, then 2 sets of summary data will be printed during
476:   PetscFinalize().
477: .vb
478:       PetscInitialize(int *argc,char ***args,0,0);
479:       [stage 0 of code]
480:       PetscLogStagePush(1);
481:       [stage 1 of code]
482:       PetscLogStagePop();
483:       PetscBarrier(...);
484:       [more stage 0 of code]
485:       PetscFinalize();
486: .ve

488:   Notes:
489:   Use PetscLogStageRegister() to register a stage.

491:   Level: intermediate

493: .seealso: PetscLogStagePush(), PetscLogStageRegister(), PetscBarrier()
494: @*/
495: PetscErrorCode  PetscLogStagePop(void)
496: {
497:   PetscStageLog  stageLog;

501:   PetscLogGetStageLog(&stageLog);
502:   PetscStageLogPop(stageLog);
503:   return(0);
504: }

506: /*@
507:   PetscLogStageSetActive - Determines stage activity for PetscLogEventBegin() and PetscLogEventEnd().

509:   Not Collective

511:   Input Parameters:
512: + stage    - The stage
513: - isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE)

515:   Level: intermediate

517: .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage()
518: @*/
519: PetscErrorCode  PetscLogStageSetActive(PetscLogStage stage, PetscBool isActive)
520: {
521:   PetscStageLog  stageLog;

525:   PetscLogGetStageLog(&stageLog);
526:   PetscStageLogSetActive(stageLog, stage, isActive);
527:   return(0);
528: }

530: /*@
531:   PetscLogStageGetActive - Returns stage activity for PetscLogEventBegin() and PetscLogEventEnd().

533:   Not Collective

535:   Input Parameter:
536: . stage    - The stage

538:   Output Parameter:
539: . isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE)

541:   Level: intermediate

543: .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage()
544: @*/
545: PetscErrorCode  PetscLogStageGetActive(PetscLogStage stage, PetscBool  *isActive)
546: {
547:   PetscStageLog  stageLog;

551:   PetscLogGetStageLog(&stageLog);
552:   PetscStageLogGetActive(stageLog, stage, isActive);
553:   return(0);
554: }

556: /*@
557:   PetscLogStageSetVisible - Determines stage visibility in PetscLogView()

559:   Not Collective

561:   Input Parameters:
562: + stage     - The stage
563: - isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE)

565:   Level: intermediate

567: .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView()
568: @*/
569: PetscErrorCode  PetscLogStageSetVisible(PetscLogStage stage, PetscBool isVisible)
570: {
571:   PetscStageLog  stageLog;

575:   PetscLogGetStageLog(&stageLog);
576:   PetscStageLogSetVisible(stageLog, stage, isVisible);
577:   return(0);
578: }

580: /*@
581:   PetscLogStageGetVisible - Returns stage visibility in PetscLogView()

583:   Not Collective

585:   Input Parameter:
586: . stage     - The stage

588:   Output Parameter:
589: . isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE)

591:   Level: intermediate

593: .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView()
594: @*/
595: PetscErrorCode  PetscLogStageGetVisible(PetscLogStage stage, PetscBool  *isVisible)
596: {
597:   PetscStageLog  stageLog;

601:   PetscLogGetStageLog(&stageLog);
602:   PetscStageLogGetVisible(stageLog, stage, isVisible);
603:   return(0);
604: }

606: /*@C
607:   PetscLogStageGetId - Returns the stage id when given the stage name.

609:   Not Collective

611:   Input Parameter:
612: . name  - The stage name

614:   Output Parameter:
615: . stage - The stage, , or -1 if no stage with that name exists

617:   Level: intermediate

619: .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage()
620: @*/
621: PetscErrorCode  PetscLogStageGetId(const char name[], PetscLogStage *stage)
622: {
623:   PetscStageLog  stageLog;

627:   PetscLogGetStageLog(&stageLog);
628:   PetscStageLogGetStage(stageLog, name, stage);
629:   return(0);
630: }

632: /*------------------------------------------------ Event Functions --------------------------------------------------*/
633: /*@C
634:   PetscLogEventRegister - Registers an event name for logging operations in an application code.

636:   Not Collective

638:   Input Parameters:
639: + name   - The name associated with the event
640: - classid - The classid associated to the class for this event, obtain either with
641:            PetscClassIdRegister() or use a predefined one such as KSP_CLASSID, SNES_CLASSID, the predefined ones
642:            are only available in C code

644:   Output Parameter:
645: . event - The event id for use with PetscLogEventBegin() and PetscLogEventEnd().

647:   Example of Usage:
648: .vb
649:       PetscLogEvent USER_EVENT;
650:       PetscClassId classid;
651:       PetscLogDouble user_event_flops;
652:       PetscClassIdRegister("class name",&classid);
653:       PetscLogEventRegister("User event name",classid,&USER_EVENT);
654:       PetscLogEventBegin(USER_EVENT,0,0,0,0);
655:          [code segment to monitor]
656:          PetscLogFlops(user_event_flops);
657:       PetscLogEventEnd(USER_EVENT,0,0,0,0);
658: .ve

660:   Notes:
661:   PETSc automatically logs library events if the code has been
662:   configured with --with-log (which is the default) and
663:   -log_view or -log_all is specified.  PetscLogEventRegister() is
664:   intended for logging user events to supplement this PETSc
665:   information.

667:   PETSc can gather data for use with the utilities Jumpshot
668:   (part of the MPICH distribution).  If PETSc has been compiled
669:   with flag -DPETSC_HAVE_MPE (MPE is an additional utility within
670:   MPICH), the user can employ another command line option, -log_mpe,
671:   to create a logfile, "mpe.log", which can be visualized
672:   Jumpshot.

674:   The classid is associated with each event so that classes of events
675:   can be disabled simultaneously, such as all matrix events. The user
676:   can either use an existing classid, such as MAT_CLASSID, or create
677:   their own as shown in the example.

679:   If an existing event with the same name exists, its event handle is
680:   returned instead of creating a new event.

682:   Level: intermediate

684: .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogFlops(),
685:           PetscLogEventActivate(), PetscLogEventDeactivate(), PetscClassIdRegister()
686: @*/
687: PetscErrorCode  PetscLogEventRegister(const char name[],PetscClassId classid,PetscLogEvent *event)
688: {
689:   PetscStageLog  stageLog;
690:   int            stage;

694:   *event = PETSC_DECIDE;
695:   PetscLogGetStageLog(&stageLog);
696:   PetscEventRegLogGetEvent(stageLog->eventLog, name, event);
697:   if (*event > 0) return(0);
698:   PetscEventRegLogRegister(stageLog->eventLog, name, classid, event);
699:   for (stage = 0; stage < stageLog->numStages; stage++) {
700:     PetscEventPerfLogEnsureSize(stageLog->stageInfo[stage].eventLog, stageLog->eventLog->numEvents);
701:     PetscClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses);
702:   }
703:   return(0);
704: }

706: /*@
707:   PetscLogEventSetCollective - Indicates that a particular event is collective.

709:   Not Collective

711:   Input Parameters:
712: + event - The event id
713: - collective - Bolean flag indicating whether a particular event is collective

715:   Note:
716:   New events returned from PetscLogEventRegister() are collective by default.

718:   Level: developer

720: .seealso: PetscLogEventRegister()
721: @*/
722: PetscErrorCode PetscLogEventSetCollective(PetscLogEvent event,PetscBool collective)
723: {
724:   PetscStageLog    stageLog;
725:   PetscEventRegLog eventRegLog;
726:   PetscErrorCode   ierr;

729:   PetscLogGetStageLog(&stageLog);
730:   PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
731:   if (event < 0 || event > eventRegLog->numEvents) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid event id");
732:   eventRegLog->eventInfo[event].collective = collective;
733:   return(0);
734: }

736: /*@
737:   PetscLogEventIncludeClass - Activates event logging for a PETSc object class in every stage.

739:   Not Collective

741:   Input Parameter:
742: . classid - The object class, for example MAT_CLASSID, SNES_CLASSID, etc.

744:   Level: developer

746: .seealso: PetscLogEventActivateClass(),PetscLogEventDeactivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate()
747: @*/
748: PetscErrorCode  PetscLogEventIncludeClass(PetscClassId classid)
749: {
750:   PetscStageLog  stageLog;
751:   int            stage;

755:   PetscLogGetStageLog(&stageLog);
756:   for (stage = 0; stage < stageLog->numStages; stage++) {
757:     PetscEventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);
758:   }
759:   return(0);
760: }

762: /*@
763:   PetscLogEventExcludeClass - Deactivates event logging for a PETSc object class in every stage.

765:   Not Collective

767:   Input Parameter:
768: . classid - The object class, for example MAT_CLASSID, SNES_CLASSID, etc.

770:   Level: developer

772: .seealso: PetscLogEventDeactivateClass(),PetscLogEventActivateClass(),PetscLogEventDeactivate(),PetscLogEventActivate()
773: @*/
774: PetscErrorCode  PetscLogEventExcludeClass(PetscClassId classid)
775: {
776:   PetscStageLog  stageLog;
777:   int            stage;

781:   PetscLogGetStageLog(&stageLog);
782:   for (stage = 0; stage < stageLog->numStages; stage++) {
783:     PetscEventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);
784:   }
785:   return(0);
786: }

788: /*@
789:   PetscLogEventActivate - Indicates that a particular event should be logged.

791:   Not Collective

793:   Input Parameter:
794: . event - The event id

796:   Usage:
797: .vb
798:       PetscLogEventDeactivate(VEC_SetValues);
799:         [code where you do not want to log VecSetValues()]
800:       PetscLogEventActivate(VEC_SetValues);
801:         [code where you do want to log VecSetValues()]
802: .ve

804:   Note:
805:   The event may be either a pre-defined PETSc event (found in include/petsclog.h)
806:   or an event number obtained with PetscLogEventRegister().

808:   Level: advanced

810: .seealso: PlogEventDeactivate(), PlogEventDeactivatePush(), PetscLogEventDeactivatePop()
811: @*/
812: PetscErrorCode  PetscLogEventActivate(PetscLogEvent event)
813: {
814:   PetscStageLog  stageLog;
815:   int            stage;

819:   PetscLogGetStageLog(&stageLog);
820:   PetscStageLogGetCurrent(stageLog, &stage);
821:   PetscEventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event);
822:   return(0);
823: }

825: /*@
826:   PetscLogEventDeactivate - Indicates that a particular event should not be logged.

828:   Not Collective

830:   Input Parameter:
831: . event - The event id

833:   Usage:
834: .vb
835:       PetscLogEventDeactivate(VEC_SetValues);
836:         [code where you do not want to log VecSetValues()]
837:       PetscLogEventActivate(VEC_SetValues);
838:         [code where you do want to log VecSetValues()]
839: .ve

841:   Note:
842:   The event may be either a pre-defined PETSc event (found in
843:   include/petsclog.h) or an event number obtained with PetscLogEventRegister()).

845:   Level: advanced

847: .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePush(), PetscLogEventDeactivatePop()
848: @*/
849: PetscErrorCode  PetscLogEventDeactivate(PetscLogEvent event)
850: {
851:   PetscStageLog  stageLog;
852:   int            stage;

856:   PetscLogGetStageLog(&stageLog);
857:   PetscStageLogGetCurrent(stageLog, &stage);
858:   PetscEventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event);
859:   return(0);
860: }

862: /*@
863:   PetscLogEventDeactivatePush - Indicates that a particular event should not be logged.

865:   Not Collective

867:   Input Parameter:
868: . event - The event id

870:   Usage:
871: .vb
872:       PetscLogEventDeactivatePush(VEC_SetValues);
873:         [code where you do not want to log VecSetValues()]
874:       PetscLogEventDeactivatePop(VEC_SetValues);
875:         [code where you do want to log VecSetValues()]
876: .ve

878:   Note:
879:   The event may be either a pre-defined PETSc event (found in
880:   include/petsclog.h) or an event number obtained with PetscLogEventRegister()).

882:   Level: advanced

884: .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePop()
885: @*/
886: PetscErrorCode  PetscLogEventDeactivatePush(PetscLogEvent event)
887: {
888:   PetscStageLog  stageLog;
889:   int            stage;

893:   PetscLogGetStageLog(&stageLog);
894:   PetscStageLogGetCurrent(stageLog, &stage);
895:   PetscEventPerfLogDeactivatePush(stageLog->stageInfo[stage].eventLog, event);
896:   return(0);
897: }

899: /*@
900:   PetscLogEventDeactivatePop - Indicates that a particular event shouldbe logged.

902:   Not Collective

904:   Input Parameter:
905: . event - The event id

907:   Usage:
908: .vb
909:       PetscLogEventDeactivatePush(VEC_SetValues);
910:         [code where you do not want to log VecSetValues()]
911:       PetscLogEventDeactivatePop(VEC_SetValues);
912:         [code where you do want to log VecSetValues()]
913: .ve

915:   Note:
916:   The event may be either a pre-defined PETSc event (found in
917:   include/petsclog.h) or an event number obtained with PetscLogEventRegister()).

919:   Level: advanced

921: .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePush()
922: @*/
923: PetscErrorCode  PetscLogEventDeactivatePop(PetscLogEvent event)
924: {
925:   PetscStageLog  stageLog;
926:   int            stage;

930:   PetscLogGetStageLog(&stageLog);
931:   PetscStageLogGetCurrent(stageLog, &stage);
932:   PetscEventPerfLogDeactivatePop(stageLog->stageInfo[stage].eventLog, event);
933:   return(0);
934: }

936: /*@
937:   PetscLogEventSetActiveAll - Sets the event activity in every stage.

939:   Not Collective

941:   Input Parameters:
942: + event    - The event id
943: - isActive - The activity flag determining whether the event is logged

945:   Level: advanced

947: .seealso: PlogEventActivate(),PlogEventDeactivate()
948: @*/
949: PetscErrorCode  PetscLogEventSetActiveAll(PetscLogEvent event, PetscBool isActive)
950: {
951:   PetscStageLog  stageLog;
952:   int            stage;

956:   PetscLogGetStageLog(&stageLog);
957:   for (stage = 0; stage < stageLog->numStages; stage++) {
958:     if (isActive) {
959:       PetscEventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event);
960:     } else {
961:       PetscEventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event);
962:     }
963:   }
964:   return(0);
965: }

967: /*@
968:   PetscLogEventActivateClass - Activates event logging for a PETSc object class.

970:   Not Collective

972:   Input Parameter:
973: . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc.

975:   Level: developer

977: .seealso: PetscLogEventDeactivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate()
978: @*/
979: PetscErrorCode  PetscLogEventActivateClass(PetscClassId classid)
980: {
981:   PetscStageLog  stageLog;
982:   int            stage;

986:   PetscLogGetStageLog(&stageLog);
987:   PetscStageLogGetCurrent(stageLog, &stage);
988:   PetscEventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);
989:   return(0);
990: }

992: /*@
993:   PetscLogEventDeactivateClass - Deactivates event logging for a PETSc object class.

995:   Not Collective

997:   Input Parameter:
998: . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc.

1000:   Level: developer

1002: .seealso: PetscLogEventActivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate()
1003: @*/
1004: PetscErrorCode  PetscLogEventDeactivateClass(PetscClassId classid)
1005: {
1006:   PetscStageLog  stageLog;
1007:   int            stage;

1011:   PetscLogGetStageLog(&stageLog);
1012:   PetscStageLogGetCurrent(stageLog, &stage);
1013:   PetscEventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);
1014:   return(0);
1015: }

1017: /*MC
1018:    PetscLogEventSync - Synchronizes the beginning of a user event.

1020:    Synopsis:
1021: #include <petsclog.h>
1022:    PetscErrorCode PetscLogEventSync(int e,MPI_Comm comm)

1024:    Collective

1026:    Input Parameters:
1027: +  e - integer associated with the event obtained from PetscLogEventRegister()
1028: -  comm - an MPI communicator

1030:    Usage:
1031: .vb
1032:      PetscLogEvent USER_EVENT;
1033:      PetscLogEventRegister("User event",0,&USER_EVENT);
1034:      PetscLogEventSync(USER_EVENT,PETSC_COMM_WORLD);
1035:      PetscLogEventBegin(USER_EVENT,0,0,0,0);
1036:         [code segment to monitor]
1037:      PetscLogEventEnd(USER_EVENT,0,0,0,0);
1038: .ve

1040:    Notes:
1041:    This routine should be called only if there is not a
1042:    PetscObject available to pass to PetscLogEventBegin().

1044:    Level: developer

1046: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd()

1048: M*/

1050: /*MC
1051:    PetscLogEventBegin - Logs the beginning of a user event.

1053:    Synopsis:
1054: #include <petsclog.h>
1055:    PetscErrorCode PetscLogEventBegin(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)

1057:    Not Collective

1059:    Input Parameters:
1060: +  e - integer associated with the event obtained from PetscLogEventRegister()
1061: -  o1,o2,o3,o4 - objects associated with the event, or 0

1063:    Fortran Synopsis:
1064:    void PetscLogEventBegin(int e,PetscErrorCode ierr)

1066:    Usage:
1067: .vb
1068:      PetscLogEvent USER_EVENT;
1069:      PetscLogDouble user_event_flops;
1070:      PetscLogEventRegister("User event",0,&USER_EVENT);
1071:      PetscLogEventBegin(USER_EVENT,0,0,0,0);
1072:         [code segment to monitor]
1073:         PetscLogFlops(user_event_flops);
1074:      PetscLogEventEnd(USER_EVENT,0,0,0,0);
1075: .ve

1077:    Notes:
1078:    You need to register each integer event with the command
1079:    PetscLogEventRegister().

1081:    Level: intermediate

1083: .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops()

1085: M*/

1087: /*MC
1088:    PetscLogEventEnd - Log the end of a user event.

1090:    Synopsis:
1091: #include <petsclog.h>
1092:    PetscErrorCode PetscLogEventEnd(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)

1094:    Not Collective

1096:    Input Parameters:
1097: +  e - integer associated with the event obtained with PetscLogEventRegister()
1098: -  o1,o2,o3,o4 - objects associated with the event, or 0

1100:    Fortran Synopsis:
1101:    void PetscLogEventEnd(int e,PetscErrorCode ierr)

1103:    Usage:
1104: .vb
1105:      PetscLogEvent USER_EVENT;
1106:      PetscLogDouble user_event_flops;
1107:      PetscLogEventRegister("User event",0,&USER_EVENT,);
1108:      PetscLogEventBegin(USER_EVENT,0,0,0,0);
1109:         [code segment to monitor]
1110:         PetscLogFlops(user_event_flops);
1111:      PetscLogEventEnd(USER_EVENT,0,0,0,0);
1112: .ve

1114:    Notes:
1115:    You should also register each additional integer event with the command
1116:    PetscLogEventRegister().

1118:    Level: intermediate

1120: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogFlops()

1122: M*/

1124: /*@C
1125:   PetscLogEventGetId - Returns the event id when given the event name.

1127:   Not Collective

1129:   Input Parameter:
1130: . name  - The event name

1132:   Output Parameter:
1133: . event - The event, or -1 if no event with that name exists

1135:   Level: intermediate

1137: .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogStageGetId()
1138: @*/
1139: PetscErrorCode  PetscLogEventGetId(const char name[], PetscLogEvent *event)
1140: {
1141:   PetscStageLog  stageLog;

1145:   PetscLogGetStageLog(&stageLog);
1146:   PetscEventRegLogGetEvent(stageLog->eventLog, name, event);
1147:   return(0);
1148: }

1150: /*------------------------------------------------ Output Functions -------------------------------------------------*/
1151: /*@C
1152:   PetscLogDump - Dumps logs of objects to a file. This file is intended to
1153:   be read by bin/petscview. This program no longer exists.

1155:   Collective on PETSC_COMM_WORLD

1157:   Input Parameter:
1158: . name - an optional file name

1160:   Usage:
1161: .vb
1162:      PetscInitialize(...);
1163:      PetscLogDefaultBegin(); or PetscLogAllBegin();
1164:      ... code ...
1165:      PetscLogDump(filename);
1166:      PetscFinalize();
1167: .ve

1169:   Notes:
1170:   The default file name is
1171: $    Log.<rank>
1172:   where <rank> is the processor number. If no name is specified,
1173:   this file will be used.

1175:   Level: advanced

1177: .seealso: PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogView()
1178: @*/
1179: PetscErrorCode  PetscLogDump(const char sname[])
1180: {
1181:   PetscStageLog      stageLog;
1182:   PetscEventPerfInfo *eventInfo;
1183:   FILE               *fd;
1184:   char               file[PETSC_MAX_PATH_LEN], fname[PETSC_MAX_PATH_LEN];
1185:   PetscLogDouble     flops, _TotalTime;
1186:   PetscMPIInt        rank;
1187:   int                action, object, curStage;
1188:   PetscLogEvent      event;
1189:   PetscErrorCode     ierr;

1192:   /* Calculate the total elapsed time */
1193:   PetscTime(&_TotalTime);
1194:   _TotalTime -= petsc_BaseTime;
1195:   /* Open log file */
1196:   MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
1197:   if (sname && sname[0]) sprintf(file, "%s.%d", sname, rank);
1198:   else sprintf(file, "Log.%d", rank);
1199:   PetscFixFilename(file, fname);
1200:   PetscFOpen(PETSC_COMM_WORLD, fname, "w", &fd);
1201:   if ((rank == 0) && (!fd)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN, "Cannot open file: %s", fname);
1202:   /* Output totals */
1203:   PetscFPrintf(PETSC_COMM_WORLD, fd, "Total Flop %14e %16.8e\n", petsc_TotalFlops, _TotalTime);
1204:   PetscFPrintf(PETSC_COMM_WORLD, fd, "Clock Resolution %g\n", 0.0);
1205:   /* Output actions */
1206:   if (petsc_logActions) {
1207:     PetscFPrintf(PETSC_COMM_WORLD, fd, "Actions accomplished %d\n", petsc_numActions);
1208:     for (action = 0; action < petsc_numActions; action++) {
1209:       PetscFPrintf(PETSC_COMM_WORLD, fd, "%g %d %d %d %d %d %d %g %g %g\n",
1210:                           petsc_actions[action].time, petsc_actions[action].action, (int)petsc_actions[action].event, (int)petsc_actions[action].classid, petsc_actions[action].id1,
1211:                           petsc_actions[action].id2, petsc_actions[action].id3, petsc_actions[action].flops, petsc_actions[action].mem, petsc_actions[action].maxmem);
1212:     }
1213:   }
1214:   /* Output objects */
1215:   if (petsc_logObjects) {
1216:     PetscFPrintf(PETSC_COMM_WORLD, fd, "Objects created %d destroyed %d\n", petsc_numObjects, petsc_numObjectsDestroyed);
1217:     for (object = 0; object < petsc_numObjects; object++) {
1218:       PetscFPrintf(PETSC_COMM_WORLD, fd, "Parent ID: %d Memory: %d\n", petsc_objects[object].parent, (int) petsc_objects[object].mem);
1219:       if (!petsc_objects[object].name[0]) {
1220:         PetscFPrintf(PETSC_COMM_WORLD, fd,"No Name\n");
1221:       } else {
1222:         PetscFPrintf(PETSC_COMM_WORLD, fd, "Name: %s\n", petsc_objects[object].name);
1223:       }
1224:       if (petsc_objects[object].info[0] != 0) {
1225:         PetscFPrintf(PETSC_COMM_WORLD, fd, "No Info\n");
1226:       } else {
1227:         PetscFPrintf(PETSC_COMM_WORLD, fd, "Info: %s\n", petsc_objects[object].info);
1228:       }
1229:     }
1230:   }
1231:   /* Output events */
1232:   PetscFPrintf(PETSC_COMM_WORLD, fd, "Event log:\n");
1233:   PetscLogGetStageLog(&stageLog);
1234:   PetscIntStackTop(stageLog->stack, &curStage);
1235:   eventInfo = stageLog->stageInfo[curStage].eventLog->eventInfo;
1236:   for (event = 0; event < stageLog->stageInfo[curStage].eventLog->numEvents; event++) {
1237:     if (eventInfo[event].time != 0.0) flops = eventInfo[event].flops/eventInfo[event].time;
1238:     else flops = 0.0;
1239:     PetscFPrintf(PETSC_COMM_WORLD, fd, "%d %16d %16g %16g %16g\n", event, eventInfo[event].count,
1240:                         eventInfo[event].flops, eventInfo[event].time, flops);
1241:   }
1242:   PetscFClose(PETSC_COMM_WORLD, fd);
1243:   return(0);
1244: }

1246: /*
1247:   PetscLogView_Detailed - Each process prints the times for its own events

1249: */
1250: PetscErrorCode  PetscLogView_Detailed(PetscViewer viewer)
1251: {
1252:   PetscStageLog      stageLog;
1253:   PetscEventPerfInfo *eventInfo = NULL, *stageInfo = NULL;
1254:   PetscLogDouble     locTotalTime, numRed, maxMem;
1255:   int                numStages,numEvents,stage,event;
1256:   MPI_Comm           comm = PetscObjectComm((PetscObject) viewer);
1257:   PetscMPIInt        rank,size;
1258:   PetscErrorCode     ierr;

1261:   MPI_Comm_size(comm, &size);
1262:   MPI_Comm_rank(comm, &rank);
1263:   /* Must preserve reduction count before we go on */
1264:   numRed = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
1265:   /* Get the total elapsed time */
1266:   PetscTime(&locTotalTime);  locTotalTime -= petsc_BaseTime;
1267:   PetscViewerASCIIPrintf(viewer,"size = %d\n",size);
1268:   PetscViewerASCIIPrintf(viewer,"LocalTimes = {}\n");
1269:   PetscViewerASCIIPrintf(viewer,"LocalMessages = {}\n");
1270:   PetscViewerASCIIPrintf(viewer,"LocalMessageLens = {}\n");
1271:   PetscViewerASCIIPrintf(viewer,"LocalReductions = {}\n");
1272:   PetscViewerASCIIPrintf(viewer,"LocalFlop = {}\n");
1273:   PetscViewerASCIIPrintf(viewer,"LocalObjects = {}\n");
1274:   PetscViewerASCIIPrintf(viewer,"LocalMemory = {}\n");
1275:   PetscLogGetStageLog(&stageLog);
1276:   MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);
1277:   PetscViewerASCIIPrintf(viewer,"Stages = {}\n");
1278:   for (stage=0; stage<numStages; stage++) {
1279:     PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"] = {}\n",stageLog->stageInfo[stage].name);
1280:     PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"summary\"] = {}\n",stageLog->stageInfo[stage].name);
1281:     MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);
1282:     for (event = 0; event < numEvents; event++) {
1283:       PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"%s\"] = {}\n",stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name);
1284:     }
1285:   }
1286:   PetscMallocGetMaximumUsage(&maxMem);
1287:   PetscViewerASCIIPushSynchronized(viewer);
1288:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalTimes[%d] = %g\n",rank,locTotalTime);
1289:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessages[%d] = %g\n",rank,(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct));
1290:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessageLens[%d] = %g\n",rank,(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len));
1291:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalReductions[%d] = %g\n",rank,numRed);
1292:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalFlop[%d] = %g\n",rank,petsc_TotalFlops);
1293:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalObjects[%d] = %d\n",rank,petsc_numObjects);
1294:   PetscViewerASCIISynchronizedPrintf(viewer,"LocalMemory[%d] = %g\n",rank,maxMem);
1295:   PetscViewerFlush(viewer);
1296:   for (stage=0; stage<numStages; stage++) {
1297:     stageInfo = &stageLog->stageInfo[stage].perfInfo;
1298:     PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"summary\"][%d] = {\"time\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g}\n",
1299:                                               stageLog->stageInfo[stage].name,rank,
1300:                                               stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops);
1301:     MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);
1302:     for (event = 0; event < numEvents; event++) {
1303:       eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event];
1304:       PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"%s\"][%d] = {\"count\" : %D, \"time\" : %g, \"syncTime\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g",
1305:                                                 stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name,rank,
1306:                                                 eventInfo->count,eventInfo->time,eventInfo->syncTime,eventInfo->numMessages,eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops);
1307:       if (eventInfo->dof[0] >= 0.) {
1308:         PetscInt d, e;

1310:         PetscViewerASCIISynchronizedPrintf(viewer, ", \"dof\" : [");
1311:         for (d = 0; d < 8; ++d) {
1312:           if (d > 0) {PetscViewerASCIISynchronizedPrintf(viewer, ", ");}
1313:           PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->dof[d]);
1314:         }
1315:         PetscViewerASCIISynchronizedPrintf(viewer, "]");
1316:         PetscViewerASCIISynchronizedPrintf(viewer, ", \"error\" : [");
1317:         for (e = 0; e < 8; ++e) {
1318:           if (e > 0) {PetscViewerASCIISynchronizedPrintf(viewer, ", ");}
1319:           PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->errors[e]);
1320:         }
1321:         PetscViewerASCIISynchronizedPrintf(viewer, "]");
1322:       }
1323:       PetscViewerASCIISynchronizedPrintf(viewer,"}\n");
1324:     }
1325:   }
1326:   PetscViewerFlush(viewer);
1327:   PetscViewerASCIIPopSynchronized(viewer);
1328:   return(0);
1329: }

1331: /*
1332:   PetscLogView_CSV - Each process prints the times for its own events in Comma-Separated Value Format
1333: */
1334: PetscErrorCode  PetscLogView_CSV(PetscViewer viewer)
1335: {
1336:   PetscStageLog      stageLog;
1337:   PetscEventPerfInfo *eventInfo = NULL;
1338:   PetscLogDouble     locTotalTime, maxMem;
1339:   int                numStages,numEvents,stage,event;
1340:   MPI_Comm           comm = PetscObjectComm((PetscObject) viewer);
1341:   PetscMPIInt        rank,size;
1342:   PetscErrorCode     ierr;

1345:   MPI_Comm_size(comm, &size);
1346:   MPI_Comm_rank(comm, &rank);
1347:   /* Must preserve reduction count before we go on */
1348:   /* Get the total elapsed time */
1349:   PetscTime(&locTotalTime);  locTotalTime -= petsc_BaseTime;
1350:   PetscLogGetStageLog(&stageLog);
1351:   MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);
1352:   PetscMallocGetMaximumUsage(&maxMem);
1353:   PetscViewerASCIIPushSynchronized(viewer);
1354:   PetscViewerASCIIPrintf(viewer,"Stage Name,Event Name,Rank,Count,Time,Num Messages,Message Length,Num Reductions,FLOP,dof0,dof1,dof2,dof3,dof4,dof5,dof6,dof7,e0,e1,e2,e3,e4,e5,e6,e7,%d\n", size);
1355:   PetscViewerFlush(viewer);
1356:   for (stage=0; stage<numStages; stage++) {
1357:     PetscEventPerfInfo *stageInfo = &stageLog->stageInfo[stage].perfInfo;

1359:     PetscViewerASCIISynchronizedPrintf(viewer,"%s,summary,%d,1,%g,%g,%g,%g,%g\n",
1360:                                               stageLog->stageInfo[stage].name,rank,stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops);
1361:     MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);
1362:     for (event = 0; event < numEvents; event++) {
1363:       eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event];
1364:       PetscViewerASCIISynchronizedPrintf(viewer,"%s,%s,%d,%D,%g,%g,%g,%g,%g",stageLog->stageInfo[stage].name,
1365:                                                 stageLog->eventLog->eventInfo[event].name,rank,eventInfo->count,eventInfo->time,eventInfo->numMessages,
1366:                                                 eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops);
1367:       if (eventInfo->dof[0] >= 0.) {
1368:         PetscInt d, e;

1370:         for (d = 0; d < 8; ++d) {
1371:           PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->dof[d]);
1372:         }
1373:         for (e = 0; e < 8; ++e) {
1374:           PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->errors[e]);
1375:         }
1376:       }
1377:       PetscViewerASCIISynchronizedPrintf(viewer,"\n");
1378:     }
1379:   }
1380:   PetscViewerFlush(viewer);
1381:   PetscViewerASCIIPopSynchronized(viewer);
1382:   return(0);
1383: }

1385: static PetscErrorCode PetscLogViewWarnSync(MPI_Comm comm,FILE *fd)
1386: {
1389:   if (!PetscLogSyncOn) return(0);
1390:   PetscFPrintf(comm, fd, "\n\n");
1391:   PetscFPrintf(comm, fd, "      ##########################################################\n");
1392:   PetscFPrintf(comm, fd, "      #                                                        #\n");
1393:   PetscFPrintf(comm, fd, "      #                       WARNING!!!                       #\n");
1394:   PetscFPrintf(comm, fd, "      #                                                        #\n");
1395:   PetscFPrintf(comm, fd, "      #   This program was run with logging synchronization.   #\n");
1396:   PetscFPrintf(comm, fd, "      #   This option provides more meaningful imbalance       #\n");
1397:   PetscFPrintf(comm, fd, "      #   figures at the expense of slowing things down and    #\n");
1398:   PetscFPrintf(comm, fd, "      #   providing a distorted view of the overall runtime.   #\n");
1399:   PetscFPrintf(comm, fd, "      #                                                        #\n");
1400:   PetscFPrintf(comm, fd, "      ##########################################################\n\n\n");
1401:   return(0);
1402: }

1404: static PetscErrorCode PetscLogViewWarnDebugging(MPI_Comm comm,FILE *fd)
1405: {

1409:   if (PetscDefined(USE_DEBUG)) {
1410:     PetscFPrintf(comm, fd, "\n\n");
1411:     PetscFPrintf(comm, fd, "      ##########################################################\n");
1412:     PetscFPrintf(comm, fd, "      #                                                        #\n");
1413:     PetscFPrintf(comm, fd, "      #                       WARNING!!!                       #\n");
1414:     PetscFPrintf(comm, fd, "      #                                                        #\n");
1415:     PetscFPrintf(comm, fd, "      #   This code was compiled with a debugging option.      #\n");
1416:     PetscFPrintf(comm, fd, "      #   To get timing results run ./configure                #\n");
1417:     PetscFPrintf(comm, fd, "      #   using --with-debugging=no, the performance will      #\n");
1418:     PetscFPrintf(comm, fd, "      #   be generally two or three times faster.              #\n");
1419:     PetscFPrintf(comm, fd, "      #                                                        #\n");
1420:     PetscFPrintf(comm, fd, "      ##########################################################\n\n\n");
1421:   }
1422:   return(0);
1423: }

1425: static PetscErrorCode PetscLogViewWarnNoGpuAwareMpi(MPI_Comm comm,FILE *fd)
1426: {
1427: #if defined(PETSC_HAVE_DEVICE)
1429:   PetscMPIInt    size;

1432:   MPI_Comm_size(comm, &size);
1433:   if (use_gpu_aware_mpi || size == 1 || !PetscCreatedGpuObjects) return(0);
1434:   PetscFPrintf(comm, fd, "\n\n");
1435:   PetscFPrintf(comm, fd, "      ##########################################################\n");
1436:   PetscFPrintf(comm, fd, "      #                                                        #\n");
1437:   PetscFPrintf(comm, fd, "      #                       WARNING!!!                       #\n");
1438:   PetscFPrintf(comm, fd, "      #                                                        #\n");
1439:   PetscFPrintf(comm, fd, "      #   This code was compiled with GPU support and you've   #\n");
1440:   PetscFPrintf(comm, fd, "      #   created PETSc/GPU objects, but you intentionally     #\n");
1441:   PetscFPrintf(comm, fd, "      #   used -use_gpu_aware_mpi 0, requiring PETSc to copy   #\n");
1442:   PetscFPrintf(comm, fd, "      #   additional data between the GPU and CPU. To obtain   #\n");
1443:   PetscFPrintf(comm, fd, "      #   meaningful timing results on multi-rank runs, use    #\n");
1444:   PetscFPrintf(comm, fd, "      #   GPU-aware MPI instead.                               #\n");
1445:   PetscFPrintf(comm, fd, "      #                                                        #\n");
1446:   PetscFPrintf(comm, fd, "      ##########################################################\n\n\n");
1447:   return(0);
1448: #else
1449:   return 0;
1450: #endif
1451: }

1453: PetscErrorCode  PetscLogView_Default(PetscViewer viewer)
1454: {
1455:   FILE               *fd;
1456:   PetscLogDouble     zero       = 0.0;
1457:   PetscStageLog      stageLog;
1458:   PetscStageInfo     *stageInfo = NULL;
1459:   PetscEventPerfInfo *eventInfo = NULL;
1460:   PetscClassPerfInfo *classInfo;
1461:   char               arch[128],hostname[128],username[128],pname[PETSC_MAX_PATH_LEN],date[128];
1462:   const char         *name;
1463:   PetscLogDouble     locTotalTime, TotalTime, TotalFlops;
1464:   PetscLogDouble     numMessages, messageLength, avgMessLen, numReductions;
1465:   PetscLogDouble     stageTime, flops, flopr, mem, mess, messLen, red;
1466:   PetscLogDouble     fracTime, fracFlops, fracMessages, fracLength, fracReductions, fracMess, fracMessLen, fracRed;
1467:   PetscLogDouble     fracStageTime, fracStageFlops, fracStageMess, fracStageMessLen, fracStageRed;
1468:   PetscLogDouble     min, max, tot, ratio, avg, x, y;
1469:   PetscLogDouble     minf, maxf, totf, ratf, mint, maxt, tott, ratt, ratC, totm, totml, totr, mal, malmax, emalmax;
1470:   #if defined(PETSC_HAVE_DEVICE)
1471:   PetscLogDouble     cct, gct, csz, gsz, gmaxt, gflops, gflopr, fracgflops;
1472:   #endif
1473:   PetscMPIInt        minC, maxC;
1474:   PetscMPIInt        size, rank;
1475:   PetscBool          *localStageUsed,    *stageUsed;
1476:   PetscBool          *localStageVisible, *stageVisible;
1477:   int                numStages, localNumEvents, numEvents;
1478:   int                stage, oclass;
1479:   PetscLogEvent      event;
1480:   PetscErrorCode     ierr;
1481:   char               version[256];
1482:   MPI_Comm           comm;

1485:   PetscObjectGetComm((PetscObject)viewer,&comm);
1486:   PetscViewerASCIIGetPointer(viewer,&fd);
1487:   MPI_Comm_size(comm, &size);
1488:   MPI_Comm_rank(comm, &rank);
1489:   /* Get the total elapsed time */
1490:   PetscTime(&locTotalTime);  locTotalTime -= petsc_BaseTime;

1492:   PetscFPrintf(comm, fd, "************************************************************************************************************************\n");
1493:   PetscFPrintf(comm, fd, "***             WIDEN YOUR WINDOW TO 120 CHARACTERS.  Use 'enscript -r -fCourier9' to print this document            ***\n");
1494:   PetscFPrintf(comm, fd, "************************************************************************************************************************\n");
1495:   PetscFPrintf(comm, fd, "\n---------------------------------------------- PETSc Performance Summary: ----------------------------------------------\n\n");
1496:   PetscLogViewWarnSync(comm,fd);
1497:   PetscLogViewWarnDebugging(comm,fd);
1498:   PetscLogViewWarnNoGpuAwareMpi(comm,fd);
1499:   PetscGetArchType(arch,sizeof(arch));
1500:   PetscGetHostName(hostname,sizeof(hostname));
1501:   PetscGetUserName(username,sizeof(username));
1502:   PetscGetProgramName(pname,sizeof(pname));
1503:   PetscGetDate(date,sizeof(date));
1504:   PetscGetVersion(version,sizeof(version));
1505:   if (size == 1) {
1506:     PetscFPrintf(comm,fd,"%s on a %s named %s with %d processor, by %s %s\n", pname, arch, hostname, size, username, date);
1507:   } else {
1508:     PetscFPrintf(comm,fd,"%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date);
1509:   }
1510: #if defined(PETSC_HAVE_OPENMP)
1511:   PetscFPrintf(comm,fd,"Using %D OpenMP threads\n", PetscNumOMPThreads);
1512: #endif
1513:   PetscFPrintf(comm, fd, "Using %s\n", version);

1515:   /* Must preserve reduction count before we go on */
1516:   red = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;

1518:   /* Calculate summary information */
1519:   PetscFPrintf(comm, fd, "\n                         Max       Max/Min     Avg       Total\n");
1520:   /*   Time */
1521:   MPI_Allreduce(&locTotalTime, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1522:   MPI_Allreduce(&locTotalTime, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1523:   MPI_Allreduce(&locTotalTime, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1524:   avg  = tot/((PetscLogDouble) size);
1525:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1526:   PetscFPrintf(comm, fd, "Time (sec):           %5.3e   %7.3f   %5.3e\n", max, ratio, avg);
1527:   TotalTime = tot;
1528:   /*   Objects */
1529:   avg  = (PetscLogDouble) petsc_numObjects;
1530:   MPI_Allreduce(&avg,          &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1531:   MPI_Allreduce(&avg,          &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1532:   MPI_Allreduce(&avg,          &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1533:   avg  = tot/((PetscLogDouble) size);
1534:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1535:   PetscFPrintf(comm, fd, "Objects:              %5.3e   %7.3f   %5.3e\n", max, ratio, avg);
1536:   /*   Flops */
1537:   MPI_Allreduce(&petsc_TotalFlops,  &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1538:   MPI_Allreduce(&petsc_TotalFlops,  &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1539:   MPI_Allreduce(&petsc_TotalFlops,  &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1540:   avg  = tot/((PetscLogDouble) size);
1541:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1542:   PetscFPrintf(comm, fd, "Flop:                 %5.3e   %7.3f   %5.3e  %5.3e\n", max, ratio, avg, tot);
1543:   TotalFlops = tot;
1544:   /*   Flops/sec -- Must talk to Barry here */
1545:   if (locTotalTime != 0.0) flops = petsc_TotalFlops/locTotalTime; else flops = 0.0;
1546:   MPI_Allreduce(&flops,        &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1547:   MPI_Allreduce(&flops,        &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1548:   MPI_Allreduce(&flops,        &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1549:   avg  = tot/((PetscLogDouble) size);
1550:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1551:   PetscFPrintf(comm, fd, "Flop/sec:             %5.3e   %7.3f   %5.3e  %5.3e\n", max, ratio, avg, tot);
1552:   /*   Memory */
1553:   PetscMallocGetMaximumUsage(&mem);
1554:   if (mem > 0.0) {
1555:     MPI_Allreduce(&mem,          &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1556:     MPI_Allreduce(&mem,          &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1557:     MPI_Allreduce(&mem,          &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1558:     avg  = tot/((PetscLogDouble) size);
1559:     if (min != 0.0) ratio = max/min; else ratio = 0.0;
1560:     PetscFPrintf(comm, fd, "Memory:               %5.3e   %7.3f   %5.3e  %5.3e\n", max, ratio, avg, tot);
1561:   }
1562:   /*   Messages */
1563:   mess = 0.5*(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct);
1564:   MPI_Allreduce(&mess,         &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1565:   MPI_Allreduce(&mess,         &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1566:   MPI_Allreduce(&mess,         &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1567:   avg  = tot/((PetscLogDouble) size);
1568:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1569:   PetscFPrintf(comm, fd, "MPI Messages:         %5.3e   %7.3f   %5.3e  %5.3e\n", max, ratio, avg, tot);
1570:   numMessages = tot;
1571:   /*   Message Lengths */
1572:   mess = 0.5*(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len);
1573:   MPI_Allreduce(&mess,         &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1574:   MPI_Allreduce(&mess,         &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1575:   MPI_Allreduce(&mess,         &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1576:   if (numMessages != 0) avg = tot/numMessages; else avg = 0.0;
1577:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1578:   PetscFPrintf(comm, fd, "MPI Message Lengths:  %5.3e   %7.3f   %5.3e  %5.3e\n", max, ratio, avg, tot);
1579:   messageLength = tot;
1580:   /*   Reductions */
1581:   MPI_Allreduce(&red,          &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1582:   MPI_Allreduce(&red,          &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1583:   MPI_Allreduce(&red,          &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1584:   if (min != 0.0) ratio = max/min; else ratio = 0.0;
1585:   PetscFPrintf(comm, fd, "MPI Reductions:       %5.3e   %7.3f\n", max, ratio);
1586:   numReductions = red; /* wrong because uses count from process zero */
1587:   PetscFPrintf(comm, fd, "\nFlop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)\n");
1588:   PetscFPrintf(comm, fd, "                            e.g., VecAXPY() for real vectors of length N --> 2N flop\n");
1589:   PetscFPrintf(comm, fd, "                            and VecAXPY() for complex vectors of length N --> 8N flop\n");

1591:   /* Get total number of stages --
1592:        Currently, a single processor can register more stages than another, but stages must all be registered in order.
1593:        We can removed this requirement if necessary by having a global stage numbering and indirection on the stage ID.
1594:        This seems best accomplished by assoicating a communicator with each stage.
1595:   */
1596:   PetscLogGetStageLog(&stageLog);
1597:   MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);
1598:   PetscMalloc1(numStages, &localStageUsed);
1599:   PetscMalloc1(numStages, &stageUsed);
1600:   PetscMalloc1(numStages, &localStageVisible);
1601:   PetscMalloc1(numStages, &stageVisible);
1602:   if (numStages > 0) {
1603:     stageInfo = stageLog->stageInfo;
1604:     for (stage = 0; stage < numStages; stage++) {
1605:       if (stage < stageLog->numStages) {
1606:         localStageUsed[stage]    = stageInfo[stage].used;
1607:         localStageVisible[stage] = stageInfo[stage].perfInfo.visible;
1608:       } else {
1609:         localStageUsed[stage]    = PETSC_FALSE;
1610:         localStageVisible[stage] = PETSC_TRUE;
1611:       }
1612:     }
1613:     MPI_Allreduce(localStageUsed,    stageUsed,    numStages, MPIU_BOOL, MPI_LOR,  comm);
1614:     MPI_Allreduce(localStageVisible, stageVisible, numStages, MPIU_BOOL, MPI_LAND, comm);
1615:     for (stage = 0; stage < numStages; stage++) {
1616:       if (stageUsed[stage]) {
1617:         PetscFPrintf(comm, fd, "\nSummary of Stages:   ----- Time ------  ----- Flop ------  --- Messages ---  -- Message Lengths --  -- Reductions --\n");
1618:         PetscFPrintf(comm, fd, "                        Avg     %%Total     Avg     %%Total    Count   %%Total     Avg         %%Total    Count   %%Total\n");
1619:         break;
1620:       }
1621:     }
1622:     for (stage = 0; stage < numStages; stage++) {
1623:       if (!stageUsed[stage]) continue;
1624:       /* CANNOT use MPI_Allreduce() since it might fail the line number check */
1625:       if (localStageUsed[stage]) {
1626:         MPI_Allreduce(&stageInfo[stage].perfInfo.time,          &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1627:         MPI_Allreduce(&stageInfo[stage].perfInfo.flops,         &flops,     1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1628:         MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages,   &mess,      1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1629:         MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1630:         MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red,       1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1631:         name = stageInfo[stage].name;
1632:       } else {
1633:         MPI_Allreduce(&zero,                           &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1634:         MPI_Allreduce(&zero,                           &flops,     1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1635:         MPI_Allreduce(&zero,                           &mess,      1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1636:         MPI_Allreduce(&zero,                           &messLen,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1637:         MPI_Allreduce(&zero,                           &red,       1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1638:         name = "";
1639:       }
1640:       mess *= 0.5; messLen *= 0.5; red /= size;
1641:       if (TotalTime     != 0.0) fracTime       = stageTime/TotalTime;    else fracTime       = 0.0;
1642:       if (TotalFlops    != 0.0) fracFlops      = flops/TotalFlops;       else fracFlops      = 0.0;
1643:       /* Talk to Barry if (stageTime     != 0.0) flops          = (size*flops)/stageTime; else flops          = 0.0; */
1644:       if (numMessages   != 0.0) fracMessages   = mess/numMessages;       else fracMessages   = 0.0;
1645:       if (mess          != 0.0) avgMessLen     = messLen/mess;           else avgMessLen     = 0.0;
1646:       if (messageLength != 0.0) fracLength     = messLen/messageLength;  else fracLength     = 0.0;
1647:       if (numReductions != 0.0) fracReductions = red/numReductions;      else fracReductions = 0.0;
1648:       PetscFPrintf(comm, fd, "%2d: %15s: %6.4e %5.1f%%  %6.4e %5.1f%%  %5.3e %5.1f%%  %5.3e      %5.1f%%  %5.3e %5.1f%%\n",
1649:                           stage, name, stageTime/size, 100.0*fracTime, flops, 100.0*fracFlops,
1650:                           mess, 100.0*fracMessages, avgMessLen, 100.0*fracLength, red, 100.0*fracReductions);
1651:     }
1652:   }

1654:   PetscFPrintf(comm, fd,"\n------------------------------------------------------------------------------------------------------------------------\n");
1655:   PetscFPrintf(comm, fd, "See the 'Profiling' chapter of the users' manual for details on interpreting output.\n");
1656:   PetscFPrintf(comm, fd, "Phase summary info:\n");
1657:   PetscFPrintf(comm, fd, "   Count: number of times phase was executed\n");
1658:   PetscFPrintf(comm, fd, "   Time and Flop: Max - maximum over all processors\n");
1659:   PetscFPrintf(comm, fd, "                  Ratio - ratio of maximum to minimum over all processors\n");
1660:   PetscFPrintf(comm, fd, "   Mess: number of messages sent\n");
1661:   PetscFPrintf(comm, fd, "   AvgLen: average message length (bytes)\n");
1662:   PetscFPrintf(comm, fd, "   Reduct: number of global reductions\n");
1663:   PetscFPrintf(comm, fd, "   Global: entire computation\n");
1664:   PetscFPrintf(comm, fd, "   Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().\n");
1665:   PetscFPrintf(comm, fd, "      %%T - percent time in this phase         %%F - percent flop in this phase\n");
1666:   PetscFPrintf(comm, fd, "      %%M - percent messages in this phase     %%L - percent message lengths in this phase\n");
1667:   PetscFPrintf(comm, fd, "      %%R - percent reductions in this phase\n");
1668:   PetscFPrintf(comm, fd, "   Total Mflop/s: 10e-6 * (sum of flop over all processors)/(max time over all processors)\n");
1669:   if (PetscLogMemory) {
1670:     PetscFPrintf(comm, fd, "   Malloc Mbytes: Memory allocated and kept during event (sum over all calls to event)\n");
1671:     PetscFPrintf(comm, fd, "   EMalloc Mbytes: extra memory allocated during event and then freed (maximum over all calls to events)\n");
1672:     PetscFPrintf(comm, fd, "   MMalloc Mbytes: Increase in high water mark of allocated memory (sum over all calls to event)\n");
1673:     PetscFPrintf(comm, fd, "   RMI Mbytes: Increase in resident memory (sum over all calls to event)\n");
1674:   }
1675:   #if defined(PETSC_HAVE_DEVICE)
1676:   PetscFPrintf(comm, fd, "   GPU Mflop/s: 10e-6 * (sum of flop on GPU over all processors)/(max GPU time over all processors)\n");
1677:   PetscFPrintf(comm, fd, "   CpuToGpu Count: total number of CPU to GPU copies per processor\n");
1678:   PetscFPrintf(comm, fd, "   CpuToGpu Size (Mbytes): 10e-6 * (total size of CPU to GPU copies per processor)\n");
1679:   PetscFPrintf(comm, fd, "   GpuToCpu Count: total number of GPU to CPU copies per processor\n");
1680:   PetscFPrintf(comm, fd, "   GpuToCpu Size (Mbytes): 10e-6 * (total size of GPU to CPU copies per processor)\n");
1681:   PetscFPrintf(comm, fd, "   GPU %%F: percent flops on GPU in this event\n");
1682:   #endif
1683:   PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------\n");

1685:   PetscLogViewWarnDebugging(comm,fd);

1687:   /* Report events */
1688:   PetscFPrintf(comm, fd,"Event                Count      Time (sec)     Flop                              --- Global ---  --- Stage ----  Total");
1689:   if (PetscLogMemory) {
1690:     PetscFPrintf(comm, fd,"  Malloc EMalloc MMalloc RMI");
1691:   }
1692:   #if defined(PETSC_HAVE_DEVICE)
1693:   PetscFPrintf(comm, fd,"   GPU    - CpuToGpu -   - GpuToCpu - GPU");
1694:   #endif
1695:   PetscFPrintf(comm, fd,"\n");
1696:   PetscFPrintf(comm, fd,"                   Max Ratio  Max     Ratio   Max  Ratio  Mess   AvgLen  Reduct  %%T %%F %%M %%L %%R  %%T %%F %%M %%L %%R Mflop/s");
1697:   if (PetscLogMemory) {
1698:     PetscFPrintf(comm, fd," Mbytes Mbytes Mbytes Mbytes");
1699:   }
1700:   #if defined(PETSC_HAVE_DEVICE)
1701:   PetscFPrintf(comm, fd," Mflop/s Count   Size   Count   Size  %%F");
1702:   #endif
1703:   PetscFPrintf(comm, fd,"\n");
1704:   PetscFPrintf(comm, fd,"------------------------------------------------------------------------------------------------------------------------");
1705:   if (PetscLogMemory) {
1706:     PetscFPrintf(comm, fd,"-----------------------------");
1707:   }
1708:   #if defined(PETSC_HAVE_DEVICE)
1709:   PetscFPrintf(comm, fd,"---------------------------------------");
1710:   #endif
1711:   PetscFPrintf(comm, fd,"\n");

1713:   /* Problem: The stage name will not show up unless the stage executed on proc 1 */
1714:   for (stage = 0; stage < numStages; stage++) {
1715:     if (!stageVisible[stage]) continue;
1716:     /* CANNOT use MPI_Allreduce() since it might fail the line number check */
1717:     if (localStageUsed[stage]) {
1718:       PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name);
1719:       MPI_Allreduce(&stageInfo[stage].perfInfo.time,          &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1720:       MPI_Allreduce(&stageInfo[stage].perfInfo.flops,         &flops,     1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1721:       MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages,   &mess,      1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1722:       MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1723:       MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red,       1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1724:     } else {
1725:       PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage);
1726:       MPI_Allreduce(&zero,                           &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1727:       MPI_Allreduce(&zero,                           &flops,     1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1728:       MPI_Allreduce(&zero,                           &mess,      1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1729:       MPI_Allreduce(&zero,                           &messLen,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1730:       MPI_Allreduce(&zero,                           &red,       1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1731:     }
1732:     mess *= 0.5; messLen *= 0.5; red /= size;

1734:     /* Get total number of events in this stage --
1735:        Currently, a single processor can register more events than another, but events must all be registered in order,
1736:        just like stages. We can removed this requirement if necessary by having a global event numbering and indirection
1737:        on the event ID. This seems best accomplished by associating a communicator with each stage.

1739:        Problem: If the event did not happen on proc 1, its name will not be available.
1740:        Problem: Event visibility is not implemented
1741:     */
1742:     if (localStageUsed[stage]) {
1743:       eventInfo      = stageLog->stageInfo[stage].eventLog->eventInfo;
1744:       localNumEvents = stageLog->stageInfo[stage].eventLog->numEvents;
1745:     } else localNumEvents = 0;
1746:     MPI_Allreduce(&localNumEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);
1747:     for (event = 0; event < numEvents; event++) {
1748:       /* CANNOT use MPI_Allreduce() since it might fail the line number check */
1749:       if (localStageUsed[stage] && (event < stageLog->stageInfo[stage].eventLog->numEvents) && (eventInfo[event].depth == 0)) {
1750:         if ((eventInfo[event].count > 0) && (eventInfo[event].time > 0.0)) flopr = eventInfo[event].flops; else flopr = 0.0;
1751:         MPI_Allreduce(&flopr,                          &minf,  1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1752:         MPI_Allreduce(&flopr,                          &maxf,  1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1753:         MPI_Allreduce(&eventInfo[event].flops,         &totf,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1754:         MPI_Allreduce(&eventInfo[event].time,          &mint,  1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1755:         MPI_Allreduce(&eventInfo[event].time,          &maxt,  1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1756:         MPI_Allreduce(&eventInfo[event].time,          &tott,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1757:         MPI_Allreduce(&eventInfo[event].numMessages,   &totm,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1758:         MPI_Allreduce(&eventInfo[event].messageLength, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1759:         MPI_Allreduce(&eventInfo[event].numReductions, &totr,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1760:         MPI_Allreduce(&eventInfo[event].count,         &minC,  1, MPI_INT,             MPI_MIN, comm);
1761:         MPI_Allreduce(&eventInfo[event].count,         &maxC,  1, MPI_INT,             MPI_MAX, comm);
1762:         if (PetscLogMemory) {
1763:           MPI_Allreduce(&eventInfo[event].memIncrease,    &mem,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1764:           MPI_Allreduce(&eventInfo[event].mallocSpace,    &mal,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1765:           MPI_Allreduce(&eventInfo[event].mallocIncrease, &malmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1766:           MPI_Allreduce(&eventInfo[event].mallocIncreaseEvent, &emalmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1767:         }
1768:         #if defined(PETSC_HAVE_DEVICE)
1769:         MPI_Allreduce(&eventInfo[event].CpuToGpuCount,    &cct,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1770:         MPI_Allreduce(&eventInfo[event].GpuToCpuCount,    &gct,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1771:         MPI_Allreduce(&eventInfo[event].CpuToGpuSize,     &csz,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1772:         MPI_Allreduce(&eventInfo[event].GpuToCpuSize,     &gsz,   1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1773:         MPI_Allreduce(&eventInfo[event].GpuFlops,         &gflops,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1774:         MPI_Allreduce(&eventInfo[event].GpuTime,          &gmaxt ,1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1775:         #endif
1776:         name = stageLog->eventLog->eventInfo[event].name;
1777:       } else {
1778:         flopr = 0.0;
1779:         MPI_Allreduce(&flopr,                         &minf,  1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1780:         MPI_Allreduce(&flopr,                         &maxf,  1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1781:         MPI_Allreduce(&zero,                          &totf,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1782:         MPI_Allreduce(&zero,                          &mint,  1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);
1783:         MPI_Allreduce(&zero,                          &maxt,  1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1784:         MPI_Allreduce(&zero,                          &tott,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1785:         MPI_Allreduce(&zero,                          &totm,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1786:         MPI_Allreduce(&zero,                          &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1787:         MPI_Allreduce(&zero,                          &totr,  1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1788:         MPI_Allreduce(&ierr,                          &minC,  1, MPI_INT,             MPI_MIN, comm);
1789:         MPI_Allreduce(&ierr,                          &maxC,  1, MPI_INT,             MPI_MAX, comm);
1790:         if (PetscLogMemory) {
1791:           MPI_Allreduce(&zero,                        &mem,    1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1792:           MPI_Allreduce(&zero,                        &mal,    1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1793:           MPI_Allreduce(&zero,                        &malmax, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1794:           MPI_Allreduce(&zero,                        &emalmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1795:         }
1796:         #if defined(PETSC_HAVE_DEVICE)
1797:         MPI_Allreduce(&zero,                          &cct,    1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1798:         MPI_Allreduce(&zero,                          &gct,    1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1799:         MPI_Allreduce(&zero,                          &csz,    1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1800:         MPI_Allreduce(&zero,                          &gsz,    1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1801:         MPI_Allreduce(&zero,                          &gflops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);
1802:         MPI_Allreduce(&zero,                          &gmaxt , 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);
1803:         #endif
1804:         name  = "";
1805:       }
1806:       if (mint < 0.0) {
1807:         PetscFPrintf(comm, fd, "WARNING!!! Minimum time %g over all processors for %s is negative! This happens\n on some machines whose times cannot handle too rapid calls.!\n artificially changing minimum to zero.\n",mint,name);
1808:         mint = 0;
1809:       }
1810:       if (minf < 0.0) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Minimum flop %g over all processors for %s is negative! Not possible!",minf,name);
1811:       totm *= 0.5; totml *= 0.5; totr /= size;

1813:       if (maxC != 0) {
1814:         if (minC          != 0)   ratC             = ((PetscLogDouble)maxC)/minC;else ratC             = 0.0;
1815:         if (mint          != 0.0) ratt             = maxt/mint;                  else ratt             = 0.0;
1816:         if (minf          != 0.0) ratf             = maxf/minf;                  else ratf             = 0.0;
1817:         if (TotalTime     != 0.0) fracTime         = tott/TotalTime;             else fracTime         = 0.0;
1818:         if (TotalFlops    != 0.0) fracFlops        = totf/TotalFlops;            else fracFlops        = 0.0;
1819:         if (stageTime     != 0.0) fracStageTime    = tott/stageTime;             else fracStageTime    = 0.0;
1820:         if (flops         != 0.0) fracStageFlops   = totf/flops;                 else fracStageFlops   = 0.0;
1821:         if (numMessages   != 0.0) fracMess         = totm/numMessages;           else fracMess         = 0.0;
1822:         if (messageLength != 0.0) fracMessLen      = totml/messageLength;        else fracMessLen      = 0.0;
1823:         if (numReductions != 0.0) fracRed          = totr/numReductions;         else fracRed          = 0.0;
1824:         if (mess          != 0.0) fracStageMess    = totm/mess;                  else fracStageMess    = 0.0;
1825:         if (messLen       != 0.0) fracStageMessLen = totml/messLen;              else fracStageMessLen = 0.0;
1826:         if (red           != 0.0) fracStageRed     = totr/red;                   else fracStageRed     = 0.0;
1827:         if (totm          != 0.0) totml           /= totm;                       else totml            = 0.0;
1828:         if (maxt          != 0.0) flopr            = totf/maxt;                  else flopr            = 0.0;
1829:         if (fracStageTime > 1.00) {PetscFPrintf(comm, fd,"Warning -- total time of event greater than time of entire stage -- something is wrong with the timer\n");}
1830:         PetscFPrintf(comm, fd,
1831:                             "%-16s %7d%4.1f %5.4e%4.1f %3.2e%4.1f %2.1e %2.1e %2.1e%3.0f%3.0f%3.0f%3.0f%3.0f %3.0f%3.0f%3.0f%3.0f%3.0f %5.0f",
1832:                             name, maxC, ratC, maxt, ratt, maxf, ratf, totm, totml, totr,
1833:                             100.0*fracTime, 100.0*fracFlops, 100.0*fracMess, 100.0*fracMessLen, 100.0*fracRed,
1834:                             100.0*fracStageTime, 100.0*fracStageFlops, 100.0*fracStageMess, 100.0*fracStageMessLen, 100.0*fracStageRed,
1835:                             PetscAbs(flopr)/1.0e6);
1836:         if (PetscLogMemory) {
1837:           PetscFPrintf(comm, fd," %5.0f   %5.0f   %5.0f   %5.0f",mal/1.0e6,emalmax/1.0e6,malmax/1.0e6,mem/1.0e6);
1838:         }
1839:         #if defined(PETSC_HAVE_DEVICE)
1840:         if (totf  != 0.0) fracgflops = gflops/totf;  else fracgflops = 0.0;
1841:         if (gmaxt != 0.0) gflopr     = gflops/gmaxt; else gflopr     = 0.0;
1842:         PetscFPrintf(comm, fd,"   %5.0f   %4.0f %3.2e %4.0f %3.2e% 3.0f",PetscAbs(gflopr)/1.0e6,cct/size,csz/(1.0e6*size),gct/size,gsz/(1.0e6*size),100.0*fracgflops);
1843:         #endif
1844:         PetscFPrintf(comm, fd,"\n");
1845:       }
1846:     }
1847:   }

1849:   /* Memory usage and object creation */
1850:   PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------");
1851:   if (PetscLogMemory) {
1852:     PetscFPrintf(comm, fd, "-----------------------------");
1853:   }
1854:   #if defined(PETSC_HAVE_DEVICE)
1855:   PetscFPrintf(comm, fd, "---------------------------------------");
1856:   #endif
1857:   PetscFPrintf(comm, fd, "\n");
1858:   PetscFPrintf(comm, fd, "\n");
1859:   PetscFPrintf(comm, fd, "Memory usage is given in bytes:\n\n");

1861:   /* Right now, only stages on the first processor are reported here, meaning only objects associated with
1862:      the global communicator, or MPI_COMM_SELF for proc 1. We really should report global stats and then
1863:      stats for stages local to processor sets.
1864:   */
1865:   /* We should figure out the longest object name here (now 20 characters) */
1866:   PetscFPrintf(comm, fd, "Object Type          Creations   Destructions     Memory  Descendants' Mem.\n");
1867:   PetscFPrintf(comm, fd, "Reports information only for process 0.\n");
1868:   for (stage = 0; stage < numStages; stage++) {
1869:     if (localStageUsed[stage]) {
1870:       classInfo = stageLog->stageInfo[stage].classLog->classInfo;
1871:       PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name);
1872:       for (oclass = 0; oclass < stageLog->stageInfo[stage].classLog->numClasses; oclass++) {
1873:         if ((classInfo[oclass].creations > 0) || (classInfo[oclass].destructions > 0)) {
1874:           PetscFPrintf(comm, fd, "%20s %5d          %5d  %11.0f     %g\n", stageLog->classLog->classInfo[oclass].name,
1875:                               classInfo[oclass].creations, classInfo[oclass].destructions, classInfo[oclass].mem,
1876:                               classInfo[oclass].descMem);
1877:         }
1878:       }
1879:     } else {
1880:       if (!localStageVisible[stage]) continue;
1881:       PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage);
1882:     }
1883:   }

1885:   PetscFree(localStageUsed);
1886:   PetscFree(stageUsed);
1887:   PetscFree(localStageVisible);
1888:   PetscFree(stageVisible);

1890:   /* Information unrelated to this particular run */
1891:   PetscFPrintf(comm, fd, "========================================================================================================================\n");
1892:   PetscTime(&y);
1893:   PetscTime(&x);
1894:   PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y);
1895:   PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y);
1896:   PetscFPrintf(comm,fd,"Average time to get PetscTime(): %g\n", (y-x)/10.0);
1897:   /* MPI information */
1898:   if (size > 1) {
1899:     MPI_Status  status;
1900:     PetscMPIInt tag;
1901:     MPI_Comm    newcomm;

1903:     MPI_Barrier(comm);
1904:     PetscTime(&x);
1905:     MPI_Barrier(comm);
1906:     MPI_Barrier(comm);
1907:     MPI_Barrier(comm);
1908:     MPI_Barrier(comm);
1909:     MPI_Barrier(comm);
1910:     PetscTime(&y);
1911:     PetscFPrintf(comm, fd, "Average time for MPI_Barrier(): %g\n", (y-x)/5.0);
1912:     PetscCommDuplicate(comm,&newcomm, &tag);
1913:     MPI_Barrier(comm);
1914:     if (rank) {
1915:       MPI_Recv(NULL, 0, MPI_INT, rank-1,            tag, newcomm, &status);
1916:       MPI_Send(NULL, 0, MPI_INT, (rank+1)%size, tag, newcomm);
1917:     } else {
1918:       PetscTime(&x);
1919:       MPI_Send(NULL, 0, MPI_INT, 1,          tag, newcomm);
1920:       MPI_Recv(NULL, 0, MPI_INT, size-1, tag, newcomm, &status);
1921:       PetscTime(&y);
1922:       PetscFPrintf(comm,fd,"Average time for zero size MPI_Send(): %g\n", (y-x)/size);
1923:     }
1924:     PetscCommDestroy(&newcomm);
1925:   }
1926:   PetscOptionsView(NULL,viewer);

1928:   /* Machine and compile information */
1929: #if defined(PETSC_USE_FORTRAN_KERNELS)
1930:   PetscFPrintf(comm, fd, "Compiled with FORTRAN kernels\n");
1931: #else
1932:   PetscFPrintf(comm, fd, "Compiled without FORTRAN kernels\n");
1933: #endif
1934: #if defined(PETSC_USE_64BIT_INDICES)
1935:   PetscFPrintf(comm, fd, "Compiled with 64 bit PetscInt\n");
1936: #elif defined(PETSC_USE___FLOAT128)
1937:   PetscFPrintf(comm, fd, "Compiled with 32 bit PetscInt\n");
1938: #endif
1939: #if defined(PETSC_USE_REAL_SINGLE)
1940:   PetscFPrintf(comm, fd, "Compiled with single precision PetscScalar and PetscReal\n");
1941: #elif defined(PETSC_USE___FLOAT128)
1942:   PetscFPrintf(comm, fd, "Compiled with 128 bit precision PetscScalar and PetscReal\n");
1943: #endif
1944: #if defined(PETSC_USE_REAL_MAT_SINGLE)
1945:   PetscFPrintf(comm, fd, "Compiled with single precision matrices\n");
1946: #else
1947:   PetscFPrintf(comm, fd, "Compiled with full precision matrices (default)\n");
1948: #endif
1949:   PetscFPrintf(comm, fd, "sizeof(short) %d sizeof(int) %d sizeof(long) %d sizeof(void*) %d sizeof(PetscScalar) %d sizeof(PetscInt) %d\n",
1950:                       (int) sizeof(short), (int) sizeof(int), (int) sizeof(long), (int) sizeof(void*),(int) sizeof(PetscScalar),(int) sizeof(PetscInt));

1952:   PetscFPrintf(comm, fd, "Configure options: %s",petscconfigureoptions);
1953:   PetscFPrintf(comm, fd, "%s", petscmachineinfo);
1954:   PetscFPrintf(comm, fd, "%s", petsccompilerinfo);
1955:   PetscFPrintf(comm, fd, "%s", petsccompilerflagsinfo);
1956:   PetscFPrintf(comm, fd, "%s", petsclinkerinfo);

1958:   /* Cleanup */
1959:   PetscFPrintf(comm, fd, "\n");
1960:   PetscLogViewWarnNoGpuAwareMpi(comm,fd);
1961:   PetscLogViewWarnDebugging(comm,fd);
1962:   return(0);
1963: }

1965: /*@C
1966:   PetscLogView - Prints a summary of the logging.

1968:   Collective over MPI_Comm

1970:   Input Parameter:
1971: .  viewer - an ASCII viewer

1973:   Options Database Keys:
1974: +  -log_view [:filename] - Prints summary of log information
1975: .  -log_view :filename.py:ascii_info_detail - Saves logging information from each process as a Python file
1976: .  -log_view :filename.xml:ascii_xml - Saves a summary of the logging information in a nested format (see below for how to view it)
1977: .  -log_view :filename.txt:ascii_flamegraph - Saves logging information in a format suitable for visualising as a Flame Graph (see below for how to view it)
1978: .  -log_all - Saves a file Log.rank for each MPI process with details of each step of the computation
1979: -  -log_trace [filename] - Displays a trace of what each process is doing

1981:   Notes:
1982:   It is possible to control the logging programatically but we recommend using the options database approach whenever possible
1983:   By default the summary is printed to stdout.

1985:   Before calling this routine you must have called either PetscLogDefaultBegin() or PetscLogNestedBegin()

1987:   If PETSc is configured with --with-logging=0 then this functionality is not available

1989:   To view the nested XML format filename.xml first copy  ${PETSC_DIR}/share/petsc/xml/performance_xml2html.xsl to the current
1990:   directory then open filename.xml with your browser. Specific notes for certain browsers
1991: $    Firefox and Internet explorer - simply open the file
1992: $    Google Chrome - you must start up Chrome with the option --allow-file-access-from-files
1993: $    Safari - see https://ccm.net/faq/36342-safari-how-to-enable-local-file-access
1994:   or one can use the package http://xmlsoft.org/XSLT/xsltproc2.html to translate the xml file to html and then open it with
1995:   your browser.
1996:   Alternatively, use the script ${PETSC_DIR}/lib/petsc/bin/petsc-performance-view to automatically open a new browser
1997:   window and render the XML log file contents.

1999:   The nested XML format was kindly donated by Koos Huijssen and Christiaan M. Klaij  MARITIME  RESEARCH  INSTITUTE  NETHERLANDS

2001:   The Flame Graph output can be visualised using either the original Flame Graph script (https://github.com/brendangregg/FlameGraph)
2002:   or using speedscope (https://www.speedscope.app).
2003:   Old XML profiles may be converted into this format using the script ${PETSC_DIR}/lib/petsc/bin/xml2flamegraph.py.

2005:   Level: beginner

2007: .seealso: PetscLogDefaultBegin(), PetscLogDump()
2008: @*/
2009: PetscErrorCode  PetscLogView(PetscViewer viewer)
2010: {
2011:   PetscErrorCode    ierr;
2012:   PetscBool         isascii;
2013:   PetscViewerFormat format;
2014:   int               stage, lastStage;
2015:   PetscStageLog     stageLog;

2018:   if (!PetscLogPLB) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Must use -log_view or PetscLogDefaultBegin() before calling this routine");
2019:   /* Pop off any stages the user forgot to remove */
2020:   lastStage = 0;
2021:   PetscLogGetStageLog(&stageLog);
2022:   PetscStageLogGetCurrent(stageLog, &stage);
2023:   while (stage >= 0) {
2024:     lastStage = stage;
2025:     PetscStageLogPop(stageLog);
2026:     PetscStageLogGetCurrent(stageLog, &stage);
2027:   }
2028:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
2029:   if (!isascii) SETERRQ(PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Currently can only view logging to ASCII");
2030:   PetscViewerGetFormat(viewer,&format);
2031:   if (format == PETSC_VIEWER_DEFAULT || format == PETSC_VIEWER_ASCII_INFO) {
2032:     PetscLogView_Default(viewer);
2033:   } else if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
2034:     PetscLogView_Detailed(viewer);
2035:   } else if (format == PETSC_VIEWER_ASCII_CSV) {
2036:     PetscLogView_CSV(viewer);
2037:   } else if (format == PETSC_VIEWER_ASCII_XML) {
2038:     PetscLogView_Nested(viewer);
2039:   } else if (format == PETSC_VIEWER_ASCII_FLAMEGRAPH) {
2040:     PetscLogView_Flamegraph(viewer);
2041:   }
2042:   PetscStageLogPush(stageLog, lastStage);
2043:   return(0);
2044: }

2046: /*@C
2047:   PetscLogViewFromOptions - Processes command line options to determine if/how a PetscLog is to be viewed.

2049:   Collective on PETSC_COMM_WORLD

2051:   Not normally called by user

2053:   Level: intermediate

2055: @*/
2056: PetscErrorCode PetscLogViewFromOptions(void)
2057: {
2058:   PetscErrorCode    ierr;
2059:   PetscViewer       viewer;
2060:   PetscBool         flg;
2061:   PetscViewerFormat format;

2064:   PetscOptionsGetViewer(PETSC_COMM_WORLD,NULL,NULL,"-log_view",&viewer,&format,&flg);
2065:   if (flg) {
2066:     PetscViewerPushFormat(viewer,format);
2067:     PetscLogView(viewer);
2068:     PetscViewerPopFormat(viewer);
2069:     PetscViewerDestroy(&viewer);
2070:   }
2071:   return(0);
2072: }

2074: /*----------------------------------------------- Counter Functions -------------------------------------------------*/
2075: /*@C
2076:    PetscGetFlops - Returns the number of flops used on this processor
2077:    since the program began.

2079:    Not Collective

2081:    Output Parameter:
2082:    flops - number of floating point operations

2084:    Notes:
2085:    A global counter logs all PETSc flop counts.  The user can use
2086:    PetscLogFlops() to increment this counter to include flops for the
2087:    application code.

2089:    Level: intermediate

2091: .seealso: PetscTime(), PetscLogFlops()
2092: @*/
2093: PetscErrorCode  PetscGetFlops(PetscLogDouble *flops)
2094: {
2096:   *flops = petsc_TotalFlops;
2097:   return(0);
2098: }

2100: PetscErrorCode  PetscLogObjectState(PetscObject obj, const char format[], ...)
2101: {
2103:   size_t         fullLength;
2104:   va_list        Argp;

2107:   if (!petsc_logObjects) return(0);
2108:   va_start(Argp, format);
2109:   PetscVSNPrintf(petsc_objects[obj->id].info, 64,format,&fullLength, Argp);
2110:   va_end(Argp);
2111:   return(0);
2112: }

2114: /*MC
2115:    PetscLogFlops - Adds floating point operations to the global counter.

2117:    Synopsis:
2118: #include <petsclog.h>
2119:    PetscErrorCode PetscLogFlops(PetscLogDouble f)

2121:    Not Collective

2123:    Input Parameter:
2124: .  f - flop counter

2126:    Usage:
2127: .vb
2128:      PetscLogEvent USER_EVENT;
2129:      PetscLogEventRegister("User event",0,&USER_EVENT);
2130:      PetscLogEventBegin(USER_EVENT,0,0,0,0);
2131:         [code segment to monitor]
2132:         PetscLogFlops(user_flops)
2133:      PetscLogEventEnd(USER_EVENT,0,0,0,0);
2134: .ve

2136:    Notes:
2137:    A global counter logs all PETSc flop counts.  The user can use
2138:    PetscLogFlops() to increment this counter to include flops for the
2139:    application code.

2141:    Level: intermediate

2143: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscGetFlops()

2145: M*/

2147: /*MC
2148:    PetscPreLoadBegin - Begin a segment of code that may be preloaded (run twice)
2149:     to get accurate timings

2151:    Synopsis:
2152: #include <petsclog.h>
2153:    void PetscPreLoadBegin(PetscBool  flag,char *name);

2155:    Not Collective

2157:    Input Parameters:
2158: +   flag - PETSC_TRUE to run twice, PETSC_FALSE to run once, may be overridden
2159:            with command line option -preload true or -preload false
2160: -   name - name of first stage (lines of code timed separately with -log_view) to
2161:            be preloaded

2163:    Usage:
2164: .vb
2165:      PetscPreLoadBegin(PETSC_TRUE,"first stage);
2166:        lines of code
2167:        PetscPreLoadStage("second stage");
2168:        lines of code
2169:      PetscPreLoadEnd();
2170: .ve

2172:    Notes:
2173:     Only works in C/C++, not Fortran

2175:      Flags available within the macro.
2176: +    PetscPreLoadingUsed - true if we are or have done preloading
2177: .    PetscPreLoadingOn - true if it is CURRENTLY doing preload
2178: .    PetscPreLoadIt - 0 for the first computation (with preloading turned off it is only 0) 1 for the second
2179: -    PetscPreLoadMax - number of times it will do the computation, only one when preloading is turned on
2180:      The first two variables are available throughout the program, the second two only between the PetscPreLoadBegin()
2181:      and PetscPreLoadEnd()

2183:    Level: intermediate

2185: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadEnd(), PetscPreLoadStage()

2187: M*/

2189: /*MC
2190:    PetscPreLoadEnd - End a segment of code that may be preloaded (run twice)
2191:     to get accurate timings

2193:    Synopsis:
2194: #include <petsclog.h>
2195:    void PetscPreLoadEnd(void);

2197:    Not Collective

2199:    Usage:
2200: .vb
2201:      PetscPreLoadBegin(PETSC_TRUE,"first stage);
2202:        lines of code
2203:        PetscPreLoadStage("second stage");
2204:        lines of code
2205:      PetscPreLoadEnd();
2206: .ve

2208:    Notes:
2209:     only works in C/C++ not fortran

2211:    Level: intermediate

2213: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadStage()

2215: M*/

2217: /*MC
2218:    PetscPreLoadStage - Start a new segment of code to be timed separately.
2219:     to get accurate timings

2221:    Synopsis:
2222: #include <petsclog.h>
2223:    void PetscPreLoadStage(char *name);

2225:    Not Collective

2227:    Usage:
2228: .vb
2229:      PetscPreLoadBegin(PETSC_TRUE,"first stage);
2230:        lines of code
2231:        PetscPreLoadStage("second stage");
2232:        lines of code
2233:      PetscPreLoadEnd();
2234: .ve

2236:    Notes:
2237:     only works in C/C++ not fortran

2239:    Level: intermediate

2241: .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd()

2243: M*/

2245: #if defined(PETSC_HAVE_DEVICE)

2247: #if defined(PETSC_HAVE_CUDA)
2248: #include <petscdevice.h>
2249: PETSC_EXTERN cudaEvent_t petsc_gputimer_begin;
2250: PETSC_EXTERN cudaEvent_t petsc_gputimer_end;
2251: #endif

2253: #if defined(PETSC_HAVE_HIP)
2254: #include <petscdevice.h>
2255: PETSC_EXTERN hipEvent_t petsc_gputimer_begin;
2256: PETSC_EXTERN hipEvent_t petsc_gputimer_end;
2257: #endif

2259: /*-------------------------------------------- GPU event Functions ----------------------------------------------*/
2260: /*@C
2261:   PetscLogGpuTimeBegin - Start timer for device

2263:   Notes:
2264:     When CUDA or HIP is enabled, the timer is run on the GPU, it is a separate logging of time devoted to GPU computations (excluding kernel launch times).
2265:     When CUDA or HIP is not available, the timer is run on the CPU, it is a separate logging of time devoted to GPU computations (including kernel launch times).
2266:     There is no need to call WaitForCUDA() or WaitForHIP() between PetscLogGpuTimeBegin and PetscLogGpuTimeEnd
2267:     This timer should NOT include times for data transfers between the GPU and CPU, nor setup actions such as allocating space.
2268:     The regular logging captures the time for data transfers and any CPU activites during the event
2269:     It is used to compute the flop rate on the GPU as it is actively engaged in running a kernel.

2271:   Developer Notes:
2272:     The GPU event timer captures the execution time of all the kernels launched in the default stream by the CPU between PetscLogGpuTimeBegin() and PetsLogGpuTimeEnd().
2273:     PetscLogGpuTimeBegin() and PetsLogGpuTimeEnd() insert the begin and end events into the default stream (stream 0). The device will record a time stamp for the event when it reaches that event in the stream. The function xxxEventSynchronize() is called in PetsLogGpuTimeEnd() to block CPU execution, but not continued GPU excution, until the timer event is recorded.

2275:   Level: intermediate

2277: .seealso:  PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeEnd()
2278: @*/
2279: PetscErrorCode PetscLogGpuTimeBegin(void)
2280: {
2281: #if defined(PETSC_HAVE_CUDA)
2282:   cudaError_t    cerr;
2283: #elif defined(PETSC_HAVE_HIP)
2284:   hipError_t     cerr;
2285: #else
2287: #endif
2289: #if defined(PETSC_USE_DEBUG)
2290:   if (petsc_gtime_inuse) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Forgot to call PetscLogGpuTimeEnd()?");
2291:   petsc_gtime_inuse = PETSC_TRUE;
2292: #endif
2293:   if (!PetscLogPLB) return(0);
2294: #if defined(PETSC_HAVE_CUDA)
2295:   cerr = cudaEventRecord(petsc_gputimer_begin,PetscDefaultCudaStream);CHKERRCUDA(cerr);
2296: #elif defined(PETSC_HAVE_HIP)
2297:   cerr = hipEventRecord(petsc_gputimer_begin,PetscDefaultHipStream);CHKERRHIP(cerr);
2298: #else
2299:   PetscTimeSubtract(&petsc_gtime);
2300: #endif
2301:   return(0);
2302: }

2304: /*@C
2305:   PetscLogGpuTimeEnd - Stop timer for device

2307:   Level: intermediate

2309: .seealso:  PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeBegin()
2310: @*/
2311: PetscErrorCode PetscLogGpuTimeEnd(void)
2312: {
2313: #if defined(PETSC_HAVE_CUDA)
2314:   float          gtime;
2315:   cudaError_t    cerr;
2316: #elif defined(PETSC_HAVE_HIP)
2317:   float          gtime;
2318:   hipError_t     cerr;
2319: #else
2321: #endif
2323: #if defined(PETSC_USE_DEBUG)
2324:   if (!petsc_gtime_inuse) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Forgot to call PetscLogGpuTimeBegin()?");
2325:   petsc_gtime_inuse = PETSC_FALSE;
2326: #endif
2327:   if (!PetscLogPLE) return(0);
2328: #if defined(PETSC_HAVE_CUDA)
2329:   cerr = cudaEventRecord(petsc_gputimer_end,PetscDefaultCudaStream);CHKERRCUDA(cerr);
2330:   cerr = cudaEventSynchronize(petsc_gputimer_end);CHKERRCUDA(cerr);
2331:   cerr = cudaEventElapsedTime(&gtime,petsc_gputimer_begin,petsc_gputimer_end);CHKERRCUDA(cerr);
2332:   petsc_gtime += (PetscLogDouble)gtime/1000.0; /* convert milliseconds to seconds */
2333: #elif defined(PETSC_HAVE_HIP)
2334:   cerr = hipEventRecord(petsc_gputimer_end,PetscDefaultHipStream);CHKERRHIP(cerr);
2335:   cerr = hipEventSynchronize(petsc_gputimer_end);CHKERRHIP(cerr);
2336:   cerr = hipEventElapsedTime(&gtime,petsc_gputimer_begin,petsc_gputimer_end);CHKERRHIP(cerr);
2337:   petsc_gtime += (PetscLogDouble)gtime/1000.0; /* convert milliseconds to seconds */
2338: #else
2339:   PetscTimeAdd(&petsc_gtime);
2340: #endif
2341:   return(0);
2342: }
2343: #endif /* end of PETSC_HAVE_DEVICE */

2345: #else /* end of -DPETSC_USE_LOG section */

2347: PetscErrorCode  PetscLogObjectState(PetscObject obj, const char format[], ...)
2348: {
2350:   return(0);
2351: }

2353: #endif /* PETSC_USE_LOG*/

2355: PetscClassId PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID;
2356: PetscClassId PETSC_OBJECT_CLASSID  = 0;

2358: /*@C
2359:   PetscClassIdRegister - Registers a new class name for objects and logging operations in an application code.

2361:   Not Collective

2363:   Input Parameter:
2364: . name   - The class name

2366:   Output Parameter:
2367: . oclass - The class id or classid

2369:   Level: developer

2371: @*/
2372: PetscErrorCode  PetscClassIdRegister(const char name[],PetscClassId *oclass)
2373: {
2374: #if defined(PETSC_USE_LOG)
2375:   PetscStageLog  stageLog;
2376:   PetscInt       stage;
2378: #endif

2381:   *oclass = ++PETSC_LARGEST_CLASSID;
2382: #if defined(PETSC_USE_LOG)
2383:   PetscLogGetStageLog(&stageLog);
2384:   PetscClassRegLogRegister(stageLog->classLog, name, *oclass);
2385:   for (stage = 0; stage < stageLog->numStages; stage++) {
2386:     PetscClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses);
2387:   }
2388: #endif
2389:   return(0);
2390: }

2392: #if defined(PETSC_USE_LOG) && defined(PETSC_HAVE_MPE)
2393: #include <mpe.h>

2395: PetscBool PetscBeganMPE = PETSC_FALSE;

2397: PETSC_INTERN PetscErrorCode PetscLogEventBeginMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject);
2398: PETSC_INTERN PetscErrorCode PetscLogEventEndMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject);

2400: /*@C
2401:    PetscLogMPEBegin - Turns on MPE logging of events. This creates large log files
2402:    and slows the program down.

2404:    Collective over PETSC_COMM_WORLD

2406:    Options Database Keys:
2407: . -log_mpe - Prints extensive log information

2409:    Notes:
2410:    A related routine is PetscLogDefaultBegin() (with the options key -log_view), which is
2411:    intended for production runs since it logs only flop rates and object
2412:    creation (and should not significantly slow the programs).

2414:    Level: advanced

2416: .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogEventActivate(),
2417:           PetscLogEventDeactivate()
2418: @*/
2419: PetscErrorCode  PetscLogMPEBegin(void)
2420: {

2424:   /* Do MPE initialization */
2425:   if (!MPE_Initialized_logging()) { /* This function exists in mpich 1.1.2 and higher */
2426:     PetscInfo(0,"Initializing MPE.\n");
2427:     MPE_Init_log();

2429:     PetscBeganMPE = PETSC_TRUE;
2430:   } else {
2431:     PetscInfo(0,"MPE already initialized. Not attempting to reinitialize.\n");
2432:   }
2433:   PetscLogSet(PetscLogEventBeginMPE, PetscLogEventEndMPE);
2434:   return(0);
2435: }

2437: /*@C
2438:    PetscLogMPEDump - Dumps the MPE logging info to file for later use with Jumpshot.

2440:    Collective over PETSC_COMM_WORLD

2442:    Level: advanced

2444: .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogMPEBegin()
2445: @*/
2446: PetscErrorCode  PetscLogMPEDump(const char sname[])
2447: {
2448:   char           name[PETSC_MAX_PATH_LEN];

2452:   if (PetscBeganMPE) {
2453:     PetscInfo(0,"Finalizing MPE.\n");
2454:     if (sname) {
2455:       PetscStrcpy(name,sname);
2456:     } else {
2457:       PetscGetProgramName(name,sizeof(name));
2458:     }
2459:     MPE_Finish_log(name);
2460:   } else {
2461:     PetscInfo(0,"Not finalizing MPE (not started by PETSc).\n");
2462:   }
2463:   return(0);
2464: }

2466: #define PETSC_RGB_COLORS_MAX 39
2467: static const char *PetscLogMPERGBColors[PETSC_RGB_COLORS_MAX] = {
2468:   "OliveDrab:      ",
2469:   "BlueViolet:     ",
2470:   "CadetBlue:      ",
2471:   "CornflowerBlue: ",
2472:   "DarkGoldenrod:  ",
2473:   "DarkGreen:      ",
2474:   "DarkKhaki:      ",
2475:   "DarkOliveGreen: ",
2476:   "DarkOrange:     ",
2477:   "DarkOrchid:     ",
2478:   "DarkSeaGreen:   ",
2479:   "DarkSlateGray:  ",
2480:   "DarkTurquoise:  ",
2481:   "DeepPink:       ",
2482:   "DarkKhaki:      ",
2483:   "DimGray:        ",
2484:   "DodgerBlue:     ",
2485:   "GreenYellow:    ",
2486:   "HotPink:        ",
2487:   "IndianRed:      ",
2488:   "LavenderBlush:  ",
2489:   "LawnGreen:      ",
2490:   "LemonChiffon:   ",
2491:   "LightCoral:     ",
2492:   "LightCyan:      ",
2493:   "LightPink:      ",
2494:   "LightSalmon:    ",
2495:   "LightSlateGray: ",
2496:   "LightYellow:    ",
2497:   "LimeGreen:      ",
2498:   "MediumPurple:   ",
2499:   "MediumSeaGreen: ",
2500:   "MediumSlateBlue:",
2501:   "MidnightBlue:   ",
2502:   "MintCream:      ",
2503:   "MistyRose:      ",
2504:   "NavajoWhite:    ",
2505:   "NavyBlue:       ",
2506:   "OliveDrab:      "
2507: };

2509: /*@C
2510:   PetscLogMPEGetRGBColor - This routine returns a rgb color useable with PetscLogEventRegister()

2512:   Not collective. Maybe it should be?

2514:   Output Parameter:
2515: . str - character string representing the color

2517:   Level: developer

2519: .seealso: PetscLogEventRegister
2520: @*/
2521: PetscErrorCode  PetscLogMPEGetRGBColor(const char *str[])
2522: {
2523:   static int idx = 0;

2526:   *str = PetscLogMPERGBColors[idx];
2527:   idx  = (idx + 1)% PETSC_RGB_COLORS_MAX;
2528:   return(0);
2529: }

2531: #endif /* PETSC_USE_LOG && PETSC_HAVE_MPE */