Actual source code: plog.c

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

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

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

  9:       ***

 11:       This file, and only this file, is for functions that interact with the global logging state
 12: */
 13: #include <petsc/private/logimpl.h>
 14: #include <petsc/private/loghandlerimpl.h>
 15: #include <petsctime.h>
 16: #include <petscviewer.h>
 17: #include <petscdevice.h>
 18: #include <petsc/private/deviceimpl.h>

 20: #if defined(PETSC_HAVE_THREADSAFETY)

 22: PetscInt           petsc_log_gid = -1; /* Global threadId counter */
 23: PETSC_TLS PetscInt petsc_log_tid = -1; /* Local threadId */

 25: /* shared variables */
 26: PetscSpinlock PetscLogSpinLock;

 28: PetscInt PetscLogGetTid(void)
 29: {
 30:   if (petsc_log_tid < 0) {
 31:     PetscCall(PetscSpinlockLock(&PetscLogSpinLock));
 32:     petsc_log_tid = ++petsc_log_gid;
 33:     PetscCall(PetscSpinlockUnlock(&PetscLogSpinLock));
 34:   }
 35:   return petsc_log_tid;
 36: }

 38: #endif

 40: /* Global counters */
 41: PetscLogDouble petsc_BaseTime        = 0.0;
 42: PetscLogDouble petsc_TotalFlops      = 0.0; /* The number of flops */
 43: PetscLogDouble petsc_send_ct         = 0.0; /* The number of sends */
 44: PetscLogDouble petsc_recv_ct         = 0.0; /* The number of receives */
 45: PetscLogDouble petsc_send_len        = 0.0; /* The total length of all sent messages */
 46: PetscLogDouble petsc_recv_len        = 0.0; /* The total length of all received messages */
 47: PetscLogDouble petsc_isend_ct        = 0.0; /* The number of immediate sends */
 48: PetscLogDouble petsc_irecv_ct        = 0.0; /* The number of immediate receives */
 49: PetscLogDouble petsc_isend_len       = 0.0; /* The total length of all immediate send messages */
 50: PetscLogDouble petsc_irecv_len       = 0.0; /* The total length of all immediate receive messages */
 51: PetscLogDouble petsc_wait_ct         = 0.0; /* The number of waits */
 52: PetscLogDouble petsc_wait_any_ct     = 0.0; /* The number of anywaits */
 53: PetscLogDouble petsc_wait_all_ct     = 0.0; /* The number of waitalls */
 54: PetscLogDouble petsc_sum_of_waits_ct = 0.0; /* The total number of waits */
 55: PetscLogDouble petsc_allreduce_ct    = 0.0; /* The number of reductions */
 56: PetscLogDouble petsc_gather_ct       = 0.0; /* The number of gathers and gathervs */
 57: PetscLogDouble petsc_scatter_ct      = 0.0; /* The number of scatters and scattervs */

 59: /* Thread Local storage */
 60: PETSC_TLS PetscLogDouble petsc_TotalFlops_th      = 0.0;
 61: PETSC_TLS PetscLogDouble petsc_send_ct_th         = 0.0;
 62: PETSC_TLS PetscLogDouble petsc_recv_ct_th         = 0.0;
 63: PETSC_TLS PetscLogDouble petsc_send_len_th        = 0.0;
 64: PETSC_TLS PetscLogDouble petsc_recv_len_th        = 0.0;
 65: PETSC_TLS PetscLogDouble petsc_isend_ct_th        = 0.0;
 66: PETSC_TLS PetscLogDouble petsc_irecv_ct_th        = 0.0;
 67: PETSC_TLS PetscLogDouble petsc_isend_len_th       = 0.0;
 68: PETSC_TLS PetscLogDouble petsc_irecv_len_th       = 0.0;
 69: PETSC_TLS PetscLogDouble petsc_wait_ct_th         = 0.0;
 70: PETSC_TLS PetscLogDouble petsc_wait_any_ct_th     = 0.0;
 71: PETSC_TLS PetscLogDouble petsc_wait_all_ct_th     = 0.0;
 72: PETSC_TLS PetscLogDouble petsc_sum_of_waits_ct_th = 0.0;
 73: PETSC_TLS PetscLogDouble petsc_allreduce_ct_th    = 0.0;
 74: PETSC_TLS PetscLogDouble petsc_gather_ct_th       = 0.0;
 75: PETSC_TLS PetscLogDouble petsc_scatter_ct_th      = 0.0;

 77: PetscLogDouble petsc_ctog_ct        = 0.0; /* The total number of CPU to GPU copies */
 78: PetscLogDouble petsc_gtoc_ct        = 0.0; /* The total number of GPU to CPU copies */
 79: PetscLogDouble petsc_ctog_sz        = 0.0; /* The total size of CPU to GPU copies */
 80: PetscLogDouble petsc_gtoc_sz        = 0.0; /* The total size of GPU to CPU copies */
 81: PetscLogDouble petsc_ctog_ct_scalar = 0.0; /* The total number of CPU to GPU copies */
 82: PetscLogDouble petsc_gtoc_ct_scalar = 0.0; /* The total number of GPU to CPU copies */
 83: PetscLogDouble petsc_ctog_sz_scalar = 0.0; /* The total size of CPU to GPU copies */
 84: PetscLogDouble petsc_gtoc_sz_scalar = 0.0; /* The total size of GPU to CPU copies */
 85: PetscLogDouble petsc_gflops         = 0.0; /* The flops done on a GPU */
 86: PetscLogDouble petsc_gtime          = 0.0; /* The time spent on a GPU */
 87: PetscLogDouble petsc_genergy        = 0.0; /* The energy (estimated with power*gtime) consumed on a GPU */
 88: PetscLogDouble petsc_genergy_meter  = 0.0; /* Readings from the energy meter on a GPU */

 90: PETSC_TLS PetscLogDouble petsc_ctog_ct_th        = 0.0;
 91: PETSC_TLS PetscLogDouble petsc_gtoc_ct_th        = 0.0;
 92: PETSC_TLS PetscLogDouble petsc_ctog_sz_th        = 0.0;
 93: PETSC_TLS PetscLogDouble petsc_gtoc_sz_th        = 0.0;
 94: PETSC_TLS PetscLogDouble petsc_ctog_ct_scalar_th = 0.0;
 95: PETSC_TLS PetscLogDouble petsc_gtoc_ct_scalar_th = 0.0;
 96: PETSC_TLS PetscLogDouble petsc_ctog_sz_scalar_th = 0.0;
 97: PETSC_TLS PetscLogDouble petsc_gtoc_sz_scalar_th = 0.0;
 98: PETSC_TLS PetscLogDouble petsc_gflops_th         = 0.0;
 99: PETSC_TLS PetscLogDouble petsc_gtime_th          = 0.0;

101: PetscBool PetscLogMemory = PETSC_FALSE;
102: PetscBool PetscLogSyncOn = PETSC_FALSE;

104: PetscBool PetscLogGpuTimeFlag        = PETSC_FALSE;
105: PetscBool PetscLogGpuEnergyFlag      = PETSC_FALSE;
106: PetscBool PetscLogGpuEnergyMeterFlag = PETSC_FALSE;

108: PetscInt PetscLogNumViewersCreated   = 0;
109: PetscInt PetscLogNumViewersDestroyed = 0;

111: PetscLogState petsc_log_state = NULL;

113: #define PETSC_LOG_HANDLER_HOT_BLANK {NULL, NULL, NULL, NULL, NULL, NULL}

115: PetscLogHandlerHot PetscLogHandlers[PETSC_LOG_HANDLER_MAX] = {
116:   PETSC_LOG_HANDLER_HOT_BLANK,
117:   PETSC_LOG_HANDLER_HOT_BLANK,
118:   PETSC_LOG_HANDLER_HOT_BLANK,
119:   PETSC_LOG_HANDLER_HOT_BLANK,
120: };

122: #undef PETSC_LOG_HANDLERS_HOT_BLANK

124: #if defined(PETSC_USE_LOG)
125: #include <../src/sys/logging/handler/impls/default/logdefault.h>

127:   #if defined(PETSC_HAVE_THREADSAFETY)
128: PetscErrorCode PetscAddLogDouble(PetscLogDouble *tot, PetscLogDouble *tot_th, PetscLogDouble tmp)
129: {
130:   *tot_th += tmp;
131:   PetscCall(PetscSpinlockLock(&PetscLogSpinLock));
132:   *tot += tmp;
133:   PetscCall(PetscSpinlockUnlock(&PetscLogSpinLock));
134:   return PETSC_SUCCESS;
135: }

137: PetscErrorCode PetscAddLogDoubleCnt(PetscLogDouble *cnt, PetscLogDouble *tot, PetscLogDouble *cnt_th, PetscLogDouble *tot_th, PetscLogDouble tmp)
138: {
139:   *cnt_th = *cnt_th + 1;
140:   *tot_th += tmp;
141:   PetscCall(PetscSpinlockLock(&PetscLogSpinLock));
142:   *tot += (PetscLogDouble)tmp;
143:   *cnt += *cnt + 1;
144:   PetscCall(PetscSpinlockUnlock(&PetscLogSpinLock));
145:   return PETSC_SUCCESS;
146: }

148:   #endif

150: static PetscErrorCode PetscLogTryGetHandler(PetscLogHandlerType type, PetscLogHandler *handler)
151: {
152:   PetscFunctionBegin;
153:   PetscAssertPointer(handler, 2);
154:   *handler = NULL;
155:   for (int i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
156:     PetscLogHandler h = PetscLogHandlers[i].handler;
157:     if (h) {
158:       PetscBool match;

160:       PetscCall(PetscObjectTypeCompare((PetscObject)h, type, &match));
161:       if (match) {
162:         *handler = PetscLogHandlers[i].handler;
163:         PetscFunctionReturn(PETSC_SUCCESS);
164:       }
165:     }
166:   }
167:   PetscFunctionReturn(PETSC_SUCCESS);
168: }

170: /*@
171:   PetscLogGetDefaultHandler - Get the default log handler if it is running.

173:   Not collective

175:   Output Parameter:
176: . handler - the default `PetscLogHandler`, or `NULL` if it is not running.

178:   Level: developer

180:   Notes:
181:   The default handler is started with `PetscLogDefaultBegin()`,
182:   if the options flags `-log_all` or `-log_view` is given without arguments,
183:   or for `-log_view :output:format` if `format` is not `ascii_xml` or `ascii_flamegraph`.

185: .seealso: [](ch_profiling)
186: @*/
187: PetscErrorCode PetscLogGetDefaultHandler(PetscLogHandler *handler)
188: {
189:   PetscFunctionBegin;
190:   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, handler));
191:   PetscFunctionReturn(PETSC_SUCCESS);
192: }

194: static PetscErrorCode PetscLogGetHandler(PetscLogHandlerType type, PetscLogHandler *handler)
195: {
196:   PetscFunctionBegin;
197:   PetscAssertPointer(handler, 2);
198:   PetscCall(PetscLogTryGetHandler(type, handler));
199:   PetscCheck(*handler != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "A PetscLogHandler of type %s has not been started.", type);
200:   PetscFunctionReturn(PETSC_SUCCESS);
201: }

203: /*@
204:   PetscLogGetState - Get the `PetscLogState` for PETSc's global logging, used
205:   by all default log handlers (`PetscLogDefaultBegin()`,
206:   `PetscLogNestedBegin()`, `PetscLogTraceBegin()`, `PetscLogMPEBegin()`,
207:   `PetscLogPerfstubsBegin()`).

209:   Collective on `PETSC_COMM_WORLD`

211:   Output Parameter:
212: . state - The `PetscLogState` changed by registrations (such as
213:           `PetscLogEventRegister()`) and actions (such as `PetscLogEventBegin()` or
214:           `PetscLogStagePush()`), or `NULL` if logging is not active

216:   Level: developer

218: .seealso: [](ch_profiling), `PetscLogState`
219: @*/
220: PetscErrorCode PetscLogGetState(PetscLogState *state)
221: {
222:   PetscFunctionBegin;
223:   PetscAssertPointer(state, 1);
224:   *state = petsc_log_state;
225:   PetscFunctionReturn(PETSC_SUCCESS);
226: }

228: static PetscErrorCode PetscLogHandlerCopyToHot(PetscLogHandler h, PetscLogHandlerHot *hot)
229: {
230:   PetscFunctionBegin;
231:   hot->handler       = h;
232:   hot->eventBegin    = h->ops->eventbegin;
233:   hot->eventEnd      = h->ops->eventend;
234:   hot->eventSync     = h->ops->eventsync;
235:   hot->objectCreate  = h->ops->objectcreate;
236:   hot->objectDestroy = h->ops->objectdestroy;
237:   PetscFunctionReturn(PETSC_SUCCESS);
238: }

240: /*@
241:   PetscLogHandlerStart - Connect a log handler to PETSc's global logging stream and state.

243:   Logically collective

245:   Input Parameters:
246: . h - a `PetscLogHandler`

248:   Level: developer

250:   Notes:
251:   Users should only need this if they create their own log handlers: handlers that are started
252:   from the command line (such as `-log_view` and `-log_trace`) or from a function like
253:   `PetscLogNestedBegin()` will automatically be started.

255:   There is a limit of `PESC_LOG_HANDLER_MAX` handlers that can be active at one time.

257:   To disconnect a handler from the global stream call `PetscLogHandlerStop()`.

259:   When a log handler is started, stages that have already been pushed with `PetscLogStagePush()`,
260:   will be pushed for the new log handler, but it will not be informed of any events that are
261:   in progress.  It is recommended to start any user-defined log handlers immediately following
262:   `PetscInitialize()`  before any user-defined stages are pushed.

264: .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogState`, `PetscLogHandlerStop()`, `PetscInitialize()`
265: @*/
266: PetscErrorCode PetscLogHandlerStart(PetscLogHandler h)
267: {
268:   PetscFunctionBegin;
269:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
270:     if (PetscLogHandlers[i].handler == h) PetscFunctionReturn(PETSC_SUCCESS);
271:   }
272:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
273:     if (PetscLogHandlers[i].handler == NULL) {
274:       PetscCall(PetscObjectReference((PetscObject)h));
275:       PetscCall(PetscLogHandlerCopyToHot(h, &PetscLogHandlers[i]));
276:       if (petsc_log_state) {
277:         PetscLogStage stack_height;
278:         PetscIntStack orig_stack, temp_stack;

280:         PetscCall(PetscLogHandlerSetState(h, petsc_log_state));
281:         stack_height = petsc_log_state->stage_stack->top + 1;
282:         PetscCall(PetscIntStackCreate(&temp_stack));
283:         orig_stack                     = petsc_log_state->stage_stack;
284:         petsc_log_state->stage_stack   = temp_stack;
285:         petsc_log_state->current_stage = -1;
286:         for (int s = 0; s < stack_height; s++) {
287:           PetscLogStage stage = orig_stack->stack[s];
288:           PetscCall(PetscLogHandlerStagePush(h, stage));
289:           PetscCall(PetscIntStackPush(temp_stack, stage));
290:           petsc_log_state->current_stage = stage;
291:         }
292:         PetscCall(PetscIntStackDestroy(temp_stack));
293:         petsc_log_state->stage_stack = orig_stack;
294:       }
295:       PetscFunctionReturn(PETSC_SUCCESS);
296:     }
297:   }
298:   SETERRQ(PetscObjectComm((PetscObject)h), PETSC_ERR_ARG_WRONGSTATE, "%d log handlers already started, cannot start another", PETSC_LOG_HANDLER_MAX);
299:   PetscFunctionReturn(PETSC_SUCCESS);
300: }

302: /*@
303:   PetscLogHandlerStop - Disconnect a log handler from PETSc's global logging stream.

305:   Logically collective

307:   Input Parameters:
308: . h - a `PetscLogHandler`

310:   Level: developer

312:   Note:
313:   After `PetscLogHandlerStop()`, the handler can still access the global logging state
314:   with `PetscLogHandlerGetState()`, so that it can access the registry when post-processing
315:   (for instance, in `PetscLogHandlerView()`),

317:   When a log handler is stopped, the remaining stages will be popped before it is
318:   disconnected from the log stream.

320: .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogState`, `PetscLogHandlerStart()`
321: @*/
322: PetscErrorCode PetscLogHandlerStop(PetscLogHandler h)
323: {
324:   PetscFunctionBegin;
325:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
326:     if (PetscLogHandlers[i].handler == h) {
327:       if (petsc_log_state) {
328:         PetscLogState state;
329:         PetscLogStage stack_height;
330:         PetscIntStack orig_stack, temp_stack;

332:         PetscCall(PetscLogHandlerGetState(h, &state));
333:         PetscCheck(state == petsc_log_state, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONGSTATE, "Called PetscLogHandlerStop() for a PetscLogHander that was not started.");
334:         stack_height = petsc_log_state->stage_stack->top + 1;
335:         PetscCall(PetscIntStackCreate(&temp_stack));
336:         orig_stack                   = petsc_log_state->stage_stack;
337:         petsc_log_state->stage_stack = temp_stack;
338:         for (int s = 0; s < stack_height; s++) {
339:           PetscLogStage stage = orig_stack->stack[s];

341:           PetscCall(PetscIntStackPush(temp_stack, stage));
342:         }
343:         for (int s = 0; s < stack_height; s++) {
344:           PetscLogStage stage;
345:           PetscBool     empty;

347:           PetscCall(PetscIntStackPop(temp_stack, &stage));
348:           PetscCall(PetscIntStackEmpty(temp_stack, &empty));
349:           if (!empty) PetscCall(PetscIntStackTop(temp_stack, &petsc_log_state->current_stage));
350:           else petsc_log_state->current_stage = -1;
351:           PetscCall(PetscLogHandlerStagePop(h, stage));
352:         }
353:         PetscCall(PetscIntStackDestroy(temp_stack));
354:         petsc_log_state->stage_stack = orig_stack;
355:         PetscCall(PetscIntStackTop(petsc_log_state->stage_stack, &petsc_log_state->current_stage));
356:       }
357:       PetscCall(PetscArrayzero(&PetscLogHandlers[i], 1));
358:       PetscCall(PetscObjectDereference((PetscObject)h));
359:     }
360:   }
361:   PetscFunctionReturn(PETSC_SUCCESS);
362: }

364: /*@
365:   PetscLogIsActive - Check if logging (profiling) is currently in progress.

367:   Not Collective

369:   Output Parameter:
370: . isActive - `PETSC_TRUE` if logging is in progress, `PETSC_FALSE` otherwise

372:   Level: beginner

374: .seealso: [](ch_profiling), `PetscLogDefaultBegin()`
375: @*/
376: PetscErrorCode PetscLogIsActive(PetscBool *isActive)
377: {
378:   PetscFunctionBegin;
379:   *isActive = PETSC_FALSE;
380:   if (petsc_log_state) {
381:     for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
382:       if (PetscLogHandlers[i].handler) {
383:         *isActive = PETSC_TRUE;
384:         PetscFunctionReturn(PETSC_SUCCESS);
385:       }
386:     }
387:   }
388:   PetscFunctionReturn(PETSC_SUCCESS);
389: }

391: PETSC_UNUSED static PetscErrorCode PetscLogEventBeginIsActive(PetscBool *isActive)
392: {
393:   PetscFunctionBegin;
394:   *isActive = PETSC_FALSE;
395:   if (petsc_log_state) {
396:     for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
397:       if (PetscLogHandlers[i].eventBegin) {
398:         *isActive = PETSC_TRUE;
399:         PetscFunctionReturn(PETSC_SUCCESS);
400:       }
401:     }
402:   }
403:   PetscFunctionReturn(PETSC_SUCCESS);
404: }

406: PETSC_UNUSED static PetscErrorCode PetscLogEventEndIsActive(PetscBool *isActive)
407: {
408:   PetscFunctionBegin;
409:   *isActive = PETSC_FALSE;
410:   if (petsc_log_state) {
411:     for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
412:       if (PetscLogHandlers[i].eventEnd) {
413:         *isActive = PETSC_TRUE;
414:         PetscFunctionReturn(PETSC_SUCCESS);
415:       }
416:     }
417:   }
418:   PetscFunctionReturn(PETSC_SUCCESS);
419: }

421: PETSC_INTERN PetscErrorCode PetscLogTypeBegin(PetscLogHandlerType type)
422: {
423:   PetscLogHandler handler;

425:   PetscFunctionBegin;
426:   PetscCall(PetscLogTryGetHandler(type, &handler));
427:   if (handler) PetscFunctionReturn(PETSC_SUCCESS);
428:   PetscCall(PetscLogHandlerCreate(PETSC_COMM_WORLD, &handler));
429:   PetscCall(PetscLogHandlerSetType(handler, type));
430:   PetscCall(PetscLogHandlerStart(handler));
431:   PetscCall(PetscLogHandlerDestroy(&handler));
432:   PetscFunctionReturn(PETSC_SUCCESS);
433: }

435: /*@
436:   PetscLogDefaultBegin - Turns on logging (profiling) of PETSc code using the default log handler (profiler). This logs time, flop
437:   rates, and object creation and should not slow programs down too much.

439:   Logically Collective on `PETSC_COMM_WORLD`

441:   Options Database Key:
442: . -log_view [viewertype:filename:viewerformat] - Prints summary of flop and timing (profiling) information to the
443:                                                  screen (for PETSc configured with `--with-log=1` (which is the default)).
444:                                                  This option must be provided before `PetscInitialize()`.

446:   Example Usage:
447: .vb
448:       PetscInitialize(...);
449:       PetscLogDefaultBegin();
450:        ... code ...
451:       PetscLogView(viewer); or PetscLogDump();
452:       PetscFinalize();
453: .ve

455:   Level: advanced

457:   Notes:
458:   `PetscLogView()` or `PetscLogDump()` actually cause the printing of
459:   the logging information.

461:   This routine may be called more than once.

463:   To provide the `-log_view` option in your source code you must call  PetscCall(PetscOptionsSetValue(NULL, "-log_view", NULL));
464:   before you call `PetscInitialize()`

466: .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogView()`, `PetscLogTraceBegin()`
467: @*/
468: PetscErrorCode PetscLogDefaultBegin(void)
469: {
470:   PetscFunctionBegin;
471:   PetscCall(PetscLogTypeBegin(PETSCLOGHANDLERDEFAULT));
472:   PetscFunctionReturn(PETSC_SUCCESS);
473: }

475: /*@C
476:   PetscLogTraceBegin - Begins trace logging.  Every time a PETSc event
477:   begins or ends, the event name is printed.

479:   Logically Collective on `PETSC_COMM_WORLD`, No Fortran Support

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

484:   Options Database Key:
485: . -log_trace [filename] - Begins `PetscLogTraceBegin()`

487:   Level: intermediate

489:   Notes:
490:   `PetscLogTraceBegin()` prints the processor number, the execution time (sec),
491:   then "Event begin:" or "Event end:" followed by the event name.

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

497: .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogView()`, `PetscLogDefaultBegin()`
498: @*/
499: PetscErrorCode PetscLogTraceBegin(FILE *file)
500: {
501:   PetscLogHandler handler;

503:   PetscFunctionBegin;
504:   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERTRACE, &handler));
505:   if (handler) PetscFunctionReturn(PETSC_SUCCESS);
506:   PetscCall(PetscLogHandlerCreateTrace(PETSC_COMM_WORLD, file, &handler));
507:   PetscCall(PetscLogHandlerStart(handler));
508:   PetscCall(PetscLogHandlerDestroy(&handler));
509:   PetscFunctionReturn(PETSC_SUCCESS);
510: }

512: PETSC_INTERN PetscErrorCode PetscLogHandlerCreate_Nested(MPI_Comm, PetscLogHandler *);

514: /*@
515:   PetscLogNestedBegin - Turns on nested logging of objects and events. This logs flop
516:   rates and object creation and should not slow programs down too much.

518:   Logically Collective on `PETSC_COMM_WORLD`, No Fortran Support

520:   Options Database Keys:
521: . -log_view :filename.xml:ascii_xml - Prints an XML summary of flop and timing information to the file

523:   Example Usage:
524: .vb
525:       PetscInitialize(...);
526:       PetscLogNestedBegin();
527:        ... code ...
528:       PetscLogView(viewer);
529:       PetscFinalize();
530: .ve

532:   Level: advanced

534: .seealso: `PetscLogDump()`, `PetscLogView()`, `PetscLogTraceBegin()`, `PetscLogDefaultBegin()`
535: @*/
536: PetscErrorCode PetscLogNestedBegin(void)
537: {
538:   PetscFunctionBegin;
539:   PetscCall(PetscLogTypeBegin(PETSCLOGHANDLERNESTED));
540:   PetscFunctionReturn(PETSC_SUCCESS);
541: }

543: /*@C
544:   PetscLogLegacyCallbacksBegin - Create and start a log handler from callbacks
545:   matching the now deprecated function pointers `PetscLogPLB`, `PetscLogPLE`,
546:   `PetscLogPHC`, `PetscLogPHD`.

548:   Logically Collective on `PETSC_COMM_WORLD`

550:   Input Parameters:
551: + PetscLogPLB - A callback that will be executed by `PetscLogEventBegin()` (or `NULL`)
552: . PetscLogPLE - A callback that will be executed by `PetscLogEventEnd()` (or `NULL`)
553: . PetscLogPHC - A callback that will be executed by `PetscLogObjectCreate()` (or `NULL`)
554: - PetscLogPHD - A callback that will be executed by `PetscLogObjectCreate()` (or `NULL`)

556:   Calling sequence of `PetscLogPLB`:
557: + e  - a `PetscLogEvent` that is beginning
558: . _i - deprecated, unused
559: . o1 - a `PetscObject` associated with `e` (or `NULL`)
560: . o2 - a `PetscObject` associated with `e` (or `NULL`)
561: . o3 - a `PetscObject` associated with `e` (or `NULL`)
562: - o4 - a `PetscObject` associated with `e` (or `NULL`)

564:   Calling sequence of `PetscLogPLE`:
565: + e  - a `PetscLogEvent` that is beginning
566: . _i - deprecated, unused
567: . o1 - a `PetscObject` associated with `e` (or `NULL`)
568: . o2 - a `PetscObject` associated with `e` (or `NULL`)
569: . o3 - a `PetscObject` associated with `e` (or `NULL`)
570: - o4 - a `PetscObject` associated with `e` (or `NULL`)

572:   Calling sequence of `PetscLogPHC`:
573: . o - a `PetscObject` that has just been created

575:   Calling sequence of `PetscLogPHD`:
576: . o - a `PetscObject` that is about to be destroyed

578:   Level: advanced

580:   Notes:
581:   This is for transitioning from the deprecated function `PetscLogSet()` and should not be used in new code.

583:   This should help migrate external log handlers to use `PetscLogHandler`, but
584:   callbacks that depend on the deprecated `PetscLogStage` datatype will have to be
585:   updated.

587: .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogHandlerStart()`, `PetscLogState`
588: @*/
589: PetscErrorCode PetscLogLegacyCallbacksBegin(PetscErrorCode (*PetscLogPLB)(PetscLogEvent e, int _i, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4), PetscErrorCode (*PetscLogPLE)(PetscLogEvent e, int _i, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4), PetscErrorCode (*PetscLogPHC)(PetscObject o), PetscErrorCode (*PetscLogPHD)(PetscObject o))
590: {
591:   PetscLogHandler handler;

593:   PetscFunctionBegin;
594:   PetscCall(PetscLogHandlerCreateLegacy(PETSC_COMM_WORLD, PetscLogPLB, PetscLogPLE, PetscLogPHC, PetscLogPHD, &handler));
595:   PetscCall(PetscLogHandlerStart(handler));
596:   PetscCall(PetscLogHandlerDestroy(&handler));
597:   PetscFunctionReturn(PETSC_SUCCESS);
598: }

600:   #if defined(PETSC_HAVE_MPE)
601:     #include <mpe.h>
602: static PetscBool PetscBeganMPE = PETSC_FALSE;
603:   #endif

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

609:   Collective on `PETSC_COMM_WORLD`, No Fortran Support

611:   Options Database Key:
612: . -log_mpe - Prints extensive log information

614:   Level: advanced

616:   Note:
617:   A related routine is `PetscLogDefaultBegin()` (with the options key `-log_view`), which is
618:   intended for production runs since it logs only flop rates and object creation (and should
619:   not significantly slow the programs).

621: .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogDefaultBegin()`, `PetscLogEventActivate()`,
622:           `PetscLogEventDeactivate()`
623: @*/
624: PetscErrorCode PetscLogMPEBegin(void)
625: {
626:   PetscFunctionBegin;
627:   #if defined(PETSC_HAVE_MPE)
628:   /* Do MPE initialization */
629:   if (!MPE_Initialized_logging()) { /* This function exists in mpich 1.1.2 and higher */
630:     PetscCall(PetscInfo(0, "Initializing MPE.\n"));
631:     PetscCall(MPE_Init_log());

633:     PetscBeganMPE = PETSC_TRUE;
634:   } else {
635:     PetscCall(PetscInfo(0, "MPE already initialized. Not attempting to reinitialize.\n"));
636:   }
637:   PetscCall(PetscLogTypeBegin(PETSCLOGHANDLERMPE));
638:   #else
639:   SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "PETSc was configured without MPE support, reconfigure with --with-mpe or --download-mpe");
640:   #endif
641:   PetscFunctionReturn(PETSC_SUCCESS);
642: }

644:   #if defined(PETSC_HAVE_TAU_PERFSTUBS)
645: #include <../src/sys/perfstubs/timer.h>
646:   #endif

648: /*@C
649:   PetscLogPerfstubsBegin - Turns on logging of events using the perfstubs interface.

651:   Collective on `PETSC_COMM_WORLD`, No Fortran Support

653:   Options Database Key:
654: . -log_perfstubs - use an external log handler through the perfstubs interface

656:   Level: advanced

658: .seealso: [](ch_profiling), `PetscLogDefaultBegin()`, `PetscLogEventActivate()`
659: @*/
660: PetscErrorCode PetscLogPerfstubsBegin(void)
661: {
662:   PetscFunctionBegin;
663:   #if defined(PETSC_HAVE_TAU_PERFSTUBS)
664:   PetscCall(PetscLogTypeBegin(PETSCLOGHANDLERPERFSTUBS));
665:   #else
666:   SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "PETSc was configured without perfstubs support, reconfigure with --with-tau-perfstubs");
667:   #endif
668:   PetscFunctionReturn(PETSC_SUCCESS);
669: }

671: /*@
672:   PetscLogActions - Determines whether actions are logged for the default log handler.

674:   Not Collective

676:   Input Parameter:
677: . flag - `PETSC_TRUE` if actions are to be logged

679:   Options Database Key:
680: + -log_exclude_actions - (deprecated) Does nothing
681: - -log_include_actions - Turn on action logging

683:   Level: intermediate

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

689: .seealso: [](ch_profiling), `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogGetDefaultHandler()`
690: @*/
691: PetscErrorCode PetscLogActions(PetscBool flag)
692: {
693:   PetscFunctionBegin;
694:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
695:     PetscLogHandler h = PetscLogHandlers[i].handler;

697:     if (h) PetscCall(PetscLogHandlerSetLogActions(h, flag));
698:   }
699:   PetscFunctionReturn(PETSC_SUCCESS);
700: }

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

705:   Not Collective

707:   Input Parameter:
708: . flag - `PETSC_TRUE` if objects are to be logged

710:   Options Database Key:
711: + -log_exclude_objects - (deprecated) Does nothing
712: - -log_include_objects - Turns on object logging

714:   Level: intermediate

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

720: .seealso: [](ch_profiling), `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogGetDefaultHandler()`
721: @*/
722: PetscErrorCode PetscLogObjects(PetscBool flag)
723: {
724:   PetscFunctionBegin;
725:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
726:     PetscLogHandler h = PetscLogHandlers[i].handler;

728:     if (h) PetscCall(PetscLogHandlerSetLogObjects(h, flag));
729:   }
730:   PetscFunctionReturn(PETSC_SUCCESS);
731: }

733: /*------------------------------------------------ Stage Functions --------------------------------------------------*/
734: /*@
735:   PetscLogStageRegister - Attaches a character string name to a logging stage.

737:   Not Collective

739:   Input Parameter:
740: . sname - The name to associate with that stage

742:   Output Parameter:
743: . stage - The stage number or -1 if logging is not active (`PetscLogIsActive()`).

745:   Level: intermediate

747: .seealso: [](ch_profiling), `PetscLogStagePush()`, `PetscLogStagePop()`
748: @*/
749: PetscErrorCode PetscLogStageRegister(const char sname[], PetscLogStage *stage)
750: {
751:   PetscLogState state;

753:   PetscFunctionBegin;
754:   *stage = -1;
755:   PetscCall(PetscLogGetState(&state));
756:   if (state) PetscCall(PetscLogStateStageRegister(state, sname, stage));
757:   PetscFunctionReturn(PETSC_SUCCESS);
758: }

760: /*@
761:   PetscLogStagePush - This function pushes a stage on the logging stack. Events started and stopped until `PetscLogStagePop()` will be associated with the stage

763:   Not Collective

765:   Input Parameter:
766: . stage - The stage on which to log

768:   Example Usage:
769:   If the option -log_view is used to run the program containing the
770:   following code, then 2 sets of summary data will be printed during
771:   PetscFinalize().
772: .vb
773:       PetscInitialize(int *argc,char ***args,0,0);
774:       [stage 0 of code]
775:       PetscLogStagePush(1);
776:       [stage 1 of code]
777:       PetscLogStagePop();
778:       PetscBarrier(...);
779:       [more stage 0 of code]
780:       PetscFinalize();
781: .ve

783:   Level: intermediate

785:   Note:
786:   Use `PetscLogStageRegister()` to register a stage.

788: .seealso: [](ch_profiling), `PetscLogStagePop()`, `PetscLogStageRegister()`, `PetscBarrier()`
789: @*/
790: PetscErrorCode PetscLogStagePush(PetscLogStage stage)
791: {
792:   PetscLogState state;

794:   PetscFunctionBegin;
795:   PetscCall(PetscLogGetState(&state));
796:   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
797:   for (int i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
798:     PetscLogHandler h = PetscLogHandlers[i].handler;
799:     if (h) PetscCall(PetscLogHandlerStagePush(h, stage));
800:   }
801:   PetscCall(PetscLogStateStagePush(state, stage));
802:   PetscFunctionReturn(PETSC_SUCCESS);
803: }

805: /*@
806:   PetscLogStagePop - This function pops a stage from the logging stack that was pushed with `PetscLogStagePush()`

808:   Not Collective

810:   Example Usage:
811:   If the option -log_view is used to run the program containing the
812:   following code, then 2 sets of summary data will be printed during
813:   PetscFinalize().
814: .vb
815:       PetscInitialize(int *argc,char ***args,0,0);
816:       [stage 0 of code]
817:       PetscLogStagePush(1);
818:       [stage 1 of code]
819:       PetscLogStagePop();
820:       PetscBarrier(...);
821:       [more stage 0 of code]
822:       PetscFinalize();
823: .ve

825:   Level: intermediate

827: .seealso: [](ch_profiling), `PetscLogStagePush()`, `PetscLogStageRegister()`, `PetscBarrier()`
828: @*/
829: PetscErrorCode PetscLogStagePop(void)
830: {
831:   PetscLogState state;
832:   PetscLogStage current_stage;

834:   PetscFunctionBegin;
835:   PetscCall(PetscLogGetState(&state));
836:   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
837:   current_stage = state->current_stage;
838:   PetscCall(PetscLogStateStagePop(state));
839:   for (int i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
840:     PetscLogHandler h = PetscLogHandlers[i].handler;
841:     if (h) PetscCall(PetscLogHandlerStagePop(h, current_stage));
842:   }
843:   PetscFunctionReturn(PETSC_SUCCESS);
844: }

846: /*@
847:   PetscLogStageSetActive - Sets if a stage is used for `PetscLogEventBegin()` and `PetscLogEventEnd()`.

849:   Not Collective

851:   Input Parameters:
852: + stage    - The stage
853: - isActive - The activity flag, `PETSC_TRUE` for logging, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)

855:   Level: intermediate

857:   Note:
858:   If this is set to `PETSC_FALSE` the logging acts as if the stage did not exist

860: .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
861: @*/
862: PetscErrorCode PetscLogStageSetActive(PetscLogStage stage, PetscBool isActive)
863: {
864:   PetscLogState state;

866:   PetscFunctionBegin;
867:   PetscCall(PetscLogGetState(&state));
868:   if (state) PetscCall(PetscLogStateStageSetActive(state, stage, isActive));
869:   PetscFunctionReturn(PETSC_SUCCESS);
870: }

872: /*@
873:   PetscLogStageGetActive - Checks if a stage is used for `PetscLogEventBegin()` and `PetscLogEventEnd()`.

875:   Not Collective

877:   Input Parameter:
878: . stage - The stage

880:   Output Parameter:
881: . isActive - The activity flag, `PETSC_TRUE` for logging, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)

883:   Level: intermediate

885: .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
886: @*/
887: PetscErrorCode PetscLogStageGetActive(PetscLogStage stage, PetscBool *isActive)
888: {
889:   PetscLogState state;

891:   PetscFunctionBegin;
892:   *isActive = PETSC_FALSE;
893:   PetscCall(PetscLogGetState(&state));
894:   if (state) PetscCall(PetscLogStateStageGetActive(state, stage, isActive));
895:   PetscFunctionReturn(PETSC_SUCCESS);
896: }

898: /*@
899:   PetscLogStageSetVisible - Determines stage visibility in `PetscLogView()`

901:   Not Collective

903:   Input Parameters:
904: + stage     - The stage
905: - isVisible - The visibility flag, `PETSC_TRUE` to print, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)

907:   Level: intermediate

909:   Developer Notes:
910:   Visibility only affects the default log handler in `PetscLogView()`: stages that are
911:   set to invisible are suppressed from output.

913: .seealso: [](ch_profiling), `PetscLogStageGetVisible()`, `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogView()`, `PetscLogGetDefaultHandler()`
914: @*/
915: PetscErrorCode PetscLogStageSetVisible(PetscLogStage stage, PetscBool isVisible)
916: {
917:   PetscFunctionBegin;
918:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
919:     PetscLogHandler h = PetscLogHandlers[i].handler;

921:     if (h) PetscCall(PetscLogHandlerStageSetVisible(h, stage, isVisible));
922:   }
923:   PetscFunctionReturn(PETSC_SUCCESS);
924: }

926: /*@
927:   PetscLogStageGetVisible - Returns stage visibility in `PetscLogView()`

929:   Not Collective

931:   Input Parameter:
932: . stage - The stage

934:   Output Parameter:
935: . isVisible - The visibility flag, `PETSC_TRUE` to print, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)

937:   Level: intermediate

939: .seealso: [](ch_profiling), `PetscLogStageSetVisible()`, `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogView()`, `PetscLogGetDefaultHandler()`
940: @*/
941: PetscErrorCode PetscLogStageGetVisible(PetscLogStage stage, PetscBool *isVisible)
942: {
943:   PetscLogHandler handler;

945:   PetscFunctionBegin;
946:   *isVisible = PETSC_FALSE;
947:   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
948:   if (handler) PetscCall(PetscLogHandlerStageGetVisible(handler, stage, isVisible));
949:   PetscFunctionReturn(PETSC_SUCCESS);
950: }

952: /*@
953:   PetscLogStageGetId - Returns the stage id when given the stage name.

955:   Not Collective

957:   Input Parameter:
958: . name - The stage name

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

963:   Level: intermediate

965: .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
966: @*/
967: PetscErrorCode PetscLogStageGetId(const char name[], PetscLogStage *stage)
968: {
969:   PetscLogState state;

971:   PetscFunctionBegin;
972:   *stage = -1;
973:   PetscCall(PetscLogGetState(&state));
974:   if (state) PetscCall(PetscLogStateGetStageFromName(state, name, stage));
975:   PetscFunctionReturn(PETSC_SUCCESS);
976: }

978: /*@
979:   PetscLogStageGetName - Returns the stage name when given the stage id.

981:   Not Collective

983:   Input Parameter:
984: . stage - The stage

986:   Output Parameter:
987: . name - The stage name

989:   Level: intermediate

991: .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
992: @*/
993: PetscErrorCode PetscLogStageGetName(PetscLogStage stage, const char *name[])
994: {
995:   PetscLogStageInfo stage_info;
996:   PetscLogState     state;

998:   PetscFunctionBegin;
999:   *name = NULL;
1000:   PetscCall(PetscLogGetState(&state));
1001:   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
1002:   PetscCall(PetscLogStateStageGetInfo(state, stage, &stage_info));
1003:   *name = stage_info.name;
1004:   PetscFunctionReturn(PETSC_SUCCESS);
1005: }

1007: /*------------------------------------------------ Event Functions --------------------------------------------------*/

1009: /*@
1010:   PetscLogEventRegister - Registers an event name for logging operations

1012:   Not Collective

1014:   Input Parameters:
1015: + name    - The name associated with the event
1016: - classid - The classid associated to the class for this event, obtain either with
1017:            `PetscClassIdRegister()` or use a predefined one such as `KSP_CLASSID`, `SNES_CLASSID`, the predefined ones
1018:            are only available in C code

1020:   Output Parameter:
1021: . event - The event id for use with `PetscLogEventBegin()` and `PetscLogEventEnd()`.

1023:   Example Usage:
1024: .vb
1025:       PetscLogEvent USER_EVENT;
1026:       PetscClassId classid;
1027:       PetscLogDouble user_event_flops;
1028:       PetscClassIdRegister("class name",&classid);
1029:       PetscLogEventRegister("User event name",classid,&USER_EVENT);
1030:       PetscLogEventBegin(USER_EVENT,0,0,0,0);
1031:          [code segment to monitor]
1032:          PetscLogFlops(user_event_flops);
1033:       PetscLogEventEnd(USER_EVENT,0,0,0,0);
1034: .ve

1036:   Level: intermediate

1038:   Notes:
1039:   PETSc automatically logs library events if the code has been
1040:   configured with --with-log (which is the default) and
1041:   -log_view or -log_all is specified.  `PetscLogEventRegister()` is
1042:   intended for logging user events to supplement this PETSc
1043:   information.

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

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

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

1060: .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogFlops()`,
1061:           `PetscLogEventActivate()`, `PetscLogEventDeactivate()`, `PetscClassIdRegister()`
1062: @*/
1063: PetscErrorCode PetscLogEventRegister(const char name[], PetscClassId classid, PetscLogEvent *event)
1064: {
1065:   PetscLogState state;

1067:   PetscFunctionBegin;
1068:   *event = -1;
1069:   PetscCall(PetscLogGetState(&state));
1070:   if (state) PetscCall(PetscLogStateEventRegister(state, name, classid, event));
1071:   PetscFunctionReturn(PETSC_SUCCESS);
1072: }

1074: /*@
1075:   PetscLogEventSetCollective - Indicates that a particular event is collective.

1077:   Logically Collective

1079:   Input Parameters:
1080: + event      - The event id
1081: - collective - `PetscBool` indicating whether a particular event is collective

1083:   Level: developer

1085:   Notes:
1086:   New events returned from `PetscLogEventRegister()` are collective by default.

1088:   Collective events are handled specially if the command line option `-log_sync` is used. In that case the logging saves information about
1089:   two parts of the event; the time for all the MPI ranks to synchronize and then the time for the actual computation/communication
1090:   to be performed. This option is useful to debug imbalance within the computations or communications.

1092: .seealso: [](ch_profiling), `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogEventRegister()`
1093: @*/
1094: PetscErrorCode PetscLogEventSetCollective(PetscLogEvent event, PetscBool collective)
1095: {
1096:   PetscLogState state;

1098:   PetscFunctionBegin;
1099:   PetscCall(PetscLogGetState(&state));
1100:   if (state) PetscCall(PetscLogStateEventSetCollective(state, event, collective));
1101:   PetscFunctionReturn(PETSC_SUCCESS);
1102: }

1104: /*
1105:   PetscLogClassSetActiveAll - Activate or inactivate logging for all events associated with a PETSc object class in every stage.

1107:   Not Collective

1109:   Input Parameters:
1110: + classid - The object class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1111: - isActive - if `PETSC_FALSE`, events associated with this class will not be send to log handlers.

1113:   Level: developer

1115: .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventActivateAll()`, `PetscLogStageSetActive()`, `PetscLogEventActivateClass()`
1116: */
1117: static PetscErrorCode PetscLogClassSetActiveAll(PetscClassId classid, PetscBool isActive)
1118: {
1119:   PetscLogState state;

1121:   PetscFunctionBegin;
1122:   PetscCall(PetscLogGetState(&state));
1123:   if (state) PetscCall(PetscLogStateClassSetActiveAll(state, classid, isActive));
1124:   PetscFunctionReturn(PETSC_SUCCESS);
1125: }

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

1130:   Not Collective

1132:   Input Parameter:
1133: . classid - The object class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.

1135:   Level: developer

1137: .seealso: [](ch_profiling), `PetscLogEventActivateClass()`, `PetscLogEventDeactivateClass()`, `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1138: @*/
1139: PetscErrorCode PetscLogEventIncludeClass(PetscClassId classid)
1140: {
1141:   PetscFunctionBegin;
1142:   PetscCall(PetscLogClassSetActiveAll(classid, PETSC_TRUE));
1143:   PetscFunctionReturn(PETSC_SUCCESS);
1144: }

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

1149:   Not Collective

1151:   Input Parameter:
1152: . classid - The object class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.

1154:   Level: developer

1156:   Note:
1157:   If a class is excluded then events associated with that class are not logged.

1159: .seealso: [](ch_profiling), `PetscLogEventDeactivateClass()`, `PetscLogEventActivateClass()`, `PetscLogEventDeactivate()`, `PetscLogEventActivate()`
1160: @*/
1161: PetscErrorCode PetscLogEventExcludeClass(PetscClassId classid)
1162: {
1163:   PetscFunctionBegin;
1164:   PetscCall(PetscLogClassSetActiveAll(classid, PETSC_FALSE));
1165:   PetscFunctionReturn(PETSC_SUCCESS);
1166: }

1168: /*
1169:   PetscLogEventSetActive - Activate or inactivate logging for an event in a given stage

1171:   Not Collective

1173:   Input Parameters:
1174: + stage - A registered `PetscLogStage` (or `PETSC_DEFAULT` for the current stage)
1175: . event - A `PetscLogEvent`
1176: - isActive - If `PETSC_FALSE`, activity from this event (`PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogEventSync()`) will not be sent to log handlers during this stage

1178:   Usage:
1179: .vb
1180:       PetscLogEventSetActive(VEC_SetValues, PETSC_FALSE);
1181:         [code where you do not want to log VecSetValues()]
1182:       PetscLogEventSetActive(VEC_SetValues, PETSC_TRUE);
1183:         [code where you do want to log VecSetValues()]
1184: .ve

1186:   Level: advanced

1188:   Note:
1189:   The event may be either a pre-defined PETSc event (found in include/petsclog.h)
1190:   or an event number obtained with `PetscLogEventRegister()`.

1192: .seealso: [](ch_profiling), `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`
1193: */
1194: static PetscErrorCode PetscLogEventSetActive(PetscLogStage stage, PetscLogEvent event, PetscBool isActive)
1195: {
1196:   PetscLogState state;

1198:   PetscFunctionBegin;
1199:   PetscCall(PetscLogGetState(&state));
1200:   if (state) PetscCall(PetscLogStateEventSetActive(state, stage, event, isActive));
1201:   PetscFunctionReturn(PETSC_SUCCESS);
1202: }

1204: /*@
1205:   PetscLogEventActivate - Indicates that a particular event should be logged.

1207:   Not Collective

1209:   Input Parameter:
1210: . event - The event id

1212:   Example Usage:
1213: .vb
1214:       PetscLogEventDeactivate(VEC_SetValues);
1215:         [code where you do not want to log VecSetValues()]
1216:       PetscLogEventActivate(VEC_SetValues);
1217:         [code where you do want to log VecSetValues()]
1218: .ve

1220:   Level: advanced

1222:   Note:
1223:   The event may be either a pre-defined PETSc event (found in include/petsclog.h)
1224:   or an event number obtained with `PetscLogEventRegister()`.

1226: .seealso: [](ch_profiling), `PetscLogEventDeactivate()`, `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`
1227: @*/
1228: PetscErrorCode PetscLogEventActivate(PetscLogEvent event)
1229: {
1230:   PetscFunctionBegin;
1231:   PetscCall(PetscLogEventSetActive(PETSC_DEFAULT, event, PETSC_TRUE));
1232:   PetscFunctionReturn(PETSC_SUCCESS);
1233: }

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

1238:   Not Collective

1240:   Input Parameter:
1241: . event - The event id

1243:   Example Usage:
1244: .vb
1245:       PetscLogEventDeactivate(VEC_SetValues);
1246:         [code where you do not want to log VecSetValues()]
1247:       PetscLogEventActivate(VEC_SetValues);
1248:         [code where you do want to log VecSetValues()]
1249: .ve

1251:   Level: advanced

1253:   Note:
1254:   The event may be either a pre-defined PETSc event (found in
1255:   include/petsclog.h) or an event number obtained with `PetscLogEventRegister()`).

1257: .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`
1258: @*/
1259: PetscErrorCode PetscLogEventDeactivate(PetscLogEvent event)
1260: {
1261:   PetscFunctionBegin;
1262:   PetscCall(PetscLogEventSetActive(PETSC_DEFAULT, event, PETSC_FALSE));
1263:   PetscFunctionReturn(PETSC_SUCCESS);
1264: }

1266: /*@
1267:   PetscLogEventDeactivatePush - Indicates that a particular event should not be logged until `PetscLogEventDeactivatePop()` is called

1269:   Not Collective

1271:   Input Parameter:
1272: . event - The event id

1274:   Example Usage:
1275: .vb
1276:       PetscLogEventDeactivatePush(VEC_SetValues);
1277:         [code where you do not want to log VecSetValues()]
1278:       PetscLogEventDeactivatePop(VEC_SetValues);
1279:         [code where you do want to log VecSetValues()]
1280: .ve

1282:   Level: advanced

1284:   Note:
1285:   The event may be either a pre-defined PETSc event (found in
1286:   include/petsclog.h) or an event number obtained with `PetscLogEventRegister()`).

1288:   PETSc's default log handler (`PetscLogDefaultBegin()`) respects this function because it can make the output of `PetscLogView()` easier to interpret, but other handlers (such as the nested handler, `PetscLogNestedBegin()`) ignore it because suppressing events is not helpful in their output formats.

1290: .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivate()`, `PetscLogEventDeactivatePop()`
1291: @*/
1292: PetscErrorCode PetscLogEventDeactivatePush(PetscLogEvent event)
1293: {
1294:   PetscFunctionBegin;
1295:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1296:     PetscLogHandler h = PetscLogHandlers[i].handler;

1298:     if (h) PetscCall(PetscLogHandlerEventDeactivatePush(h, PETSC_DEFAULT, event));
1299:   }
1300:   PetscFunctionReturn(PETSC_SUCCESS);
1301: }

1303: /*@
1304:   PetscLogEventDeactivatePop - Indicates that a particular event should again be logged after the logging was turned off with `PetscLogEventDeactivatePush()`

1306:   Not Collective

1308:   Input Parameter:
1309: . event - The event id

1311:   Example Usage:
1312: .vb
1313:       PetscLogEventDeactivatePush(VEC_SetValues);
1314:         [code where you do not want to log VecSetValues()]
1315:       PetscLogEventDeactivatePop(VEC_SetValues);
1316:         [code where you do want to log VecSetValues()]
1317: .ve

1319:   Level: advanced

1321:   Note:
1322:   The event may be either a pre-defined PETSc event (found in
1323:   include/petsclog.h) or an event number obtained with `PetscLogEventRegister()`).

1325: .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivatePush()`
1326: @*/
1327: PetscErrorCode PetscLogEventDeactivatePop(PetscLogEvent event)
1328: {
1329:   PetscFunctionBegin;
1330:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1331:     PetscLogHandler h = PetscLogHandlers[i].handler;

1333:     if (h) PetscCall(PetscLogHandlerEventDeactivatePop(h, PETSC_DEFAULT, event));
1334:   }
1335:   PetscFunctionReturn(PETSC_SUCCESS);
1336: }

1338: /*@
1339:   PetscLogEventSetActiveAll - Turns on logging of all events

1341:   Not Collective

1343:   Input Parameters:
1344: + event    - The event id
1345: - isActive - The activity flag determining whether the event is logged

1347:   Level: advanced

1349: .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1350: @*/
1351: PetscErrorCode PetscLogEventSetActiveAll(PetscLogEvent event, PetscBool isActive)
1352: {
1353:   PetscLogState state;

1355:   PetscFunctionBegin;
1356:   PetscCall(PetscLogGetState(&state));
1357:   if (state) PetscCall(PetscLogStateEventSetActiveAll(state, event, isActive));
1358:   PetscFunctionReturn(PETSC_SUCCESS);
1359: }

1361: /*
1362:   PetscLogClassSetActive - Activates event logging for a PETSc object class for the current stage

1364:   Not Collective

1366:   Input Parameters:
1367: + stage - A registered `PetscLogStage` (or `PETSC_DEFAULT` for the current stage)
1368: . classid - The event class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1369: - isActive - If `PETSC_FALSE`, events associated with this class are not sent to log handlers.

1371:   Level: developer

1373: .seealso: [](ch_profiling), `PetscLogEventIncludeClass()`, `PetscLogEventActivate()`, `PetscLogEventActivateAll()`, `PetscLogStageSetActive()`
1374: */
1375: static PetscErrorCode PetscLogClassSetActive(PetscLogStage stage, PetscClassId classid, PetscBool isActive)
1376: {
1377:   PetscLogState state;

1379:   PetscFunctionBegin;
1380:   PetscCall(PetscLogGetState(&state));
1381:   if (state) PetscCall(PetscLogStateClassSetActive(state, stage, classid, isActive));
1382:   PetscFunctionReturn(PETSC_SUCCESS);
1383: }

1385: /*@
1386:   PetscLogEventActivateClass - Activates event logging for a PETSc object class for the current stage

1388:   Not Collective

1390:   Input Parameter:
1391: . classid - The event class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.

1393:   Level: developer

1395: .seealso: [](ch_profiling), `PetscLogEventIncludeClass()`, `PetscLogEventExcludeClass()`, `PetscLogEventDeactivateClass()`, `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1396: @*/
1397: PetscErrorCode PetscLogEventActivateClass(PetscClassId classid)
1398: {
1399:   PetscFunctionBegin;
1400:   PetscCall(PetscLogClassSetActive(PETSC_DEFAULT, classid, PETSC_TRUE));
1401:   PetscFunctionReturn(PETSC_SUCCESS);
1402: }

1404: /*@
1405:   PetscLogEventDeactivateClass - Deactivates event logging for a PETSc object class for the current stage

1407:   Not Collective

1409:   Input Parameter:
1410: . classid - The event class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.

1412:   Level: developer

1414: .seealso: [](ch_profiling), `PetscLogEventIncludeClass()`, `PetscLogEventExcludeClass()`, `PetscLogEventActivateClass()`, `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1415: @*/
1416: PetscErrorCode PetscLogEventDeactivateClass(PetscClassId classid)
1417: {
1418:   PetscFunctionBegin;
1419:   PetscCall(PetscLogClassSetActive(PETSC_DEFAULT, classid, PETSC_FALSE));
1420:   PetscFunctionReturn(PETSC_SUCCESS);
1421: }

1423: /*MC
1424:   PetscLogEventSync - Synchronizes the beginning of a user event.

1426:   Synopsis:
1427: #include <petsclog.h>
1428:   PetscErrorCode PetscLogEventSync(PetscLogEvent e, MPI_Comm comm)

1430:   Collective

1432:   Input Parameters:
1433: + e    - `PetscLogEvent` obtained from `PetscLogEventRegister()`
1434: - comm - an MPI communicator

1436:   Example Usage:
1437: .vb
1438:   PetscLogEvent USER_EVENT;

1440:   PetscLogEventRegister("User event", 0, &USER_EVENT);
1441:   PetscLogEventSync(USER_EVENT, PETSC_COMM_WORLD);
1442:   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
1443:   [code segment to monitor]
1444:   PetscLogEventEnd(USER_EVENT, 0, 0, 0 , 0);
1445: .ve

1447:   Level: developer

1449:   Note:
1450:   This routine should be called only if there is not a `PetscObject` available to pass to
1451:   `PetscLogEventBegin()`.

1453: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`
1454: M*/

1456: /*MC
1457:   PetscLogEventBegin - Logs the beginning of a user event.

1459:   Synopsis:
1460: #include <petsclog.h>
1461:   PetscErrorCode PetscLogEventBegin(PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)

1463:   Not Collective

1465:   Input Parameters:
1466: + e  - `PetscLogEvent` obtained from `PetscLogEventRegister()`
1467: . o1 - object associated with the event, or `NULL`
1468: . o2 - object associated with the event, or `NULL`
1469: . o3 - object associated with the event, or `NULL`
1470: - o4 - object associated with the event, or `NULL`

1472:   Fortran Synopsis:
1473:   void PetscLogEventBegin(int e, PetscErrorCode ierr)

1475:   Example Usage:
1476: .vb
1477:   PetscLogEvent USER_EVENT;

1479:   PetscLogDouble user_event_flops;
1480:   PetscLogEventRegister("User event",0, &USER_EVENT);
1481:   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
1482:   [code segment to monitor]
1483:   PetscLogFlops(user_event_flops);
1484:   PetscLogEventEnd(USER_EVENT, 0, 0, 0, 0);
1485: .ve

1487:   Level: intermediate

1489:   Developer Note:
1490:   `PetscLogEventBegin()` and `PetscLogEventBegin()` return error codes instead of explicitly
1491:   handling the errors that occur in the macro directly because other packages that use this
1492:   macros have used them in their own functions or methods that do not return error codes and it
1493:   would be disruptive to change the current behavior.

1495: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventEnd()`, `PetscLogFlops()`
1496: M*/

1498: /*MC
1499:   PetscLogEventEnd - Log the end of a user event.

1501:   Synopsis:
1502: #include <petsclog.h>
1503:   PetscErrorCode PetscLogEventEnd(PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)

1505:   Not Collective

1507:   Input Parameters:
1508: + e  - `PetscLogEvent` obtained from `PetscLogEventRegister()`
1509: . o1 - object associated with the event, or `NULL`
1510: . o2 - object associated with the event, or `NULL`
1511: . o3 - object associated with the event, or `NULL`
1512: - o4 - object associated with the event, or `NULL`

1514:   Fortran Synopsis:
1515:   void PetscLogEventEnd(int e, PetscErrorCode ierr)

1517:   Example Usage:
1518: .vb
1519:   PetscLogEvent USER_EVENT;

1521:   PetscLogDouble user_event_flops;
1522:   PetscLogEventRegister("User event", 0, &USER_EVENT);
1523:   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
1524:   [code segment to monitor]
1525:   PetscLogFlops(user_event_flops);
1526:   PetscLogEventEnd(USER_EVENT, 0, 0, 0, 0);
1527: .ve

1529:   Level: intermediate

1531: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogFlops()`
1532: M*/

1534: /*@C
1535:   PetscLogStageGetPerfInfo - Return the performance information about the given stage

1537:   No Fortran Support

1539:   Input Parameters:
1540: . stage - The stage number or `PETSC_DETERMINE` for the current stage

1542:   Output Parameter:
1543: . info - This structure is filled with the performance information

1545:   Level: intermediate

1547:   Notes:
1548:   This is a low level routine used by the logging functions in PETSc.

1550:   A `PETSCLOGHANDLERDEFAULT` must be running for this to work, having been started either with
1551:   `PetscLogDefaultBegin()` or from the command line with `-log_view`.  If it was not started,
1552:   all performance statistics in `info` will be zeroed.

1554: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogGetDefaultHandler()`
1555: @*/
1556: PetscErrorCode PetscLogStageGetPerfInfo(PetscLogStage stage, PetscEventPerfInfo *info)
1557: {
1558:   PetscLogHandler     handler;
1559:   PetscEventPerfInfo *event_info;

1561:   PetscFunctionBegin;
1562:   PetscAssertPointer(info, 2);
1563:   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
1564:   if (handler) {
1565:     PetscCall(PetscLogHandlerGetStagePerfInfo(handler, stage, &event_info));
1566:     *info = *event_info;
1567:   } else {
1568:     PetscCall(PetscInfo(NULL, "Default log handler is not running, PetscLogStageGetPerfInfo() returning zeros\n"));
1569:     PetscCall(PetscMemzero(info, sizeof(*info)));
1570:   }
1571:   PetscFunctionReturn(PETSC_SUCCESS);
1572: }

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

1577:   No Fortran Support

1579:   Input Parameters:
1580: + stage - The stage number or `PETSC_DETERMINE` for the current stage
1581: - event - The event number

1583:   Output Parameter:
1584: . info - This structure is filled with the performance information

1586:   Level: intermediate

1588:   Note:
1589:   This is a low level routine used by the logging functions in PETSc

1591:   A `PETSCLOGHANDLERDEFAULT` must be running for this to work, having been started either with
1592:   `PetscLogDefaultBegin()` or from the command line with `-log_view`.  If it was not started,
1593:   all performance statistics in `info` will be zeroed.

1595: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogGetDefaultHandler()`
1596: @*/
1597: PetscErrorCode PetscLogEventGetPerfInfo(PetscLogStage stage, PetscLogEvent event, PetscEventPerfInfo *info)
1598: {
1599:   PetscLogHandler     handler;
1600:   PetscEventPerfInfo *event_info;

1602:   PetscFunctionBegin;
1603:   PetscAssertPointer(info, 3);
1604:   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
1605:   if (handler) {
1606:     PetscCall(PetscLogHandlerGetEventPerfInfo(handler, stage, event, &event_info));
1607:     *info = *event_info;
1608:   } else {
1609:     PetscCall(PetscInfo(NULL, "Default log handler is not running, PetscLogEventGetPerfInfo() returning zeros\n"));
1610:     PetscCall(PetscMemzero(info, sizeof(*info)));
1611:   }
1612:   PetscFunctionReturn(PETSC_SUCCESS);
1613: }

1615: /*@
1616:   PetscLogEventSetDof - Set the nth number of degrees of freedom of a numerical problem associated with this event

1618:   Not Collective

1620:   Input Parameters:
1621: + event - The event id to log
1622: . n     - The dof index, in [0, 8)
1623: - dof   - The number of dofs

1625:   Options Database Key:
1626: . -log_view - Activates log summary

1628:   Level: developer

1630:   Note:
1631:   This is to enable logging of convergence

1633: .seealso: `PetscLogEventSetError()`, `PetscLogEventRegister()`, `PetscLogGetDefaultHandler()`
1634: @*/
1635: PetscErrorCode PetscLogEventSetDof(PetscLogEvent event, PetscInt n, PetscLogDouble dof)
1636: {
1637:   PetscFunctionBegin;
1638:   PetscCheck(!(n < 0) && !(n > 7), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Error index %" PetscInt_FMT " is not in [0, 8)", n);
1639:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1640:     PetscLogHandler h = PetscLogHandlers[i].handler;

1642:     if (h) {
1643:       PetscEventPerfInfo *event_info;

1645:       PetscCall(PetscLogHandlerGetEventPerfInfo(h, PETSC_DEFAULT, event, &event_info));
1646:       if (event_info) event_info->dof[n] = dof;
1647:     }
1648:   }
1649:   PetscFunctionReturn(PETSC_SUCCESS);
1650: }

1652: /*@
1653:   PetscLogEventSetError - Set the nth error associated with a numerical problem associated with this event

1655:   Not Collective

1657:   Input Parameters:
1658: + event - The event id to log
1659: . n     - The error index, in [0, 8)
1660: - error - The error

1662:   Options Database Key:
1663: . -log_view - Activates log summary

1665:   Level: developer

1667:   Notes:
1668:   This is to enable logging of convergence, and enable users to interpret the errors as they wish. For example,
1669:   as different norms, or as errors for different fields

1671:   This is a low level routine used by the logging functions in PETSc

1673: .seealso: `PetscLogEventSetDof()`, `PetscLogEventRegister()`, `PetscLogGetDefaultHandler()`
1674: @*/
1675: PetscErrorCode PetscLogEventSetError(PetscLogEvent event, PetscInt n, PetscLogDouble error)
1676: {
1677:   PetscFunctionBegin;
1678:   PetscCheck(!(n < 0) && !(n > 7), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Error index %" PetscInt_FMT " is not in [0, 8)", n);
1679:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1680:     PetscLogHandler h = PetscLogHandlers[i].handler;

1682:     if (h) {
1683:       PetscEventPerfInfo *event_info;

1685:       PetscCall(PetscLogHandlerGetEventPerfInfo(h, PETSC_DEFAULT, event, &event_info));
1686:       if (event_info) event_info->errors[n] = error;
1687:     }
1688:   }
1689:   PetscFunctionReturn(PETSC_SUCCESS);
1690: }

1692: /*@
1693:   PetscLogEventGetId - Returns the event id when given the event name.

1695:   Not Collective

1697:   Input Parameter:
1698: . name - The event name

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

1703:   Level: intermediate

1705: .seealso: [](ch_profiling), `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogStageGetId()`
1706: @*/
1707: PetscErrorCode PetscLogEventGetId(const char name[], PetscLogEvent *event)
1708: {
1709:   PetscLogState state;

1711:   PetscFunctionBegin;
1712:   *event = -1;
1713:   PetscCall(PetscLogGetState(&state));
1714:   if (state) PetscCall(PetscLogStateGetEventFromName(state, name, event));
1715:   PetscFunctionReturn(PETSC_SUCCESS);
1716: }

1718: /*@
1719:   PetscLogEventGetName - Returns the event name when given the event id.

1721:   Not Collective

1723:   Input Parameter:
1724: . event - The event

1726:   Output Parameter:
1727: . name - The event name

1729:   Level: intermediate

1731: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
1732: @*/
1733: PetscErrorCode PetscLogEventGetName(PetscLogEvent event, const char *name[])
1734: {
1735:   PetscLogEventInfo event_info;
1736:   PetscLogState     state;

1738:   PetscFunctionBegin;
1739:   *name = NULL;
1740:   PetscCall(PetscLogGetState(&state));
1741:   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
1742:   PetscCall(PetscLogStateEventGetInfo(state, event, &event_info));
1743:   *name = event_info.name;
1744:   PetscFunctionReturn(PETSC_SUCCESS);
1745: }

1747: /*@
1748:   PetscLogEventsPause - Put event logging into "paused" mode: timers and counters for in-progress events are paused, and any events that happen before logging is resumed with `PetscLogEventsResume()` are logged in the "Main Stage" of execution.

1750:   Not collective

1752:   Level: advanced

1754:   Notes:
1755:   When an external library or runtime has is initialized it can involve lots of setup time that skews the statistics of any unrelated running events: this function is intended to isolate such calls in the default log summary (`PetscLogDefaultBegin()`, `PetscLogView()`).

1757:   Other log handlers (such as the nested handler, `PetscLogNestedBegin()`) will ignore this function.

1759: .seealso: [](ch_profiling), `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`, `PetscLogEventsResume()`, `PetscLogGetDefaultHandler()`
1760: @*/
1761: PetscErrorCode PetscLogEventsPause(void)
1762: {
1763:   PetscFunctionBegin;
1764:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1765:     PetscLogHandler h = PetscLogHandlers[i].handler;

1767:     if (h) PetscCall(PetscLogHandlerEventsPause(h));
1768:   }
1769:   PetscFunctionReturn(PETSC_SUCCESS);
1770: }

1772: /*@
1773:   PetscLogEventsResume - Return logging to normal behavior after it was paused with `PetscLogEventsPause()`.

1775:   Not collective

1777:   Level: advanced

1779: .seealso: [](ch_profiling), `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`, `PetscLogEventsPause()`, `PetscLogGetDefaultHandler()`
1780: @*/
1781: PetscErrorCode PetscLogEventsResume(void)
1782: {
1783:   PetscFunctionBegin;
1784:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1785:     PetscLogHandler h = PetscLogHandlers[i].handler;

1787:     if (h) PetscCall(PetscLogHandlerEventsResume(h));
1788:   }
1789:   PetscFunctionReturn(PETSC_SUCCESS);
1790: }

1792: /*------------------------------------------------ Class Functions --------------------------------------------------*/

1794: /*MC
1795:    PetscLogObjectCreate - Log the creation of a `PetscObject`

1797:    Synopsis:
1798: #include <petsclog.h>
1799:    PetscErrorCode PetscLogObjectCreate(PetscObject h)

1801:    Not Collective

1803:    Input Parameters:
1804: .  h - A `PetscObject`

1806:    Level: developer

1808:    Developer Note:
1809:      Called internally by PETSc when creating objects: users do not need to call this directly.
1810:      Notification of the object creation is sent to each `PetscLogHandler` that is running.

1812: .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogObjectDestroy()`
1813: M*/

1815: /*MC
1816:    PetscLogObjectDestroy - Logs the destruction of a `PetscObject`

1818:    Synopsis:
1819: #include <petsclog.h>
1820:    PetscErrorCode PetscLogObjectDestroy(PetscObject h)

1822:    Not Collective

1824:    Input Parameters:
1825: .  h - A `PetscObject`

1827:    Level: developer

1829:    Developer Note:
1830:      Called internally by PETSc when destroying objects: users do not need to call this directly.
1831:      Notification of the object creation is sent to each `PetscLogHandler` that is running.

1833: .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogObjectCreate()`
1834: M*/

1836: /*@
1837:   PetscLogClassGetClassId - Returns the `PetscClassId` when given the class name.

1839:   Not Collective

1841:   Input Parameter:
1842: . name - The class name

1844:   Output Parameter:
1845: . classid - The `PetscClassId` id, or -1 if no class with that name exists

1847:   Level: intermediate

1849: .seealso: [](ch_profiling), `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogStageGetId()`
1850: @*/
1851: PetscErrorCode PetscLogClassGetClassId(const char name[], PetscClassId *classid)
1852: {
1853:   PetscLogClass     log_class;
1854:   PetscLogClassInfo class_info;
1855:   PetscLogState     state;

1857:   PetscFunctionBegin;
1858:   *classid = -1;
1859:   PetscCall(PetscLogGetState(&state));
1860:   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
1861:   PetscCall(PetscLogStateGetClassFromName(state, name, &log_class));
1862:   if (log_class < 0) {
1863:     *classid = -1;
1864:     PetscFunctionReturn(PETSC_SUCCESS);
1865:   }
1866:   PetscCall(PetscLogStateClassGetInfo(state, log_class, &class_info));
1867:   *classid = class_info.classid;
1868:   PetscFunctionReturn(PETSC_SUCCESS);
1869: }

1871: /*@C
1872:   PetscLogClassIdGetName - Returns a `PetscClassId`'s name.

1874:   Not Collective

1876:   Input Parameter:
1877: . classid - A `PetscClassId`

1879:   Output Parameter:
1880: . name - The class name

1882:   Level: intermediate

1884: .seealso: [](ch_profiling), `PetscLogClassRegister()`, `PetscLogClassBegin()`, `PetscLogClassEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadClass()`
1885: @*/
1886: PetscErrorCode PetscLogClassIdGetName(PetscClassId classid, const char **name)
1887: {
1888:   PetscLogClass     log_class;
1889:   PetscLogClassInfo class_info;
1890:   PetscLogState     state;

1892:   PetscFunctionBegin;
1893:   PetscCall(PetscLogGetState(&state));
1894:   PetscCall(PetscLogStateGetClassFromClassId(state, classid, &log_class));
1895:   PetscCall(PetscLogStateClassGetInfo(state, log_class, &class_info));
1896:   *name = class_info.name;
1897:   PetscFunctionReturn(PETSC_SUCCESS);
1898: }

1900: /*------------------------------------------------ Output Functions -------------------------------------------------*/
1901: /*@
1902:   PetscLogDump - Dumps logs of objects to a file. This file is intended to
1903:   be read by bin/petscview. This program no longer exists.

1905:   Collective on `PETSC_COMM_WORLD`

1907:   Input Parameter:
1908: . sname - an optional file name

1910:   Example Usage:
1911: .vb
1912:   PetscInitialize(...);
1913:   PetscLogDefaultBegin();
1914:   // ... code ...
1915:   PetscLogDump(filename);
1916:   PetscFinalize();
1917: .ve

1919:   Level: advanced

1921:   Note:
1922:   The default file name is Log.<rank> where <rank> is the MPI process rank. If no name is specified,
1923:   this file will be used.

1925: .seealso: [](ch_profiling), `PetscLogDefaultBegin()`, `PetscLogView()`, `PetscLogGetDefaultHandler()`
1926: @*/
1927: PetscErrorCode PetscLogDump(const char sname[])
1928: {
1929:   PetscLogHandler handler;

1931:   PetscFunctionBegin;
1932:   PetscCall(PetscLogGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
1933:   PetscCall(PetscLogHandlerDump(handler, sname));
1934:   PetscFunctionReturn(PETSC_SUCCESS);
1935: }

1937: /*@
1938:   PetscLogMPEDump - Dumps the MPE logging info to file for later use with Jumpshot.

1940:   Collective on `PETSC_COMM_WORLD`

1942:   Input Parameter:
1943: . sname - filename for the MPE logfile

1945:   Level: advanced

1947: .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogMPEBegin()`
1948: @*/
1949: PetscErrorCode PetscLogMPEDump(const char sname[])
1950: {
1951:   PetscFunctionBegin;
1952:   #if defined(PETSC_HAVE_MPE)
1953:   if (PetscBeganMPE) {
1954:     char name[PETSC_MAX_PATH_LEN];

1956:     PetscCall(PetscInfo(0, "Finalizing MPE.\n"));
1957:     if (sname) {
1958:       PetscCall(PetscStrncpy(name, sname, sizeof(name)));
1959:     } else {
1960:       PetscCall(PetscGetProgramName(name, sizeof(name)));
1961:     }
1962:     PetscCall(MPE_Finish_log(name));
1963:   } else {
1964:     PetscCall(PetscInfo(0, "Not finalizing MPE (not started by PETSc).\n"));
1965:   }
1966:   #else
1967:   SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "PETSc was configured without MPE support, reconfigure with --with-mpe or --download-mpe");
1968:   #endif
1969:   PetscFunctionReturn(PETSC_SUCCESS);
1970: }

1972: /*@
1973:   PetscLogView - Prints a summary of the logging.

1975:   Collective

1977:   Input Parameter:
1978: . viewer - an ASCII viewer

1980:   Options Database Keys:
1981: + -log_view [:filename]                    - Prints summary of log information
1982: . -log_view :filename.py:ascii_info_detail - Saves logging information from each process as a Python file
1983: . -log_view :filename.xml:ascii_xml        - Saves a summary of the logging information in a nested format (see below for how to view it)
1984: . -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)
1985: . -log_view_memory                         - Also display memory usage in each event
1986: . -log_view_gpu_time                       - Also display time in each event for GPU kernels (Note this may slow the computation)
1987: . -log_view_gpu_energy                     - Also display energy (estimated with power*gtime) in Joules for GPU kernels
1988: . -log_view_gpu_energy_meter               - [Experimental] Also display energy (readings from energy meters) in Joules for GPU kernels. This option is ignored if -log_view_gpu_energy is provided.
1989: . -log_all                                 - Saves a file Log.rank for each MPI rank with details of each step of the computation
1990: - -log_trace [filename]                    - Displays a trace of what each process is doing

1992:   Level: beginner

1994:   Notes:
1995:   It is possible to control the logging programmatically but we recommend using the options database approach whenever possible
1996:   By default the summary is printed to stdout.

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

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

2002:   To view the nested XML format filename.xml first copy  ${PETSC_DIR}/share/petsc/xml/performance_xml2html.xsl to the current
2003:   directory then open filename.xml with your browser. Specific notes for certain browsers
2004: .vb
2005:     Firefox and Internet explorer - simply open the file
2006:     Google Chrome - you must start up Chrome with the option --allow-file-access-from-files
2007:     Safari - see https://ccm.net/faq/36342-safari-how-to-enable-local-file-access
2008: .ve
2009:   or one can use the package <http://xmlsoft.org/XSLT/xsltproc2.html> to translate the xml file to html and then open it with
2010:   your browser.
2011:   Alternatively, use the script ${PETSC_DIR}/lib/petsc/bin/petsc-performance-view to automatically open a new browser
2012:   window and render the XML log file contents.

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

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

2020: .seealso: [](ch_profiling), `PetscLogDefaultBegin()`, `PetscLogDump()`
2021: @*/
2022: PetscErrorCode PetscLogView(PetscViewer viewer)
2023: {
2024:   PetscBool         isascii;
2025:   PetscViewerFormat format;
2026:   int               stage;
2027:   PetscLogState     state;
2028:   PetscIntStack     temp_stack;
2029:   PetscLogHandler   handler;
2030:   PetscBool         is_empty;

2032:   PetscFunctionBegin;
2033:   PetscCall(PetscLogGetState(&state));
2034:   /* Pop off any stages the user forgot to remove */
2035:   PetscCall(PetscIntStackCreate(&temp_stack));
2036:   PetscCall(PetscLogStateGetCurrentStage(state, &stage));
2037:   while (stage >= 0) {
2038:     PetscCall(PetscLogStagePop());
2039:     PetscCall(PetscIntStackPush(temp_stack, stage));
2040:     PetscCall(PetscLogStateGetCurrentStage(state, &stage));
2041:   }
2042:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
2043:   PetscCheck(isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Currently can only view logging to ASCII");
2044:   PetscCall(PetscViewerGetFormat(viewer, &format));
2045:   if (format == PETSC_VIEWER_ASCII_XML || format == PETSC_VIEWER_ASCII_FLAMEGRAPH) {
2046:     PetscCall(PetscLogGetHandler(PETSCLOGHANDLERNESTED, &handler));
2047:     PetscCall(PetscLogHandlerView(handler, viewer));
2048:   } else {
2049:     PetscCall(PetscLogGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
2050:     PetscCall(PetscLogHandlerView(handler, viewer));
2051:   }
2052:   PetscCall(PetscIntStackEmpty(temp_stack, &is_empty));
2053:   while (!is_empty) {
2054:     PetscCall(PetscIntStackPop(temp_stack, &stage));
2055:     PetscCall(PetscLogStagePush(stage));
2056:     PetscCall(PetscIntStackEmpty(temp_stack, &is_empty));
2057:   }
2058:   PetscCall(PetscIntStackDestroy(temp_stack));
2059:   PetscFunctionReturn(PETSC_SUCCESS);
2060: }

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

2065:   Collective on `PETSC_COMM_WORLD`

2067:   Level: developer

2069: .seealso: [](ch_profiling), `PetscLogView()`
2070: @*/
2071: PetscErrorCode PetscLogViewFromOptions(void)
2072: {
2073:   PetscInt          n_max = PETSC_LOG_VIEW_FROM_OPTIONS_MAX;
2074:   PetscViewer       viewers[PETSC_LOG_VIEW_FROM_OPTIONS_MAX];
2075:   PetscViewerFormat formats[PETSC_LOG_VIEW_FROM_OPTIONS_MAX];
2076:   PetscBool         flg;

2078:   PetscFunctionBegin;
2079:   PetscCall(PetscOptionsCreateViewers(PETSC_COMM_WORLD, NULL, NULL, "-log_view", &n_max, viewers, formats, &flg));
2080:   /*
2081:      PetscLogHandlerView_Default_Info() wants to be sure that the only objects still around are these viewers, so keep track of how many there are
2082:    */
2083:   PetscLogNumViewersCreated = n_max;
2084:   for (PetscInt i = 0; i < n_max; i++) {
2085:     PetscInt refct;

2087:     PetscCall(PetscViewerPushFormat(viewers[i], formats[i]));
2088:     PetscCall(PetscLogView(viewers[i]));
2089:     PetscCall(PetscViewerPopFormat(viewers[i]));
2090:     PetscCall(PetscObjectGetReference((PetscObject)viewers[i], &refct));
2091:     PetscCall(PetscViewerDestroy(&viewers[i]));
2092:     if (refct == 1) PetscLogNumViewersDestroyed++;
2093:   }
2094:   PetscLogNumViewersDestroyed = 0;
2095:   PetscFunctionReturn(PETSC_SUCCESS);
2096: }

2098: PETSC_INTERN PetscErrorCode PetscLogHandlerNestedSetThreshold(PetscLogHandler, PetscLogDouble, PetscLogDouble *);

2100: /*@
2101:   PetscLogSetThreshold - Set the threshold time for logging the events; this is a percentage out of 100, so 1. means any event
2102:   that takes 1 or more percent of the time.

2104:   Logically Collective on `PETSC_COMM_WORLD`

2106:   Input Parameter:
2107: . newThresh - the threshold to use

2109:   Output Parameter:
2110: . oldThresh - the previously set threshold value

2112:   Options Database Keys:
2113: . -log_view :filename.xml:ascii_xml - Prints an XML summary of flop and timing information to the file

2115:   Example Usage:
2116: .vb
2117:   PetscInitialize(...);
2118:   PetscLogNestedBegin();
2119:   PetscLogSetThreshold(0.1,&oldthresh);
2120:   // ... code ...
2121:   PetscLogView(viewer);
2122:   PetscFinalize();
2123: .ve

2125:   Level: advanced

2127:   Note:
2128:   This threshold is only used by the nested log handler

2130: .seealso: `PetscLogDump()`, `PetscLogView()`, `PetscLogTraceBegin()`, `PetscLogDefaultBegin()`,
2131:           `PetscLogNestedBegin()`
2132: @*/
2133: PetscErrorCode PetscLogSetThreshold(PetscLogDouble newThresh, PetscLogDouble *oldThresh)
2134: {
2135:   PetscLogHandler handler;

2137:   PetscFunctionBegin;
2138:   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERNESTED, &handler));
2139:   PetscCall(PetscLogHandlerNestedSetThreshold(handler, newThresh, oldThresh));
2140:   PetscFunctionReturn(PETSC_SUCCESS);
2141: }

2143: /*----------------------------------------------- Counter Functions -------------------------------------------------*/
2144: /*@
2145:   PetscGetFlops - Returns the number of flops used on this processor
2146:   since the program began.

2148:   Not Collective

2150:   Output Parameter:
2151: . flops - number of floating point operations

2153:   Level: intermediate

2155:   Notes:
2156:   A global counter logs all PETSc flop counts.  The user can use
2157:   `PetscLogFlops()` to increment this counter to include flops for the
2158:   application code.

2160:   A separate counter `PetscLogGpuFlops()` logs the flops that occur on any GPU associated with this MPI rank

2162: .seealso: [](ch_profiling), `PetscLogGpuFlops()`, `PetscTime()`, `PetscLogFlops()`
2163: @*/
2164: PetscErrorCode PetscGetFlops(PetscLogDouble *flops)
2165: {
2166:   PetscFunctionBegin;
2167:   *flops = petsc_TotalFlops;
2168:   PetscFunctionReturn(PETSC_SUCCESS);
2169: }

2171: /*@C
2172:   PetscLogObjectState - Record information about an object with the default log handler

2174:   Not Collective

2176:   Input Parameters:
2177: + obj    - the `PetscObject`
2178: . format - a printf-style format string
2179: - ...    - printf arguments to format

2181:   Level: developer

2183: .seealso: [](ch_profiling), `PetscLogObjectCreate()`, `PetscLogObjectDestroy()`, `PetscLogGetDefaultHandler()`
2184: @*/
2185: PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...)
2186: {
2187:   PetscFunctionBegin;
2188:   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
2189:     PetscLogHandler h = PetscLogHandlers[i].handler;

2191:     if (h) {
2192:       va_list Argp;
2193:       va_start(Argp, format);
2194:       PetscCall(PetscLogHandlerLogObjectState_Internal(h, obj, format, Argp));
2195:       va_end(Argp);
2196:     }
2197:   }
2198:   PetscFunctionReturn(PETSC_SUCCESS);
2199: }

2201: /*MC
2202:   PetscLogFlops - Adds floating point operations to the global counter.

2204:   Synopsis:
2205: #include <petsclog.h>
2206:   PetscErrorCode PetscLogFlops(PetscLogDouble f)

2208:   Not Collective

2210:   Input Parameter:
2211: . f - flop counter

2213:   Example Usage:
2214: .vb
2215:   PetscLogEvent USER_EVENT;

2217:   PetscLogEventRegister("User event", 0, &USER_EVENT);
2218:   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
2219:   [code segment to monitor]
2220:   PetscLogFlops(user_flops)
2221:   PetscLogEventEnd(USER_EVENT, 0, 0, 0, 0);
2222: .ve

2224:   Level: intermediate

2226:   Note:
2227:    A global counter logs all PETSc flop counts. The user can use PetscLogFlops() to increment
2228:    this counter to include flops for the application code.

2230: .seealso: [](ch_profiling), `PetscLogGpuFlops()`, `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscGetFlops()`
2231: M*/

2233: /*MC
2234:   PetscPreLoadBegin - Begin a segment of code that may be preloaded (run twice) to get accurate
2235:   timings

2237:   Synopsis:
2238: #include <petsclog.h>
2239:   void PetscPreLoadBegin(PetscBool flag, char *name);

2241:   Not Collective

2243:   Input Parameters:
2244: + flag - `PETSC_TRUE` to run twice, `PETSC_FALSE` to run once, may be overridden with command
2245:          line option `-preload true|false`
2246: - name - name of first stage (lines of code timed separately with `-log_view`) to be preloaded

2248:   Example Usage:
2249: .vb
2250:   PetscPreLoadBegin(PETSC_TRUE, "first stage");
2251:   // lines of code
2252:   PetscPreLoadStage("second stage");
2253:   // lines of code
2254:   PetscPreLoadEnd();
2255: .ve

2257:   Level: intermediate

2259:   Note:
2260:   Only works in C/C++, not Fortran

2262:   Flags available within the macro\:
2263: + PetscPreLoadingUsed - `PETSC_TRUE` if we are or have done preloading
2264: . PetscPreLoadingOn   - `PETSC_TRUE` if it is CURRENTLY doing preload
2265: . PetscPreLoadIt      - `0` for the first computation (with preloading turned off it is only
2266:                         `0`) `1`  for the second
2267: - PetscPreLoadMax     - number of times it will do the computation, only one when preloading is
2268:                         turned on

2270:   The first two variables are available throughout the program, the second two only between the
2271:   `PetscPreLoadBegin()` and `PetscPreLoadEnd()`

2273: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
2274: M*/

2276: /*MC
2277:   PetscPreLoadEnd - End a segment of code that may be preloaded (run twice) to get accurate
2278:   timings

2280:   Synopsis:
2281: #include <petsclog.h>
2282:   void PetscPreLoadEnd(void);

2284:   Not Collective

2286:   Example Usage:
2287: .vb
2288:   PetscPreLoadBegin(PETSC_TRUE, "first stage");
2289:   // lines of code
2290:   PetscPreLoadStage("second stage");
2291:   // lines of code
2292:   PetscPreLoadEnd();
2293: .ve

2295:   Level: intermediate

2297:   Note:
2298:   Only works in C/C++ not Fortran

2300: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadStage()`
2301: M*/

2303: /*MC
2304:   PetscPreLoadStage - Start a new segment of code to be timed separately to get accurate timings

2306:   Synopsis:
2307: #include <petsclog.h>
2308:   void PetscPreLoadStage(char *name);

2310:   Not Collective

2312:   Example Usage:
2313: .vb
2314:   PetscPreLoadBegin(PETSC_TRUE,"first stage");
2315:   // lines of code
2316:   PetscPreLoadStage("second stage");
2317:   // lines of code
2318:   PetscPreLoadEnd();
2319: .ve

2321:   Level: intermediate

2323:   Note:
2324:   Only works in C/C++ not Fortran

2326: .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`
2327: M*/

2329:   #if PetscDefined(HAVE_DEVICE)
2330: #include <petsc/private/deviceimpl.h>

2332: /*@
2333:   PetscLogGpuTime - turn on the logging of GPU time for GPU kernels

2335:   Options Database Key:
2336: . -log_view_gpu_time - provide the GPU times for all events in the `-log_view` output

2338:   Level: advanced

2340:   Notes:
2341:   Turning on the timing of the GPU kernels can slow down the entire computation and should only
2342:   be used when studying the performance of individual operations on GPU such as vector operations and
2343:   matrix-vector operations.

2345:   If this option is not used then times for most of the events in the `-log_view` output will be listed as Nan, indicating the times are not available

2347:   This routine should only be called once near the beginning of the program. Once it is started
2348:   it cannot be turned off.

2350: .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuFlops()`, `PetscLogGpuTimeEnd()`, `PetscLogGpuTimeBegin()`
2351: @*/
2352: PetscErrorCode PetscLogGpuTime(void)
2353: {
2354:   PetscFunctionBegin;
2355:   PetscCheck(petsc_gtime == 0.0, PETSC_COMM_SELF, PETSC_ERR_SUP, "GPU logging has already been turned on");
2356:   PetscLogGpuTimeFlag = PETSC_TRUE;
2357:   PetscFunctionReturn(PETSC_SUCCESS);
2358: }

2360: /*@
2361:   PetscLogGpuTimeBegin - Start timer for device

2363:   Level: intermediate

2365:   Notes:
2366:   When GPU is enabled, the timer is run on the GPU, it is a separate logging of time
2367:   devoted to GPU computations (excluding kernel launch times).

2369:   When GPU is not available, the timer is run on the CPU, it is a separate logging of
2370:   time devoted to GPU computations (including kernel launch times).

2372:   There is no need to call WaitForCUDA() or WaitForHIP() between `PetscLogGpuTimeBegin()` and
2373:   `PetscLogGpuTimeEnd()`

2375:   This timer should NOT include times for data transfers between the GPU and CPU, nor setup
2376:   actions such as allocating space.

2378:   The regular logging captures the time for data transfers and any CPU activities during the
2379:   event. It is used to compute the flop rate on the GPU as it is actively engaged in running a
2380:   kernel.

2382:   Developer Notes:
2383:   The GPU event timer captures the execution time of all the kernels launched in the default
2384:   stream by the CPU between `PetscLogGpuTimeBegin()` and `PetscLogGpuTimeEnd()`.

2386:   `PetscLogGpuTimeBegin()` and `PetscLogGpuTimeEnd()` insert the begin and end events into the
2387:   default stream (stream 0). The device will record a time stamp for the event when it reaches
2388:   that event in the stream. The function xxxEventSynchronize() is called in
2389:   `PetscLogGpuTimeEnd()` to block CPU execution, but not continued GPU execution, until the
2390:   timer event is recorded.

2392: .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuFlops()`, `PetscLogGpuTimeEnd()`, `PetscLogGpuTime()`
2393: @*/
2394: PetscErrorCode PetscLogGpuTimeBegin(void)
2395: {
2396:   PetscBool isActive;

2398:   PetscFunctionBegin;
2399:   PetscCall(PetscLogEventBeginIsActive(&isActive));
2400:   if (!isActive || !PetscLogGpuTimeFlag) PetscFunctionReturn(PETSC_SUCCESS);
2401:   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2402:     PetscDeviceContext dctx;

2404:     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2405:     PetscCall(PetscDeviceContextBeginTimer_Internal(dctx));
2406:   } else {
2407:     PetscCall(PetscTimeSubtract(&petsc_gtime));
2408:   }
2409:   PetscFunctionReturn(PETSC_SUCCESS);
2410: }

2412: /*@
2413:   PetscLogGpuTimeEnd - Stop timer for device

2415:   Level: intermediate

2417: .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuFlops()`, `PetscLogGpuTimeBegin()`
2418: @*/
2419: PetscErrorCode PetscLogGpuTimeEnd(void)
2420: {
2421:   PetscBool isActive;

2423:   PetscFunctionBegin;
2424:   PetscCall(PetscLogEventEndIsActive(&isActive));
2425:   if (!isActive || !PetscLogGpuTimeFlag) PetscFunctionReturn(PETSC_SUCCESS);
2426:   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2427:     PetscDeviceContext dctx;
2428:     PetscLogDouble     elapsed;

2430:     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2431:     PetscCall(PetscDeviceContextEndTimer_Internal(dctx, &elapsed));
2432:     petsc_gtime += (elapsed / 1000.0);
2433:     #if PetscDefined(HAVE_CUDA_VERSION_12_2PLUS)
2434:     if (PetscLogGpuEnergyFlag) {
2435:       PetscLogDouble power;
2436:       PetscCall(PetscDeviceContextGetPower_Internal(dctx, &power));
2437:       petsc_genergy += (power * elapsed / 1000000.0); // convert to Joules
2438:     }
2439:     #endif
2440:   } else {
2441:     PetscCall(PetscTimeAdd(&petsc_gtime));
2442:   }
2443:   PetscFunctionReturn(PETSC_SUCCESS);
2444: }

2446: /*@
2447:   PetscLogGpuEnergy - turn on the logging of GPU energy (estimated with power*gtime) for GPU kernels

2449:   Options Database Key:
2450: . -log_view_gpu_energy - provide the GPU energy consumption (estimated with power*gtime) for all events in the `-log_view` output

2452:   Level: advanced

2454:   Note:
2455:   This option is mutually exclusive to `-log_view_gpu_energy_meter`.

2457:   Developer Note:
2458:   This option turns on energy monitoring of GPU kernels and requires CUDA version >= 12.2. The energy consumption is estimated as
2459:   instant_power * gpu_kernel_time. Due to the delay in NVML power sampling, we read the instantaneous power draw at the end of each
2460:   event using `nvmlDeviceGetFieldValues()` with the field ID `NVML_FI_DEV_POWER_INSTANT`.

2462: .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuEnergyMeter()`
2463: @*/
2464: PetscErrorCode PetscLogGpuEnergy(void)
2465: {
2466:   PetscFunctionBegin;
2467:   PetscCheck(PetscDefined(HAVE_CUDA_VERSION_12_2PLUS), PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "-log_view_gpu_energy requires CUDA version >= 12.2");
2468:   PetscCheck(petsc_genergy == 0.0, PETSC_COMM_SELF, PETSC_ERR_SUP, "GPU energy logging has already been turned on");
2469:   PetscLogGpuEnergyFlag      = PETSC_TRUE;
2470:   PetscLogGpuEnergyMeterFlag = PETSC_FALSE;
2471:   PetscFunctionReturn(PETSC_SUCCESS);
2472: }

2474: /*@
2475:   PetscLogGpuEnergyMeter - turn on the logging of GPU energy (readings from energy meters) for GPU kernels

2477:   Options Database Key:
2478: . -log_view_gpu_energy_meter - provide the GPU energy (readings from energy meters) consumption for all events in the `-log_view` output

2480:   Level: advanced

2482:   Note:
2483:   This option is mutually exclusive to `-log_view_gpu_energy`.

2485:   Developer Note:
2486:   This option turns on energy monitoring of GPU kernels. The energy consumption is measured directly using the NVML API
2487:   `nvmlDeviceGetTotalEnergyConsumption()`, which returns the total energy used by the GPU since the driver was last initialized.
2488:   For newer GPUs, energy readings are updated every 20-100ms, so this approach may be inaccurate for short-duration GPU events.

2490: @*/
2491: PetscErrorCode PetscLogGpuEnergyMeter(void)
2492: {
2493:   PetscFunctionBegin;
2494:   PetscCheck(petsc_genergy == 0.0, PETSC_COMM_SELF, PETSC_ERR_SUP, "GPU energy logging has already been turned on");
2495:   PetscLogGpuEnergyMeterFlag = PETSC_TRUE;
2496:   PetscLogGpuEnergyFlag      = PETSC_FALSE;
2497:   PetscFunctionReturn(PETSC_SUCCESS);
2498: }

2500: /*@
2501:   PetscLogGpuEnergyMeterBegin - Start energy meter for device

2503:   Level: intermediate

2505:   Notes:
2506:   The GPU event energy meter captures the energy used by the GPU between `PetscLogGpuEnergyMeterBegin()` and `PetscLogGpuEnergyMeterEnd()`.

2508:   `PetscLogGpuEnergyMeterBegin()` and `PetscLogGpuEnergyMeterEnd()` collect the energy readings using `nvmlDeviceGetTotalEnergyConsumption()`.
2509:   The function `cupmStreamSynchronize()` is called before the energy query to ensure completion.

2511: .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuEnergyMeterEnd()`, `PetscLogGpuEnergyMeter()`
2512: @*/
2513: PetscErrorCode PetscLogGpuEnergyMeterBegin(void)
2514: {
2515:   PetscBool isActive;

2517:   PetscFunctionBegin;
2518:   PetscCall(PetscLogEventBeginIsActive(&isActive));
2519:   if (!isActive || !PetscLogGpuEnergyMeterFlag) PetscFunctionReturn(PETSC_SUCCESS);
2520:   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2521:     PetscDeviceContext dctx;

2523:     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2524:     PetscCall(PetscDeviceContextBeginEnergyMeter_Internal(dctx));
2525:   }
2526:   PetscFunctionReturn(PETSC_SUCCESS);
2527: }

2529: /*@
2530:   PetscLogGpuEnergyMeterEnd - Stop energy meter for device

2532:   Level: intermediate

2534: .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuEnergyMeterBegin()`
2535: @*/
2536: PetscErrorCode PetscLogGpuEnergyMeterEnd(void)
2537: {
2538:   PetscBool isActive;

2540:   PetscFunctionBegin;
2541:   PetscCall(PetscLogEventEndIsActive(&isActive));
2542:   if (!isActive || !PetscLogGpuEnergyMeterFlag) PetscFunctionReturn(PETSC_SUCCESS);
2543:   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2544:     PetscDeviceContext dctx;
2545:     PetscLogDouble     energy;

2547:     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2548:     PetscCall(PetscDeviceContextEndEnergyMeter_Internal(dctx, &energy));
2549:     petsc_genergy_meter += (energy / 1000.0); // convert to Joules
2550:   }
2551:   PetscFunctionReturn(PETSC_SUCCESS);
2552: }
2553:   #endif /* end of PETSC_HAVE_DEVICE */

2555: #endif /* PETSC_USE_LOG*/

2557: /* -- Utility functions for logging from Fortran -- */

2559: PETSC_EXTERN PetscErrorCode PetscASend(int count, int datatype)
2560: {
2561:   PetscFunctionBegin;
2562: #if PetscDefined(USE_LOG)
2563:   PetscCall(PetscAddLogDouble(&petsc_send_ct, &petsc_send_ct_th, 1));
2564:   #if !defined(MPIUNI_H) && !defined(PETSC_HAVE_BROKEN_RECURSIVE_MACRO)
2565:   PetscCall(PetscMPITypeSize(count, MPI_Type_f2c((MPI_Fint)datatype), &petsc_send_len, &petsc_send_len_th));
2566:   #endif
2567: #endif
2568:   PetscFunctionReturn(PETSC_SUCCESS);
2569: }

2571: PETSC_EXTERN PetscErrorCode PetscARecv(int count, int datatype)
2572: {
2573:   PetscFunctionBegin;
2574: #if PetscDefined(USE_LOG)
2575:   PetscCall(PetscAddLogDouble(&petsc_recv_ct, &petsc_recv_ct_th, 1));
2576:   #if !defined(MPIUNI_H) && !defined(PETSC_HAVE_BROKEN_RECURSIVE_MACRO)
2577:   PetscCall(PetscMPITypeSize(count, MPI_Type_f2c((MPI_Fint)datatype), &petsc_recv_len, &petsc_recv_len_th));
2578:   #endif
2579: #endif
2580:   PetscFunctionReturn(PETSC_SUCCESS);
2581: }

2583: PETSC_EXTERN PetscErrorCode PetscAReduce(void)
2584: {
2585:   PetscFunctionBegin;
2586:   if (PetscDefined(USE_LOG)) PetscCall(PetscAddLogDouble(&petsc_allreduce_ct, &petsc_allreduce_ct_th, 1));
2587:   PetscFunctionReturn(PETSC_SUCCESS);
2588: }

2590: PetscClassId PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID;
2591: PetscClassId PETSC_OBJECT_CLASSID  = 0;

2593: static PetscBool PetscLogInitializeCalled = PETSC_FALSE;

2595: PETSC_INTERN PetscErrorCode PetscLogInitialize(void)
2596: {
2597:   int stage;

2599:   PetscFunctionBegin;
2600:   if (PetscLogInitializeCalled) PetscFunctionReturn(PETSC_SUCCESS);
2601:   PetscLogInitializeCalled = PETSC_TRUE;
2602:   if (PetscDefined(USE_LOG)) {
2603:     /* Setup default logging structures */
2604:     PetscCall(PetscLogStateCreate(&petsc_log_state));
2605:     for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
2606:       if (PetscLogHandlers[i].handler) PetscCall(PetscLogHandlerSetState(PetscLogHandlers[i].handler, petsc_log_state));
2607:     }
2608:     PetscCall(PetscLogStateStageRegister(petsc_log_state, "Main Stage", &stage));
2609:     PetscCall(PetscSpinlockCreate(&PetscLogSpinLock));
2610: #if defined(PETSC_HAVE_THREADSAFETY)
2611:     petsc_log_tid = 0;
2612:     petsc_log_gid = 0;
2613: #endif

2615:     /* All processors sync here for more consistent logging */
2616:     PetscCallMPI(MPI_Barrier(PETSC_COMM_WORLD));
2617:     PetscCall(PetscTime(&petsc_BaseTime));
2618:     PetscCall(PetscLogStagePush(stage));
2619:   }
2620:   PetscFunctionReturn(PETSC_SUCCESS);
2621: }

2623: PETSC_INTERN PetscErrorCode PetscLogFinalize(void)
2624: {
2625:   PetscFunctionBegin;
2626:   if (PetscDefined(USE_LOG)) {
2627:     /* Resetting phase */
2628:     // pop remaining stages
2629:     if (petsc_log_state) {
2630:       while (petsc_log_state->current_stage >= 0) PetscCall(PetscLogStagePop());
2631:     }
2632:     for (int i = 0; i < PETSC_LOG_HANDLER_MAX; i++) PetscCall(PetscLogHandlerDestroy(&PetscLogHandlers[i].handler));
2633:     PetscCall(PetscArrayzero(PetscLogHandlers, PETSC_LOG_HANDLER_MAX));
2634:     PetscCall(PetscLogStateDestroy(&petsc_log_state));

2636:     petsc_TotalFlops         = 0.0;
2637:     petsc_BaseTime           = 0.0;
2638:     petsc_TotalFlops         = 0.0;
2639:     petsc_send_ct            = 0.0;
2640:     petsc_recv_ct            = 0.0;
2641:     petsc_send_len           = 0.0;
2642:     petsc_recv_len           = 0.0;
2643:     petsc_isend_ct           = 0.0;
2644:     petsc_irecv_ct           = 0.0;
2645:     petsc_isend_len          = 0.0;
2646:     petsc_irecv_len          = 0.0;
2647:     petsc_wait_ct            = 0.0;
2648:     petsc_wait_any_ct        = 0.0;
2649:     petsc_wait_all_ct        = 0.0;
2650:     petsc_sum_of_waits_ct    = 0.0;
2651:     petsc_allreduce_ct       = 0.0;
2652:     petsc_gather_ct          = 0.0;
2653:     petsc_scatter_ct         = 0.0;
2654:     petsc_TotalFlops_th      = 0.0;
2655:     petsc_send_ct_th         = 0.0;
2656:     petsc_recv_ct_th         = 0.0;
2657:     petsc_send_len_th        = 0.0;
2658:     petsc_recv_len_th        = 0.0;
2659:     petsc_isend_ct_th        = 0.0;
2660:     petsc_irecv_ct_th        = 0.0;
2661:     petsc_isend_len_th       = 0.0;
2662:     petsc_irecv_len_th       = 0.0;
2663:     petsc_wait_ct_th         = 0.0;
2664:     petsc_wait_any_ct_th     = 0.0;
2665:     petsc_wait_all_ct_th     = 0.0;
2666:     petsc_sum_of_waits_ct_th = 0.0;
2667:     petsc_allreduce_ct_th    = 0.0;
2668:     petsc_gather_ct_th       = 0.0;
2669:     petsc_scatter_ct_th      = 0.0;

2671:     petsc_ctog_ct       = 0.0;
2672:     petsc_gtoc_ct       = 0.0;
2673:     petsc_ctog_sz       = 0.0;
2674:     petsc_gtoc_sz       = 0.0;
2675:     petsc_gflops        = 0.0;
2676:     petsc_gtime         = 0.0;
2677:     petsc_genergy       = 0.0;
2678:     petsc_genergy_meter = 0.0;
2679:     petsc_ctog_ct_th    = 0.0;
2680:     petsc_gtoc_ct_th    = 0.0;
2681:     petsc_ctog_sz_th    = 0.0;
2682:     petsc_gtoc_sz_th    = 0.0;
2683:     petsc_gflops_th     = 0.0;
2684:     petsc_gtime_th      = 0.0;
2685:   }
2686:   PETSC_LARGEST_CLASSID    = PETSC_SMALLEST_CLASSID;
2687:   PETSC_OBJECT_CLASSID     = 0;
2688:   PetscLogInitializeCalled = PETSC_FALSE;
2689:   PetscFunctionReturn(PETSC_SUCCESS);
2690: }

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

2695:   Not Collective

2697:   Input Parameter:
2698: . name - The class name

2700:   Output Parameter:
2701: . oclass - The class id or classid

2703:   Level: developer

2705: .seealso: [](ch_profiling), `PetscLogEventRegister()`
2706: @*/
2707: PetscErrorCode PetscClassIdRegister(const char name[], PetscClassId *oclass)
2708: {
2709:   PetscFunctionBegin;
2710:   *oclass = ++PETSC_LARGEST_CLASSID;
2711: #if defined(PETSC_USE_LOG)
2712:   {
2713:     PetscLogState state;
2714:     PetscLogClass logclass;

2716:     PetscCall(PetscLogGetState(&state));
2717:     if (state) PetscCall(PetscLogStateClassRegister(state, name, *oclass, &logclass));
2718:   }
2719: #endif
2720:   PetscFunctionReturn(PETSC_SUCCESS);
2721: }