Actual source code: eventlog.c

petsc-master 2019-06-15
Report Typos and Errors

  2: /*
  3:      This defines part of the private API for logging performance information. It is intended to be used only by the
  4:    PETSc PetscLog...() interface and not elsewhere, nor by users. Hence the prototypes for these functions are NOT
  5:    in the public PETSc include files.

  7: */
  8:  #include <petsc/private/logimpl.h>

 10: PetscBool PetscLogSyncOn = PETSC_FALSE;
 11: PetscBool PetscLogMemory = PETSC_FALSE;
 12: #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) 
 13: PetscBool PetscLogGpuTraffic = PETSC_FALSE;
 14: #endif

 16: /*----------------------------------------------- Creation Functions -------------------------------------------------*/
 17: /* Note: these functions do not have prototypes in a public directory, so they are considered "internal" and not exported. */

 19: /*@C
 20:   PetscEventRegLogCreate - This creates a PetscEventRegLog object.

 22:   Not collective

 24:   Input Parameter:
 25: . eventLog - The PetscEventRegLog

 27:   Level: developer

 29: .seealso: PetscEventRegLogDestroy(), PetscStageLogCreate()
 30: @*/
 31: PetscErrorCode PetscEventRegLogCreate(PetscEventRegLog *eventLog)
 32: {
 33:   PetscEventRegLog l;
 34:   PetscErrorCode   ierr;

 37:   PetscNew(&l);
 38:   l->numEvents = 0;
 39:   l->maxEvents = 100;
 40:   PetscMalloc1(l->maxEvents,&l->eventInfo);
 41:   *eventLog    = l;
 42:   return(0);
 43: }

 45: /*@C
 46:   PetscEventRegLogDestroy - This destroys a PetscEventRegLog object.

 48:   Not collective

 50:   Input Paramter:
 51: . eventLog - The PetscEventRegLog

 53:   Level: developer

 55: .seealso: PetscEventRegLogCreate()
 56: @*/
 57: PetscErrorCode PetscEventRegLogDestroy(PetscEventRegLog eventLog)
 58: {
 59:   int            e;

 63:   for (e = 0; e < eventLog->numEvents; e++) {
 64:     PetscFree(eventLog->eventInfo[e].name);
 65:   }
 66:   PetscFree(eventLog->eventInfo);
 67:   PetscFree(eventLog);
 68:   return(0);
 69: }

 71: /*@C
 72:   PetscEventPerfLogCreate - This creates a PetscEventPerfLog object.

 74:   Not collective

 76:   Input Parameter:
 77: . eventLog - The PetscEventPerfLog

 79:   Level: developer

 81: .seealso: PetscEventPerfLogDestroy(), PetscStageLogCreate()
 82: @*/
 83: PetscErrorCode PetscEventPerfLogCreate(PetscEventPerfLog *eventLog)
 84: {
 85:   PetscEventPerfLog l;
 86:   PetscErrorCode    ierr;

 89:   PetscNew(&l);
 90:   l->numEvents = 0;
 91:   l->maxEvents = 100;
 92:   PetscMalloc1(l->maxEvents,&l->eventInfo);
 93:   *eventLog    = l;
 94:   return(0);
 95: }

 97: /*@C
 98:   PetscEventPerfLogDestroy - This destroys a PetscEventPerfLog object.

100:   Not collective

102:   Input Paramter:
103: . eventLog - The PetscEventPerfLog

105:   Level: developer

107: .seealso: PetscEventPerfLogCreate()
108: @*/
109: PetscErrorCode PetscEventPerfLogDestroy(PetscEventPerfLog eventLog)
110: {

114:   PetscFree(eventLog->eventInfo);
115:   PetscFree(eventLog);
116:   return(0);
117: }

119: /*------------------------------------------------ General Functions -------------------------------------------------*/
120: /*@C
121:   PetscEventPerfInfoClear - This clears a PetscEventPerfInfo object.

123:   Not collective

125:   Input Paramter:
126: . eventInfo - The PetscEventPerfInfo

128:   Level: developer

130: .seealso: PetscEventPerfLogCreate()
131: @*/
132: PetscErrorCode PetscEventPerfInfoClear(PetscEventPerfInfo *eventInfo)
133: {
135:   eventInfo->id            = -1;
136:   eventInfo->active        = PETSC_TRUE;
137:   eventInfo->visible       = PETSC_TRUE;
138:   eventInfo->depth         = 0;
139:   eventInfo->count         = 0;
140:   eventInfo->flops         = 0.0;
141:   eventInfo->flops2        = 0.0;
142:   eventInfo->flopsTmp      = 0.0;
143:   eventInfo->time          = 0.0;
144:   eventInfo->time2         = 0.0;
145:   eventInfo->timeTmp       = 0.0;
146:   eventInfo->syncTime      = 0.0;
147:   eventInfo->dof[0]        = -1.0;
148:   eventInfo->dof[1]        = -1.0;
149:   eventInfo->dof[2]        = -1.0;
150:   eventInfo->dof[3]        = -1.0;
151:   eventInfo->dof[4]        = -1.0;
152:   eventInfo->dof[5]        = -1.0;
153:   eventInfo->dof[6]        = -1.0;
154:   eventInfo->dof[7]        = -1.0;
155:   eventInfo->errors[0]     = -1.0;
156:   eventInfo->errors[1]     = -1.0;
157:   eventInfo->errors[2]     = -1.0;
158:   eventInfo->errors[3]     = -1.0;
159:   eventInfo->errors[4]     = -1.0;
160:   eventInfo->errors[5]     = -1.0;
161:   eventInfo->errors[6]     = -1.0;
162:   eventInfo->errors[7]     = -1.0;
163:   eventInfo->numMessages   = 0.0;
164:   eventInfo->messageLength = 0.0;
165:   eventInfo->numReductions = 0.0;
166:   #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) 
167:   eventInfo->CpuToGpuCount = 0.0;
168:   eventInfo->GpuToCpuCount = 0.0;
169:   eventInfo->CpuToGpuSize  = 0.0;
170:   eventInfo->GpuToCpuSize  = 0.0;
171:   #endif
172:   return(0);
173: }

175: /*@C
176:   PetscEventPerfInfoCopy - Copy the activity and visibility data in eventInfo to outInfo

178:   Not collective

180:   Input Paramter:
181: . eventInfo - The input PetscEventPerfInfo

183:   Output Paramter:
184: . outInfo   - The output PetscEventPerfInfo

186:   Level: developer

188: .seealso: PetscEventPerfInfoClear()
189: @*/
190: PetscErrorCode PetscEventPerfInfoCopy(PetscEventPerfInfo *eventInfo,PetscEventPerfInfo *outInfo)
191: {
193:   outInfo->id      = eventInfo->id;
194:   outInfo->active  = eventInfo->active;
195:   outInfo->visible = eventInfo->visible;
196:   return(0);
197: }

199: /*@C
200:   PetscEventPerfLogEnsureSize - This ensures that a PetscEventPerfLog is at least of a certain size.

202:   Not collective

204:   Input Paramters:
205: + eventLog - The PetscEventPerfLog
206: - size     - The size

208:   Level: developer

210: .seealso: PetscEventPerfLogCreate()
211: @*/
212: PetscErrorCode PetscEventPerfLogEnsureSize(PetscEventPerfLog eventLog,int size)
213: {
214:   PetscEventPerfInfo *eventInfo;
215:   PetscErrorCode     ierr;

218:   while (size > eventLog->maxEvents) {
219:     PetscMalloc1(eventLog->maxEvents*2,&eventInfo);
220:     PetscMemcpy(eventInfo,eventLog->eventInfo,eventLog->maxEvents * sizeof(PetscEventPerfInfo));
221:     PetscFree(eventLog->eventInfo);
222:     eventLog->eventInfo  = eventInfo;
223:     eventLog->maxEvents *= 2;
224:   }
225:   while (eventLog->numEvents < size) {
226:     PetscEventPerfInfoClear(&eventLog->eventInfo[eventLog->numEvents++]);
227:   }
228:   return(0);
229: }

231: #if defined(PETSC_HAVE_MPE)
232: #include <mpe.h>
233: PETSC_INTERN PetscErrorCode PetscLogMPEGetRGBColor(const char*[]);
234: PetscErrorCode PetscLogEventBeginMPE(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
235: {
236:   PetscErrorCode    ierr;

239:   MPE_Log_event(petsc_stageLog->eventLog->eventInfo[event].mpe_id_begin,0,NULL);
240:   return(0);
241: }

243: PetscErrorCode PetscLogEventEndMPE(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
244: {
245:   PetscErrorCode    ierr;

248:   MPE_Log_event(petsc_stageLog->eventLog->eventInfo[event].mpe_id_end,0,NULL);
249:   return(0);
250: }
251: #endif

253: /*--------------------------------------------- Registration Functions ----------------------------------------------*/
254: /*@C
255:   PetscEventRegLogRegister - Registers an event for logging operations in an application code.

257:   Not Collective

259:   Input Parameters:
260: + eventLog - The PetscEventLog
261: . ename    - The name associated with the event
262: - classid   - The classid associated to the class for this event

264:   Output Parameter:
265: . event    - The event

267:   Example of Usage:
268: .vb
269:       int USER_EVENT;
270:       PetscLogDouble user_event_flops;
271:       PetscLogEventRegister("User event name",0,&USER_EVENT);
272:       PetscLogEventBegin(USER_EVENT,0,0,0,0);
273:          [code segment to monitor]
274:          PetscLogFlops(user_event_flops);
275:       PetscLogEventEnd(USER_EVENT,0,0,0,0);
276: .ve

278:   Notes:

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

287:   Level: developer

289: .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogFlops(), 
290:           PetscEventLogActivate(), PetscEventLogDeactivate()
291: @*/
292: PetscErrorCode PetscEventRegLogRegister(PetscEventRegLog eventLog,const char ename[],PetscClassId classid,PetscLogEvent *event)
293: {
294:   PetscEventRegInfo *eventInfo;
295:   char              *str;
296:   int               e;
297:   PetscErrorCode    ierr;

302:   /* Should check classid I think */
303:   e = eventLog->numEvents++;
304:   if (eventLog->numEvents > eventLog->maxEvents) {
305:     PetscMalloc1(eventLog->maxEvents*2,&eventInfo);
306:     PetscMemcpy(eventInfo,eventLog->eventInfo,eventLog->maxEvents * sizeof(PetscEventRegInfo));
307:     PetscFree(eventLog->eventInfo);
308:     eventLog->eventInfo  = eventInfo;
309:     eventLog->maxEvents *= 2;
310:   }
311:   PetscStrallocpy(ename,&str);

313:   eventLog->eventInfo[e].name       = str;
314:   eventLog->eventInfo[e].classid    = classid;
315:   eventLog->eventInfo[e].collective = PETSC_TRUE;
316: #if defined(PETSC_HAVE_MPE)
317:   if (PetscLogPLB == PetscLogEventBeginMPE) {
318:     const char  *color;
319:     PetscMPIInt rank;
320:     int         beginID,endID;

322:     beginID = MPE_Log_get_event_number();
323:     endID   = MPE_Log_get_event_number();

325:     eventLog->eventInfo[e].mpe_id_begin = beginID;
326:     eventLog->eventInfo[e].mpe_id_end   = endID;

328:     MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
329:     if (!rank) {
330:       PetscLogMPEGetRGBColor(&color);
331:       MPE_Describe_state(beginID,endID,str,(char*)color);
332:     }
333:   }
334: #endif
335:   *event = e;
336:   return(0);
337: }

339: /*---------------------------------------------- Activation Functions -----------------------------------------------*/
340: /*@C
341:   PetscEventPerfLogActivate - Indicates that a particular event should be logged.

343:   Not Collective

345:   Input Parameters:
346: + eventLog - The PetscEventPerfLog
347: - event    - The event

349:    Usage:
350: .vb
351:       PetscEventPerfLogDeactivate(log, VEC_SetValues);
352:         [code where you do not want to log VecSetValues()]
353:       PetscEventPerfLogActivate(log, VEC_SetValues);
354:         [code where you do want to log VecSetValues()]
355: .ve

357:   Note:
358:   The event may be either a pre-defined PETSc event (found in
359:   include/petsclog.h) or an event number obtained with PetscEventRegLogRegister().

361:   Level: developer

363: .seealso: PetscEventPerfLogDeactivate()
364: @*/
365: PetscErrorCode PetscEventPerfLogActivate(PetscEventPerfLog eventLog,PetscLogEvent event)
366: {
368:   eventLog->eventInfo[event].active = PETSC_TRUE;
369:   return(0);
370: }

372: /*@C
373:   PetscEventPerfLogDeactivate - Indicates that a particular event should not be logged.

375:   Not Collective

377:   Input Parameters:
378: + eventLog - The PetscEventPerfLog
379: - event    - The event

381:    Usage:
382: .vb
383:       PetscEventPerfLogDeactivate(log, VEC_SetValues);
384:         [code where you do not want to log VecSetValues()]
385:       PetscEventPerfLogActivate(log, VEC_SetValues);
386:         [code where you do want to log VecSetValues()]
387: .ve

389:   Note:
390:   The event may be either a pre-defined PETSc event (found in
391:   include/petsclog.h) or an event number obtained with PetscEventRegLogRegister().

393:   Level: developer

395: .seealso: PetscEventPerfLogActivate()
396: @*/
397: PetscErrorCode PetscEventPerfLogDeactivate(PetscEventPerfLog eventLog,PetscLogEvent event)
398: {
400:   eventLog->eventInfo[event].active = PETSC_FALSE;
401:   return(0);
402: }

404: /*@C
405:   PetscEventPerfLogActivateClass - Activates event logging for a PETSc object class.

407:   Not Collective

409:   Input Parameters:
410: + eventLog    - The PetscEventPerfLog
411: . eventRegLog - The PetscEventRegLog
412: - classid      - The class id, for example MAT_CLASSID, SNES_CLASSID,

414:   Level: developer

416: .seealso: PetscEventPerfLogDeactivateClass(), PetscEventPerfLogActivate(), PetscEventPerfLogDeactivate()
417: @*/
418: PetscErrorCode PetscEventPerfLogActivateClass(PetscEventPerfLog eventLog,PetscEventRegLog eventRegLog,PetscClassId classid)
419: {
420:   int e;

423:   for (e = 0; e < eventLog->numEvents; e++) {
424:     int c = eventRegLog->eventInfo[e].classid;
425:     if (c == classid) eventLog->eventInfo[e].active = PETSC_TRUE;
426:   }
427:   return(0);
428: }

430: /*@C
431:   PetscEventPerfLogDeactivateClass - Deactivates event logging for a PETSc object class.

433:   Not Collective

435:   Input Parameters:
436: + eventLog    - The PetscEventPerfLog
437: . eventRegLog - The PetscEventRegLog
438: - classid - The class id, for example MAT_CLASSID, SNES_CLASSID,

440:   Level: developer

442: .seealso: PetscEventPerfLogDeactivateClass(), PetscEventPerfLogDeactivate(), PetscEventPerfLogActivate()
443: @*/
444: PetscErrorCode PetscEventPerfLogDeactivateClass(PetscEventPerfLog eventLog,PetscEventRegLog eventRegLog,PetscClassId classid)
445: {
446:   int e;

449:   for (e = 0; e < eventLog->numEvents; e++) {
450:     int c = eventRegLog->eventInfo[e].classid;
451:     if (c == classid) eventLog->eventInfo[e].active = PETSC_FALSE;
452:   }
453:   return(0);
454: }

456: /*------------------------------------------------ Query Functions --------------------------------------------------*/
457: /*@C
458:   PetscEventRegLogGetEvent - This function returns the event id given the event name.

460:   Not Collective

462:   Input Parameters:
463: + eventLog - The PetscEventRegLog
464: - name     - The stage name

466:   Output Parameter:
467: . event    - The event id, or -1 if not found

469:   Level: developer

471: .seealso: PetscEventRegLogRegister()
472: @*/
473: PetscErrorCode  PetscEventRegLogGetEvent(PetscEventRegLog eventLog,const char name[],PetscLogEvent *event)
474: {
475:   PetscBool      match;
476:   int            e;

482:   *event = -1;
483:   for (e = 0; e < eventLog->numEvents; e++) {
484:     PetscStrcasecmp(eventLog->eventInfo[e].name,name,&match);
485:     if (match) {
486:       *event = e;
487:       break;
488:     }
489:   }
490:   return(0);
491: }

493: /*@C
494:   PetscEventPerfLogSetVisible - This function determines whether an event is printed during PetscLogView()

496:   Not Collective

498:   Input Parameters:
499: + eventLog  - The PetscEventPerfLog
500: . event     - The event to log
501: - isVisible - The visibility flag, PETSC_TRUE for printing, otherwise PETSC_FALSE (default is PETSC_TRUE)

503:   Database Options:
504: . -log_view - Activates log summary

506:   Level: developer

508: .seealso: PetscEventPerfLogGetVisible(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
509: @*/
510: PetscErrorCode PetscEventPerfLogSetVisible(PetscEventPerfLog eventLog,PetscLogEvent event,PetscBool isVisible)
511: {
513:   eventLog->eventInfo[event].visible = isVisible;
514:   return(0);
515: }

517: /*@C
518:   PetscEventPerfLogGetVisible - This function returns whether an event is printed during PetscLogView()

520:   Not Collective

522:   Input Parameters:
523: + eventLog  - The PetscEventPerfLog
524: - event     - The event id to log

526:   Output Parameter:
527: . isVisible - The visibility flag, PETSC_TRUE for printing, otherwise PETSC_FALSE (default is PETSC_TRUE)

529:   Database Options:
530: . -log_view - Activates log summary

532:   Level: developer

534: .seealso: PetscEventPerfLogSetVisible(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
535: @*/
536: PetscErrorCode PetscEventPerfLogGetVisible(PetscEventPerfLog eventLog,PetscLogEvent event,PetscBool  *isVisible)
537: {
540:   *isVisible = eventLog->eventInfo[event].visible;
541:   return(0);
542: }

544: /*@C
545:   PetscLogEventGetPerfInfo - Return the performance information about the given event in the given stage

547:   Input Parameters:
548: + stage - The stage number or PETSC_DETERMINE for the current stage
549: - event - The event number

551:   Output Parameters:
552: . info - This structure is filled with the performance information

554:   Level: Intermediate

556: .seealso: PetscLogEventGetFlops()
557: @*/
558: PetscErrorCode PetscLogEventGetPerfInfo(int stage,PetscLogEvent event,PetscEventPerfInfo *info)
559: {
560:   PetscStageLog     stageLog;
561:   PetscEventPerfLog eventLog = NULL;
562:   PetscErrorCode    ierr;

566:   if (!PetscLogPLB) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Must use -log_view or PetscLogDefaultBegin() before calling this routine");
567:   PetscLogGetStageLog(&stageLog);
568:   if (stage < 0) {PetscStageLogGetCurrent(stageLog,&stage);}
569:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
570:   *info = eventLog->eventInfo[event];
571:   return(0);
572: }

574: PetscErrorCode PetscLogEventGetFlops(PetscLogEvent event,PetscLogDouble *flops)
575: {
576:   PetscStageLog     stageLog;
577:   PetscEventPerfLog eventLog = NULL;
578:   int               stage;
579:   PetscErrorCode    ierr;

582:   if (!PetscLogPLB) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Must use -log_view or PetscLogDefaultBegin() before calling this routine");
583:   PetscLogGetStageLog(&stageLog);
584:   PetscStageLogGetCurrent(stageLog,&stage);
585:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
586:   *flops = eventLog->eventInfo[event].flops;
587:   return(0);
588: }

590: PetscErrorCode PetscLogEventZeroFlops(PetscLogEvent event)
591: {
592:   PetscStageLog     stageLog;
593:   PetscEventPerfLog eventLog = NULL;
594:   int               stage;
595:   PetscErrorCode    ierr;

598:   PetscLogGetStageLog(&stageLog);
599:   PetscStageLogGetCurrent(stageLog,&stage);
600:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);

602:   eventLog->eventInfo[event].flops    = 0.0;
603:   eventLog->eventInfo[event].flops2   = 0.0;
604:   eventLog->eventInfo[event].flopsTmp = 0.0;
605:   return(0);
606: }

608: PetscErrorCode PetscLogEventSynchronize(PetscLogEvent event,MPI_Comm comm)
609: {
610:   PetscStageLog     stageLog;
611:   PetscEventRegLog  eventRegLog;
612:   PetscEventPerfLog eventLog = NULL;
613:   int               stage;
614:   PetscLogDouble    time = 0.0;
615:   PetscErrorCode    ierr;

618:   if (!PetscLogSyncOn || comm == MPI_COMM_NULL) return(0);
619:   PetscLogGetStageLog(&stageLog);
620:   PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
621:   if (!eventRegLog->eventInfo[event].collective) return(0);
622:   PetscStageLogGetCurrent(stageLog,&stage);
623:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
624:   if (eventLog->eventInfo[event].depth > 0) return(0);

626:   PetscTimeSubtract(&time);
627:   MPI_Barrier(comm);
628:   PetscTimeAdd(&time);
629:   eventLog->eventInfo[event].syncTime += time;
630:   return(0);
631: }

633: PetscErrorCode PetscLogEventBeginDefault(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
634: {
635:   PetscStageLog     stageLog;
636:   PetscEventPerfLog eventLog = NULL;
637:   int               stage;
638:   PetscErrorCode    ierr;

641:   PetscLogGetStageLog(&stageLog);
642:   PetscStageLogGetCurrent(stageLog,&stage);
643:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
644:   /* Synchronization */
645:   PetscLogEventSynchronize(event,PetscObjectComm(o1));
646:   /* Check for double counting */
647:   eventLog->eventInfo[event].depth++;
648:   if (eventLog->eventInfo[event].depth > 1) return(0);
649:   /* Log the performance info */
650:   eventLog->eventInfo[event].count++;
651:   eventLog->eventInfo[event].timeTmp = 0.0;
652:   PetscTimeSubtract(&eventLog->eventInfo[event].timeTmp);
653:   eventLog->eventInfo[event].flopsTmp       = 0.0;
654:   eventLog->eventInfo[event].flopsTmp      -= petsc_TotalFlops;
655:   eventLog->eventInfo[event].numMessages   -= petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
656:   eventLog->eventInfo[event].messageLength -= petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
657:   eventLog->eventInfo[event].numReductions -= petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
658:   if (PetscLogMemory) {
659:     PetscLogDouble usage;
660:     PetscMemoryGetCurrentUsage(&usage);
661:     eventLog->eventInfo[event].memIncrease -= usage;
662:     PetscMallocGetCurrentUsage(&usage);
663:     eventLog->eventInfo[event].mallocSpace -= usage;
664:     PetscMallocGetMaximumUsage(&usage);
665:     eventLog->eventInfo[event].mallocIncrease -= usage;
666:     PetscMallocPushMaximumUsage((int)event);
667:   }
668:   #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) 
669:   eventLog->eventInfo[event].CpuToGpuCount -= petsc_ctog_ct;
670:   eventLog->eventInfo[event].GpuToCpuCount -= petsc_gtoc_ct;
671:   eventLog->eventInfo[event].CpuToGpuSize  -= petsc_ctog_sz;
672:   eventLog->eventInfo[event].GpuToCpuSize  -= petsc_gtoc_sz;
673:   #endif
674:   return(0);
675: }

677: PetscErrorCode PetscLogEventEndDefault(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
678: {
679:   PetscStageLog     stageLog;
680:   PetscEventPerfLog eventLog = NULL;
681:   int               stage;
682:   PetscErrorCode    ierr;

685:   PetscLogGetStageLog(&stageLog);
686:   PetscStageLogGetCurrent(stageLog,&stage);
687:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
688:   /* Check for double counting */
689:   eventLog->eventInfo[event].depth--;
690:   if (eventLog->eventInfo[event].depth > 0) return(0);
691:   else if (eventLog->eventInfo[event].depth < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Logging event had unbalanced begin/end pairs");
692:   /* Log performance info */
693:   PetscTimeAdd(&eventLog->eventInfo[event].timeTmp);
694:   eventLog->eventInfo[event].time          += eventLog->eventInfo[event].timeTmp;
695:   eventLog->eventInfo[event].time2         += eventLog->eventInfo[event].timeTmp*eventLog->eventInfo[event].timeTmp;
696:   eventLog->eventInfo[event].flopsTmp      += petsc_TotalFlops;
697:   eventLog->eventInfo[event].flops         += eventLog->eventInfo[event].flopsTmp;
698:   eventLog->eventInfo[event].flops2        += eventLog->eventInfo[event].flopsTmp*eventLog->eventInfo[event].flopsTmp;
699:   eventLog->eventInfo[event].numMessages   += petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
700:   eventLog->eventInfo[event].messageLength += petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
701:   eventLog->eventInfo[event].numReductions += petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
702:   if (PetscLogMemory) {
703:     PetscLogDouble usage,musage;
704:     PetscMemoryGetCurrentUsage(&usage);
705:     eventLog->eventInfo[event].memIncrease += usage;
706:     PetscMallocGetCurrentUsage(&usage);
707:     eventLog->eventInfo[event].mallocSpace += usage;
708:     PetscMallocPopMaximumUsage((int)event,&musage);
709:     eventLog->eventInfo[event].mallocIncreaseEvent = PetscMax(musage-usage,eventLog->eventInfo[event].mallocIncreaseEvent);
710:     PetscMallocGetMaximumUsage(&usage);
711:     eventLog->eventInfo[event].mallocIncrease += usage;
712:   }
713:   #if defined(PETSC_HAVE_VIENNACL) || defined(PETSC_HAVE_CUDA) 
714:   eventLog->eventInfo[event].CpuToGpuCount += petsc_ctog_ct;
715:   eventLog->eventInfo[event].GpuToCpuCount += petsc_gtoc_ct;
716:   eventLog->eventInfo[event].CpuToGpuSize  += petsc_ctog_sz;
717:   eventLog->eventInfo[event].GpuToCpuSize  += petsc_gtoc_sz;
718:   #endif
719:   return(0);
720: }

722: PetscErrorCode PetscLogEventBeginComplete(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
723: {
724:   PetscStageLog     stageLog;
725:   PetscEventRegLog  eventRegLog;
726:   PetscEventPerfLog eventPerfLog = NULL;
727:   Action            *tmpAction;
728:   PetscLogDouble    start,end;
729:   PetscLogDouble    curTime;
730:   int               stage;
731:   PetscErrorCode    ierr;

734:   /* Dynamically enlarge logging structures */
735:   if (petsc_numActions >= petsc_maxActions) {
736:     PetscTime(&start);
737:     PetscMalloc1(petsc_maxActions*2,&tmpAction);
738:     PetscMemcpy(tmpAction,petsc_actions,petsc_maxActions * sizeof(Action));
739:     PetscFree(petsc_actions);

741:     petsc_actions     = tmpAction;
742:     petsc_maxActions *= 2;
743:     PetscTime(&end);
744:     petsc_BaseTime += (end - start);
745:   }
746:   /* Record the event */
747:   PetscLogGetStageLog(&stageLog);
748:   PetscStageLogGetCurrent(stageLog,&stage);
749:   PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
750:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
751:   PetscTime(&curTime);
752:   if (petsc_logActions) {
753:     petsc_actions[petsc_numActions].time    = curTime - petsc_BaseTime;
754:     petsc_actions[petsc_numActions].action  = ACTIONBEGIN;
755:     petsc_actions[petsc_numActions].event   = event;
756:     petsc_actions[petsc_numActions].classid = eventRegLog->eventInfo[event].classid;
757:     if (o1) petsc_actions[petsc_numActions].id1 = o1->id;
758:     else petsc_actions[petsc_numActions].id1 = -1;
759:     if (o2) petsc_actions[petsc_numActions].id2 = o2->id;
760:     else petsc_actions[petsc_numActions].id2 = -1;
761:     if (o3) petsc_actions[petsc_numActions].id3 = o3->id;
762:     else petsc_actions[petsc_numActions].id3 = -1;
763:     petsc_actions[petsc_numActions].flops = petsc_TotalFlops;

765:     PetscMallocGetCurrentUsage(&petsc_actions[petsc_numActions].mem);
766:     PetscMallocGetMaximumUsage(&petsc_actions[petsc_numActions].maxmem);
767:     petsc_numActions++;
768:   }
769:   /* Check for double counting */
770:   eventPerfLog->eventInfo[event].depth++;
771:   if (eventPerfLog->eventInfo[event].depth > 1) return(0);
772:   /* Log the performance info */
773:   eventPerfLog->eventInfo[event].count++;
774:   eventPerfLog->eventInfo[event].time          -= curTime;
775:   eventPerfLog->eventInfo[event].flops         -= petsc_TotalFlops;
776:   eventPerfLog->eventInfo[event].numMessages   -= petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
777:   eventPerfLog->eventInfo[event].messageLength -= petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
778:   eventPerfLog->eventInfo[event].numReductions -= petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
779:   return(0);
780: }

782: PetscErrorCode PetscLogEventEndComplete(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
783: {
784:   PetscStageLog     stageLog;
785:   PetscEventRegLog  eventRegLog;
786:   PetscEventPerfLog eventPerfLog = NULL;
787:   Action            *tmpAction;
788:   PetscLogDouble    start,end;
789:   PetscLogDouble    curTime;
790:   int               stage;
791:   PetscErrorCode    ierr;

794:   /* Dynamically enlarge logging structures */
795:   if (petsc_numActions >= petsc_maxActions) {
796:     PetscTime(&start);
797:     PetscMalloc1(petsc_maxActions*2,&tmpAction);
798:     PetscMemcpy(tmpAction,petsc_actions,petsc_maxActions * sizeof(Action));
799:     PetscFree(petsc_actions);

801:     petsc_actions     = tmpAction;
802:     petsc_maxActions *= 2;
803:     PetscTime(&end);
804:     petsc_BaseTime += (end - start);
805:   }
806:   /* Record the event */
807:   PetscLogGetStageLog(&stageLog);
808:   PetscStageLogGetCurrent(stageLog,&stage);
809:   PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
810:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
811:   PetscTime(&curTime);
812:   if (petsc_logActions) {
813:     petsc_actions[petsc_numActions].time    = curTime - petsc_BaseTime;
814:     petsc_actions[petsc_numActions].action  = ACTIONEND;
815:     petsc_actions[petsc_numActions].event   = event;
816:     petsc_actions[petsc_numActions].classid = eventRegLog->eventInfo[event].classid;
817:     if (o1) petsc_actions[petsc_numActions].id1 = o1->id;
818:     else petsc_actions[petsc_numActions].id1 = -1;
819:     if (o2) petsc_actions[petsc_numActions].id2 = o2->id;
820:     else petsc_actions[petsc_numActions].id2 = -1;
821:     if (o3) petsc_actions[petsc_numActions].id3 = o3->id;
822:     else petsc_actions[petsc_numActions].id3 = -1;
823:     petsc_actions[petsc_numActions].flops = petsc_TotalFlops;

825:     PetscMallocGetCurrentUsage(&petsc_actions[petsc_numActions].mem);
826:     PetscMallocGetMaximumUsage(&petsc_actions[petsc_numActions].maxmem);
827:     petsc_numActions++;
828:   }
829:   /* Check for double counting */
830:   eventPerfLog->eventInfo[event].depth--;
831:   if (eventPerfLog->eventInfo[event].depth > 0) return(0);
832:   else if (eventPerfLog->eventInfo[event].depth < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Logging event had unbalanced begin/end pairs");
833:   /* Log the performance info */
834:   eventPerfLog->eventInfo[event].count++;
835:   eventPerfLog->eventInfo[event].time          += curTime;
836:   eventPerfLog->eventInfo[event].flops         += petsc_TotalFlops;
837:   eventPerfLog->eventInfo[event].numMessages   += petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
838:   eventPerfLog->eventInfo[event].messageLength += petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
839:   eventPerfLog->eventInfo[event].numReductions += petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
840:   return(0);
841: }

843: PetscErrorCode PetscLogEventBeginTrace(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
844: {
845:   PetscStageLog     stageLog;
846:   PetscEventRegLog  eventRegLog;
847:   PetscEventPerfLog eventPerfLog = NULL;
848:   PetscLogDouble    cur_time;
849:   PetscMPIInt       rank;
850:   int               stage,err;
851:   PetscErrorCode    ierr;

854:   if (!petsc_tracetime) PetscTime(&petsc_tracetime);

856:   petsc_tracelevel++;
857:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
858:   PetscLogGetStageLog(&stageLog);
859:   PetscStageLogGetCurrent(stageLog,&stage);
860:   PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
861:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
862:   /* Check for double counting */
863:   eventPerfLog->eventInfo[event].depth++;
864:   if (eventPerfLog->eventInfo[event].depth > 1) return(0);
865:   /* Log performance info */
866:   PetscTime(&cur_time);
867:   PetscFPrintf(PETSC_COMM_SELF,petsc_tracefile,"%s[%d] %g Event begin: %s\n",petsc_tracespace,rank,cur_time-petsc_tracetime,eventRegLog->eventInfo[event].name);
868:   PetscStrncpy(petsc_tracespace,petsc_traceblanks,2*petsc_tracelevel);

870:   petsc_tracespace[2*petsc_tracelevel] = 0;

872:   err = fflush(petsc_tracefile);
873:   if (err) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fflush() failed on file");
874:   return(0);
875: }

877: PetscErrorCode PetscLogEventEndTrace(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
878: {
879:   PetscStageLog     stageLog;
880:   PetscEventRegLog  eventRegLog;
881:   PetscEventPerfLog eventPerfLog = NULL;
882:   PetscLogDouble    cur_time;
883:   int               stage,err;
884:   PetscMPIInt       rank;
885:   PetscErrorCode    ierr;

888:   petsc_tracelevel--;
889:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
890:   PetscLogGetStageLog(&stageLog);
891:   PetscStageLogGetCurrent(stageLog,&stage);
892:   PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
893:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
894:   /* Check for double counting */
895:   eventPerfLog->eventInfo[event].depth--;
896:   if (eventPerfLog->eventInfo[event].depth > 0) return(0);
897:   else if (eventPerfLog->eventInfo[event].depth < 0 || petsc_tracelevel < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Logging event had unbalanced begin/end pairs");

899:   /* Log performance info */
900:   if (petsc_tracelevel) {
901:     PetscStrncpy(petsc_tracespace,petsc_traceblanks,2*petsc_tracelevel);
902:   }
903:   petsc_tracespace[2*petsc_tracelevel] = 0;
904:   PetscTime(&cur_time);
905:   PetscFPrintf(PETSC_COMM_SELF,petsc_tracefile,"%s[%d] %g Event end: %s\n",petsc_tracespace,rank,cur_time-petsc_tracetime,eventRegLog->eventInfo[event].name);
906:   err  = fflush(petsc_tracefile);
907:   if (err) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SYS,"fflush() failed on file");
908:   return(0);
909: }

911: /*@C
912:   PetscLogEventSetDof - Set the nth number of degrees of freedom associated with this event

914:   Not Collective

916:   Input Parameters:
917: + event - The event id to log
918: . n     - The dof index, in [0, 8)
919: - dof   - The number of dofs

921:   Database Options:
922: . -log_view - Activates log summary

924:   Note: This is to enable logging of convergence

926:   Level: developer

928: .seealso: PetscLogEventSetError(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
929: @*/
930: PetscErrorCode PetscLogEventSetDof(PetscLogEvent event, PetscInt n, PetscLogDouble dof)
931: {
932:   PetscStageLog     stageLog;
933:   PetscEventPerfLog eventLog = NULL;
934:   int               stage;
935:   PetscErrorCode    ierr;

938:   if ((n < 0) || (n > 7)) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Error index %D is not in [0, 8)", n);
939:   PetscLogGetStageLog(&stageLog);
940:   PetscStageLogGetCurrent(stageLog,&stage);
941:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
942:   eventLog->eventInfo[event].dof[n] = dof;
943:   return(0);
944: }

946: /*@C
947:   PetscLogEventSetError - Set the nth error associated with this event

949:   Not Collective

951:   Input Parameters:
952: + event - The event id to log
953: . n     - The error index, in [0, 8)
954: . error - The error

956:   Database Options:
957: . -log_view - Activates log summary

959:   Note: This is to enable logging of convergence, and enable users to interpret the errors as they wish. For example,
960:   as different norms, or as errors for different fields

962:   Level: developer

964: .seealso: PetscLogEventSetDof(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
965: @*/
966: PetscErrorCode PetscLogEventSetError(PetscLogEvent event, PetscInt n, PetscLogDouble error)
967: {
968:   PetscStageLog     stageLog;
969:   PetscEventPerfLog eventLog = NULL;
970:   int               stage;
971:   PetscErrorCode    ierr;

974:   if ((n < 0) || (n > 7)) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Error index %D is not in [0, 8)", n);
975:   PetscLogGetStageLog(&stageLog);
976:   PetscStageLogGetCurrent(stageLog,&stage);
977:   PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
978:   eventLog->eventInfo[event].errors[n] = error;
979:   return(0);
980: }