Actual source code: stagelog.c
  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: PetscStageLog petsc_stageLog = NULL;
 12: /*@C
 13:   PetscLogGetStageLog - This function returns the default stage logging object.
 15:   Not collective
 17:   Output Parameter:
 18: . stageLog - The default PetscStageLog
 20:   Level: developer
 22:   Developer Notes:
 23:     Inline since called for EACH PetscEventLogBeginDefault() and PetscEventLogEndDefault()
 25: .seealso: PetscStageLogCreate()
 26: @*/
 27: PetscErrorCode PetscLogGetStageLog(PetscStageLog *stageLog)
 28: {
 31:   if (!petsc_stageLog) {
 32:     fprintf(stderr, "PETSC ERROR: Logging has not been enabled.\nYou might have forgotten to call PetscInitialize().\n");
 33:     PETSCABORT(MPI_COMM_WORLD, PETSC_ERR_SUP);
 34:   }
 35:   *stageLog = petsc_stageLog;
 36:   return(0);
 37: }
 39: /*@C
 40:   PetscStageLogGetCurrent - This function returns the stage from the top of the stack.
 42:   Not Collective
 44:   Input Parameter:
 45: . stageLog - The PetscStageLog
 47:   Output Parameter:
 48: . stage    - The current stage
 50:   Notes:
 51:   If no stage is currently active, stage is set to -1.
 53:   Level: developer
 55:   Developer Notes:
 56:     Inline since called for EACH PetscEventLogBeginDefault() and PetscEventLogEndDefault()
 58: .seealso: PetscStageLogPush(), PetscStageLogPop(), PetscLogGetStageLog()
 59: @*/
 60: PetscErrorCode  PetscStageLogGetCurrent(PetscStageLog stageLog, int *stage)
 61: {
 62:   PetscBool      empty;
 66:   PetscIntStackEmpty(stageLog->stack, &empty);
 67:   if (empty) {
 68:     *stage = -1;
 69:   } else {
 70:     PetscIntStackTop(stageLog->stack, stage);
 71:   }
 72:   if (*stage != stageLog->curStage) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB, "Inconsistency in stage log: stage %d should be %d", *stage, stageLog->curStage);
 73:   return(0);
 74: }
 76: /*@C
 77:   PetscStageLogGetEventPerfLog - This function returns the PetscEventPerfLog for the given stage.
 79:   Not Collective
 81:   Input Parameters:
 82: + stageLog - The PetscStageLog
 83: - stage    - The stage
 85:   Output Parameter:
 86: . eventLog - The PetscEventPerfLog
 88:   Level: developer
 90:   Developer Notes:
 91:     Inline since called for EACH PetscEventLogBeginDefault() and PetscEventLogEndDefault()
 93: .seealso: PetscStageLogPush(), PetscStageLogPop(), PetscLogGetStageLog()
 94: @*/
 95: PetscErrorCode  PetscStageLogGetEventPerfLog(PetscStageLog stageLog, int stage, PetscEventPerfLog *eventLog)
 96: {
 99:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
100:   *eventLog = stageLog->stageInfo[stage].eventLog;
101:   return(0);
102: }
104: /*@C
105:   PetscStageInfoDestroy - This destroys a PetscStageInfo object.
107:   Not collective
109:   Input Parameter:
110: . stageInfo - The PetscStageInfo
112:   Level: developer
114: .seealso: PetscStageLogCreate()
115: @*/
116: PetscErrorCode  PetscStageInfoDestroy(PetscStageInfo *stageInfo)
117: {
121:   PetscFree(stageInfo->name);
122:   PetscEventPerfLogDestroy(stageInfo->eventLog);
123:   PetscClassPerfLogDestroy(stageInfo->classLog);
124:   return(0);
125: }
127: /*@C
128:   PetscStageLogDestroy - This destroys a PetscStageLog object.
130:   Not collective
132:   Input Parameter:
133: . stageLog - The PetscStageLog
135:   Level: developer
137: .seealso: PetscStageLogCreate()
138: @*/
139: PetscErrorCode  PetscStageLogDestroy(PetscStageLog stageLog)
140: {
141:   int            stage;
145:   if (!stageLog) return(0);
146:   PetscIntStackDestroy(stageLog->stack);
147:   PetscEventRegLogDestroy(stageLog->eventLog);
148:   PetscClassRegLogDestroy(stageLog->classLog);
149:   for (stage = 0; stage < stageLog->numStages; stage++) {
150:     PetscStageInfoDestroy(&stageLog->stageInfo[stage]);
151:   }
152:   PetscFree(stageLog->stageInfo);
153:   PetscFree(stageLog);
154:   return(0);
155: }
157: /*@C
158:   PetscStageLogRegister - Registers a stage name for logging operations in an application code.
160:   Not Collective
162:   Input Parameter:
163: + stageLog - The PetscStageLog
164: - sname    - the name to associate with that stage
166:   Output Parameter:
167: . stage    - The stage index
169:   Level: developer
171: .seealso: PetscStageLogPush(), PetscStageLogPop(), PetscStageLogCreate()
172: @*/
173: PetscErrorCode  PetscStageLogRegister(PetscStageLog stageLog, const char sname[], int *stage)
174: {
175:   PetscStageInfo *stageInfo;
176:   int            s;
182:   /* Check stage already registered */
183:   for (s = 0; s < stageLog->numStages; ++s) {
184:     PetscBool same;
186:     PetscStrcmp(stageLog->stageInfo[s].name, sname, &same);
187:     if (same) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONG, "Duplicate stage name given: %s", sname);
188:   }
189:   /* Create new stage */
190:   s = stageLog->numStages++;
191:   if (stageLog->numStages > stageLog->maxStages) {
192:     PetscMalloc1(stageLog->maxStages*2, &stageInfo);
193:     PetscArraycpy(stageInfo, stageLog->stageInfo, stageLog->maxStages);
194:     PetscFree(stageLog->stageInfo);
195:     stageLog->stageInfo  = stageInfo;
196:     stageLog->maxStages *= 2;
197:   }
198:   /* Setup new stage info */
199:   stageInfo = &stageLog->stageInfo[s];
200:   PetscMemzero(stageInfo,sizeof(PetscStageInfo));
201:   PetscStrallocpy(sname,&stageInfo->name);
202:   stageInfo->used             = PETSC_FALSE;
203:   stageInfo->perfInfo.active  = PETSC_TRUE;
204:   stageInfo->perfInfo.visible = PETSC_TRUE;
205:   PetscEventPerfLogCreate(&stageInfo->eventLog);
206:   PetscClassPerfLogCreate(&stageInfo->classLog);
207:   *stage = s;
208:   return(0);
209: }
211: /*@C
212:   PetscStageLogPush - This function pushes a stage on the stack.
214:   Not Collective
216:   Input Parameters:
217: + stageLog   - The PetscStageLog
218: - stage - The stage to log
220:   Database Options:
221: . -log_view - Activates logging
223:   Usage:
224:   If the option -log_sumary is used to run the program containing the
225:   following code, then 2 sets of summary data will be printed during
226:   PetscFinalize().
227: .vb
228:       PetscInitialize(int *argc,char ***args,0,0);
229:       [stage 0 of code]
230:       PetscStageLogPush(stageLog,1);
231:       [stage 1 of code]
232:       PetscStageLogPop(stageLog);
233:       PetscBarrier(...);
234:       [more stage 0 of code]
235:       PetscFinalize();
236: .ve
238:   Notes:
239:   Use PetscLogStageRegister() to register a stage. All previous stages are
240:   accumulating time and flops, but events will only be logged in this stage.
242:   Level: developer
244: .seealso: PetscStageLogPop(), PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
245: @*/
246: PetscErrorCode  PetscStageLogPush(PetscStageLog stageLog, int stage)
247: {
248:   int            curStage = 0;
249:   PetscBool      empty;
253:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
255:   /* Record flops/time of previous stage */
256:   PetscIntStackEmpty(stageLog->stack, &empty);
257:   if (!empty) {
258:     PetscIntStackTop(stageLog->stack, &curStage);
259:     if (stageLog->stageInfo[curStage].perfInfo.active) {
260:       PetscTimeAdd(&stageLog->stageInfo[curStage].perfInfo.time);
261:       stageLog->stageInfo[curStage].perfInfo.flops         += petsc_TotalFlops;
262:       stageLog->stageInfo[curStage].perfInfo.numMessages   += petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
263:       stageLog->stageInfo[curStage].perfInfo.messageLength += petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
264:       stageLog->stageInfo[curStage].perfInfo.numReductions += petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
265:     }
266:   }
267:   /* Activate the stage */
268:   PetscIntStackPush(stageLog->stack, stage);
270:   stageLog->stageInfo[stage].used = PETSC_TRUE;
271:   stageLog->stageInfo[stage].perfInfo.count++;
272:   stageLog->curStage = stage;
273:   /* Subtract current quantities so that we obtain the difference when we pop */
274:   if (stageLog->stageInfo[stage].perfInfo.active) {
275:     PetscTimeSubtract(&stageLog->stageInfo[stage].perfInfo.time);
276:     stageLog->stageInfo[stage].perfInfo.flops         -= petsc_TotalFlops;
277:     stageLog->stageInfo[stage].perfInfo.numMessages   -= petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
278:     stageLog->stageInfo[stage].perfInfo.messageLength -= petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
279:     stageLog->stageInfo[stage].perfInfo.numReductions -= petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
280:   }
281:   return(0);
282: }
284: /*@C
285:   PetscStageLogPop - This function pops a stage from the stack.
287:   Not Collective
289:   Input Parameter:
290: . stageLog - The PetscStageLog
292:   Usage:
293:   If the option -log_sumary is used to run the program containing the
294:   following code, then 2 sets of summary data will be printed during
295:   PetscFinalize().
296: .vb
297:       PetscInitialize(int *argc,char ***args,0,0);
298:       [stage 0 of code]
299:       PetscStageLogPush(stageLog,1);
300:       [stage 1 of code]
301:       PetscStageLogPop(stageLog);
302:       PetscBarrier(...);
303:       [more stage 0 of code]
304:       PetscFinalize();
305: .ve
307:   Notes:
308:   Use PetscStageLogRegister() to register a stage.
310:   Level: developer
312: .seealso: PetscStageLogPush(), PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
313: @*/
314: PetscErrorCode  PetscStageLogPop(PetscStageLog stageLog)
315: {
316:   int            curStage;
317:   PetscBool      empty;
321:   /* Record flops/time of current stage */
322:   PetscIntStackPop(stageLog->stack, &curStage);
323:   if (stageLog->stageInfo[curStage].perfInfo.active) {
324:     PetscTimeAdd(&stageLog->stageInfo[curStage].perfInfo.time);
325:     stageLog->stageInfo[curStage].perfInfo.flops         += petsc_TotalFlops;
326:     stageLog->stageInfo[curStage].perfInfo.numMessages   += petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
327:     stageLog->stageInfo[curStage].perfInfo.messageLength += petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
328:     stageLog->stageInfo[curStage].perfInfo.numReductions += petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
329:   }
330:   PetscIntStackEmpty(stageLog->stack, &empty);
331:   if (!empty) {
332:     /* Subtract current quantities so that we obtain the difference when we pop */
333:     PetscIntStackTop(stageLog->stack, &curStage);
334:     if (stageLog->stageInfo[curStage].perfInfo.active) {
335:       PetscTimeSubtract(&stageLog->stageInfo[curStage].perfInfo.time);
336:       stageLog->stageInfo[curStage].perfInfo.flops         -= petsc_TotalFlops;
337:       stageLog->stageInfo[curStage].perfInfo.numMessages   -= petsc_irecv_ct  + petsc_isend_ct  + petsc_recv_ct  + petsc_send_ct;
338:       stageLog->stageInfo[curStage].perfInfo.messageLength -= petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
339:       stageLog->stageInfo[curStage].perfInfo.numReductions -= petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
340:     }
341:     stageLog->curStage = curStage;
342:   } else stageLog->curStage = -1;
343:   return(0);
344: }
347: /*@C
348:   PetscStageLogGetClassRegLog - This function returns the PetscClassRegLog for the given stage.
350:   Not Collective
352:   Input Parameters:
353: . stageLog - The PetscStageLog
355:   Output Parameter:
356: . classLog - The PetscClassRegLog
358:   Level: developer
360: .seealso: PetscStageLogPush(), PetscStageLogPop(), PetscLogGetStageLog()
361: @*/
362: PetscErrorCode  PetscStageLogGetClassRegLog(PetscStageLog stageLog, PetscClassRegLog *classLog)
363: {
366:   *classLog = stageLog->classLog;
367:   return(0);
368: }
370: /*@C
371:   PetscStageLogGetEventRegLog - This function returns the PetscEventRegLog.
373:   Not Collective
375:   Input Parameters:
376: . stageLog - The PetscStageLog
378:   Output Parameter:
379: . eventLog - The PetscEventRegLog
381:   Level: developer
383: .seealso: PetscStageLogPush(), PetscStageLogPop(), PetscLogGetStageLog()
384: @*/
385: PetscErrorCode  PetscStageLogGetEventRegLog(PetscStageLog stageLog, PetscEventRegLog *eventLog)
386: {
389:   *eventLog = stageLog->eventLog;
390:   return(0);
391: }
393: /*@C
394:   PetscStageLogGetClassPerfLog - This function returns the PetscClassPerfLog for the given stage.
396:   Not Collective
398:   Input Parameters:
399: + stageLog - The PetscStageLog
400: - stage    - The stage
402:   Output Parameter:
403: . classLog - The PetscClassPerfLog
405:   Level: developer
407: .seealso: PetscStageLogPush(), PetscStageLogPop(), PetscLogGetStageLog()
408: @*/
409: PetscErrorCode  PetscStageLogGetClassPerfLog(PetscStageLog stageLog, int stage, PetscClassPerfLog *classLog)
410: {
413:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
414:   *classLog = stageLog->stageInfo[stage].classLog;
415:   return(0);
416: }
419: /*@C
420:   PetscStageLogSetActive - This function determines whether events will be logged during this state.
422:   Not Collective
424:   Input Parameters:
425: + stageLog - The PetscStageLog
426: . stage    - The stage to log
427: - isActive - The activity flag, PETSC_TRUE for logging, otherwise PETSC_FALSE (default is PETSC_TRUE)
429:   Level: developer
431: .seealso: PetscStageLogGetActive(), PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
432: @*/
433: PetscErrorCode  PetscStageLogSetActive(PetscStageLog stageLog, int stage, PetscBool isActive)
434: {
436:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
437:   stageLog->stageInfo[stage].perfInfo.active = isActive;
438:   return(0);
439: }
441: /*@C
442:   PetscStageLogGetActive - This function returns whether events will be logged suring this stage.
444:   Not Collective
446:   Input Parameters:
447: + stageLog - The PetscStageLog
448: - stage    - The stage to log
450:   Output Parameter:
451: . isActive - The activity flag, PETSC_TRUE for logging, otherwise PETSC_FALSE (default is PETSC_TRUE)
453:   Level: developer
455: .seealso: PetscStageLogSetActive(), PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
456: @*/
457: PetscErrorCode  PetscStageLogGetActive(PetscStageLog stageLog, int stage, PetscBool  *isActive)
458: {
460:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
462:   *isActive = stageLog->stageInfo[stage].perfInfo.active;
463:   return(0);
464: }
466: /*@C
467:   PetscStageLogSetVisible - This function determines whether a stage is printed during PetscLogView()
469:   Not Collective
471:   Input Parameters:
472: + stageLog  - The PetscStageLog
473: . stage     - The stage to log
474: - isVisible - The visibility flag, PETSC_TRUE for printing, otherwise PETSC_FALSE (default is PETSC_TRUE)
476:   Database Options:
477: . -log_view - Activates log summary
479:   Level: developer
481: .seealso: PetscStageLogGetVisible(), PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
482: @*/
483: PetscErrorCode  PetscStageLogSetVisible(PetscStageLog stageLog, int stage, PetscBool isVisible)
484: {
486:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
487:   stageLog->stageInfo[stage].perfInfo.visible = isVisible;
488:   return(0);
489: }
491: /*@C
492:   PetscStageLogGetVisible - This function returns whether a stage is printed during PetscLogView()
494:   Not Collective
496:   Input Parameters:
497: + stageLog  - The PetscStageLog
498: - stage     - The stage to log
500:   Output Parameter:
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: PetscStageLogSetVisible(), PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
509: @*/
510: PetscErrorCode  PetscStageLogGetVisible(PetscStageLog stageLog, int stage, PetscBool  *isVisible)
511: {
513:   if ((stage < 0) || (stage >= stageLog->numStages)) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE, "Invalid stage %d should be in [0,%d)", stage, stageLog->numStages);
515:   *isVisible = stageLog->stageInfo[stage].perfInfo.visible;
516:   return(0);
517: }
519: /*@C
520:   PetscStageLogGetStage - This function returns the stage id given the stage name.
522:   Not Collective
524:   Input Parameters:
525: + stageLog - The PetscStageLog
526: - name     - The stage name
528:   Output Parameter:
529: . stage    - The stage id, or -1 if it does not exist
531:   Level: developer
533: .seealso: PetscStageLogGetCurrent(), PetscStageLogRegister(), PetscLogGetStageLog()
534: @*/
535: PetscErrorCode  PetscStageLogGetStage(PetscStageLog stageLog, const char name[], PetscLogStage *stage)
536: {
537:   PetscBool      match;
538:   int            s;
544:   *stage = -1;
545:   for (s = 0; s < stageLog->numStages; s++) {
546:     PetscStrcasecmp(stageLog->stageInfo[s].name, name, &match);
547:     if (match) {
548:       *stage = s;
549:       break;
550:     }
551:   }
552:   return(0);
553: }
555: /*@C
556:   PetscStageLogCreate - This creates a PetscStageLog object.
558:   Not collective
560:   Input Parameter:
561: . stageLog - The PetscStageLog
563:   Level: developer
565: .seealso: PetscStageLogCreate()
566: @*/
567: PetscErrorCode  PetscStageLogCreate(PetscStageLog *stageLog)
568: {
569:   PetscStageLog  l;
573:   PetscNew(&l);
575:   l->numStages = 0;
576:   l->maxStages = 10;
577:   l->curStage  = -1;
579:   PetscIntStackCreate(&l->stack);
580:   PetscMalloc1(l->maxStages, &l->stageInfo);
581:   PetscEventRegLogCreate(&l->eventLog);
582:   PetscClassRegLogCreate(&l->classLog);
584:   *stageLog = l;
585:   return(0);
586: }