xref: /petsc/src/sys/logging/plog.c (revision 390474f96c2cfb142235edf4f816cb7a2ce79c2a)
1 /*
2       PETSc code to log object creation and destruction and PETSc events.
3 
4       This provides the public API used by the rest of PETSc and by users.
5 
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.
8 
9       ***
10 
11       This file, and only this file, is for functions that interact with the global logging state
12 */
13 #include <petsc/private/logimpl.h> /*I    "petscsys.h"   I*/
14 #include <petsc/private/loghandlerimpl.h>
15 #include <petsctime.h>
16 #include <petscviewer.h>
17 #include <petscdevice.h>
18 #include <petsc/private/deviceimpl.h>
19 
20 #if defined(PETSC_HAVE_THREADSAFETY)
21 
22 PetscInt           petsc_log_gid = -1; /* Global threadId counter */
23 PETSC_TLS PetscInt petsc_log_tid = -1; /* Local threadId */
24 
25 /* shared variables */
26 PetscSpinlock PetscLogSpinLock;
27 
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 }
37 
38 #endif
39 
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 */
58 
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;
76 
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 */
89 
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;
100 
101 PetscBool PetscLogMemory = PETSC_FALSE;
102 PetscBool PetscLogSyncOn = PETSC_FALSE;
103 
104 PetscBool PetscLogGpuTimeFlag        = PETSC_FALSE;
105 PetscBool PetscLogGpuEnergyFlag      = PETSC_FALSE;
106 PetscBool PetscLogGpuEnergyMeterFlag = PETSC_FALSE;
107 
108 PetscInt PetscLogNumViewersCreated   = 0;
109 PetscInt PetscLogNumViewersDestroyed = 0;
110 
111 PetscLogState petsc_log_state = NULL;
112 
113 #define PETSC_LOG_HANDLER_HOT_BLANK {NULL, NULL, NULL, NULL, NULL, NULL}
114 
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 };
121 
122 #undef PETSC_LOG_HANDLERS_HOT_BLANK
123 
124 #if defined(PETSC_USE_LOG)
125   #include <../src/sys/logging/handler/impls/default/logdefault.h>
126 
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 }
136 
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 }
147 
148   #endif
149 
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;
159 
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 }
169 
170 /*@
171   PetscLogGetDefaultHandler - Get the default log handler if it is running.
172 
173   Not collective
174 
175   Output Parameter:
176 . handler - the default `PetscLogHandler`, or `NULL` if it is not running.
177 
178   Level: developer
179 
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`.
184 
185 .seealso: [](ch_profiling)
186 @*/
187 PetscErrorCode PetscLogGetDefaultHandler(PetscLogHandler *handler)
188 {
189   PetscFunctionBegin;
190   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, handler));
191   PetscFunctionReturn(PETSC_SUCCESS);
192 }
193 
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 }
202 
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()`).
208 
209   Collective on `PETSC_COMM_WORLD`
210 
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
215 
216   Level: developer
217 
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 }
227 
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 }
239 
240 /*@
241   PetscLogHandlerStart - Connect a log handler to PETSc's global logging stream and state.
242 
243   Logically collective
244 
245   Input Parameters:
246 . h - a `PetscLogHandler`
247 
248   Level: developer
249 
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.
254 
255   There is a limit of `PESC_LOG_HANDLER_MAX` handlers that can be active at one time.
256 
257   To disconnect a handler from the global stream call `PetscLogHandlerStop()`.
258 
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.
263 
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;
279 
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 }
301 
302 /*@
303   PetscLogHandlerStop - Disconnect a log handler from PETSc's global logging stream.
304 
305   Logically collective
306 
307   Input Parameters:
308 . h - a `PetscLogHandler`
309 
310   Level: developer
311 
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()`),
316 
317   When a log handler is stopped, the remaining stages will be popped before it is
318   disconnected from the log stream.
319 
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;
331 
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];
340 
341           PetscCall(PetscIntStackPush(temp_stack, stage));
342         }
343         for (int s = 0; s < stack_height; s++) {
344           PetscLogStage stage;
345           PetscBool     empty;
346 
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 }
363 
364 /*@
365   PetscLogIsActive - Check if logging (profiling) is currently in progress.
366 
367   Not Collective
368 
369   Output Parameter:
370 . isActive - `PETSC_TRUE` if logging is in progress, `PETSC_FALSE` otherwise
371 
372   Level: beginner
373 
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 }
390 
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 }
405 
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 }
420 
421 PETSC_INTERN PetscErrorCode PetscLogTypeBegin(PetscLogHandlerType type)
422 {
423   PetscLogHandler handler;
424 
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 }
434 
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.
438 
439   Logically Collective on `PETSC_COMM_WORLD`
440 
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()`.
445 
446   Example Usage:
447 .vb
448       PetscInitialize(...);
449       PetscLogDefaultBegin();
450        ... code ...
451       PetscLogView(viewer); or PetscLogDump();
452       PetscFinalize();
453 .ve
454 
455   Level: advanced
456 
457   Notes:
458   `PetscLogView()` or `PetscLogDump()` actually cause the printing of
459   the logging information.
460 
461   This routine may be called more than once.
462 
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()`
465 
466 .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogView()`, `PetscLogTraceBegin()`
467 @*/
468 PetscErrorCode PetscLogDefaultBegin(void)
469 {
470   PetscFunctionBegin;
471   PetscCall(PetscLogTypeBegin(PETSCLOGHANDLERDEFAULT));
472   PetscFunctionReturn(PETSC_SUCCESS);
473 }
474 
475 /*@C
476   PetscLogTraceBegin - Begins trace logging.  Every time a PETSc event
477   begins or ends, the event name is printed.
478 
479   Logically Collective on `PETSC_COMM_WORLD`, No Fortran Support
480 
481   Input Parameter:
482 . file - The file to print trace in (e.g. stdout)
483 
484   Options Database Key:
485 . -log_trace [filename] - Begins `PetscLogTraceBegin()`
486 
487   Level: intermediate
488 
489   Notes:
490   `PetscLogTraceBegin()` prints the processor number, the execution time (sec),
491   then "Event begin:" or "Event end:" followed by the event name.
492 
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.
496 
497 .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogView()`, `PetscLogDefaultBegin()`
498 @*/
499 PetscErrorCode PetscLogTraceBegin(FILE *file)
500 {
501   PetscLogHandler handler;
502 
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 }
511 
512 PETSC_INTERN PetscErrorCode PetscLogHandlerCreate_Nested(MPI_Comm, PetscLogHandler *);
513 
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.
517 
518   Logically Collective on `PETSC_COMM_WORLD`, No Fortran Support
519 
520   Options Database Keys:
521 . -log_view :filename.xml:ascii_xml - Prints an XML summary of flop and timing information to the file
522 
523   Example Usage:
524 .vb
525       PetscInitialize(...);
526       PetscLogNestedBegin();
527        ... code ...
528       PetscLogView(viewer);
529       PetscFinalize();
530 .ve
531 
532   Level: advanced
533 
534 .seealso: `PetscLogDump()`, `PetscLogView()`, `PetscLogTraceBegin()`, `PetscLogDefaultBegin()`
535 @*/
536 PetscErrorCode PetscLogNestedBegin(void)
537 {
538   PetscFunctionBegin;
539   PetscCall(PetscLogTypeBegin(PETSCLOGHANDLERNESTED));
540   PetscFunctionReturn(PETSC_SUCCESS);
541 }
542 
543 /*@C
544   PetscLogLegacyCallbacksBegin - Create and start a log handler from callbacks
545   matching the now deprecated function pointers `PetscLogPLB`, `PetscLogPLE`,
546   `PetscLogPHC`, `PetscLogPHD`.
547 
548   Logically Collective on `PETSC_COMM_WORLD`
549 
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`)
555 
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`)
563 
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`)
571 
572   Calling sequence of `PetscLogPHC`:
573 . o - a `PetscObject` that has just been created
574 
575   Calling sequence of `PetscLogPHD`:
576 . o - a `PetscObject` that is about to be destroyed
577 
578   Level: advanced
579 
580   Notes:
581   This is for transitioning from the deprecated function `PetscLogSet()` and should not be used in new code.
582 
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.
586 
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;
592 
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 }
599 
600   #if defined(PETSC_HAVE_MPE)
601     #include <mpe.h>
602 static PetscBool PetscBeganMPE = PETSC_FALSE;
603   #endif
604 
605 /*@C
606   PetscLogMPEBegin - Turns on MPE logging of events. This creates large log files and slows the
607   program down.
608 
609   Collective on `PETSC_COMM_WORLD`, No Fortran Support
610 
611   Options Database Key:
612 . -log_mpe - Prints extensive log information
613 
614   Level: advanced
615 
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).
620 
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());
632 
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 }
643 
644   #if defined(PETSC_HAVE_TAU_PERFSTUBS)
645     #include <../src/sys/perfstubs/timer.h>
646   #endif
647 
648 /*@C
649   PetscLogPerfstubsBegin - Turns on logging of events using the perfstubs interface.
650 
651   Collective on `PETSC_COMM_WORLD`, No Fortran Support
652 
653   Options Database Key:
654 . -log_perfstubs - use an external log handler through the perfstubs interface
655 
656   Level: advanced
657 
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 }
670 
671 /*@
672   PetscLogActions - Determines whether actions are logged for the default log handler.
673 
674   Not Collective
675 
676   Input Parameter:
677 . flag - `PETSC_TRUE` if actions are to be logged
678 
679   Options Database Key:
680 + -log_exclude_actions - (deprecated) Does nothing
681 - -log_include_actions - Turn on action logging
682 
683   Level: intermediate
684 
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.
688 
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;
696 
697     if (h) PetscCall(PetscLogHandlerSetLogActions(h, flag));
698   }
699   PetscFunctionReturn(PETSC_SUCCESS);
700 }
701 
702 /*@
703   PetscLogObjects - Determines whether objects are logged for the graphical viewer.
704 
705   Not Collective
706 
707   Input Parameter:
708 . flag - `PETSC_TRUE` if objects are to be logged
709 
710   Options Database Key:
711 + -log_exclude_objects - (deprecated) Does nothing
712 - -log_include_objects - Turns on object logging
713 
714   Level: intermediate
715 
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.
719 
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;
727 
728     if (h) PetscCall(PetscLogHandlerSetLogObjects(h, flag));
729   }
730   PetscFunctionReturn(PETSC_SUCCESS);
731 }
732 
733 /*------------------------------------------------ Stage Functions --------------------------------------------------*/
734 /*@
735   PetscLogStageRegister - Attaches a character string name to a logging stage.
736 
737   Not Collective
738 
739   Input Parameter:
740 . sname - The name to associate with that stage
741 
742   Output Parameter:
743 . stage - The stage number or -1 if logging is not active (`PetscLogIsActive()`).
744 
745   Level: intermediate
746 
747 .seealso: [](ch_profiling), `PetscLogStagePush()`, `PetscLogStagePop()`
748 @*/
749 PetscErrorCode PetscLogStageRegister(const char sname[], PetscLogStage *stage)
750 {
751   PetscLogState state;
752 
753   PetscFunctionBegin;
754   *stage = -1;
755   PetscCall(PetscLogGetState(&state));
756   if (state) PetscCall(PetscLogStateStageRegister(state, sname, stage));
757   PetscFunctionReturn(PETSC_SUCCESS);
758 }
759 
760 /*@
761   PetscLogStagePush - This function pushes a stage on the logging stack. Events started and stopped until `PetscLogStagePop()` will be associated with the stage
762 
763   Not Collective
764 
765   Input Parameter:
766 . stage - The stage on which to log
767 
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
782 
783   Level: intermediate
784 
785   Note:
786   Use `PetscLogStageRegister()` to register a stage.
787 
788 .seealso: [](ch_profiling), `PetscLogStagePop()`, `PetscLogStageRegister()`, `PetscBarrier()`
789 @*/
790 PetscErrorCode PetscLogStagePush(PetscLogStage stage)
791 {
792   PetscLogState state;
793 
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 }
804 
805 /*@
806   PetscLogStagePop - This function pops a stage from the logging stack that was pushed with `PetscLogStagePush()`
807 
808   Not Collective
809 
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
824 
825   Level: intermediate
826 
827 .seealso: [](ch_profiling), `PetscLogStagePush()`, `PetscLogStageRegister()`, `PetscBarrier()`
828 @*/
829 PetscErrorCode PetscLogStagePop(void)
830 {
831   PetscLogState state;
832   PetscLogStage current_stage;
833 
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 }
845 
846 /*@
847   PetscLogStageSetActive - Sets if a stage is used for `PetscLogEventBegin()` and `PetscLogEventEnd()`.
848 
849   Not Collective
850 
851   Input Parameters:
852 + stage    - The stage
853 - isActive - The activity flag, `PETSC_TRUE` for logging, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)
854 
855   Level: intermediate
856 
857   Note:
858   If this is set to `PETSC_FALSE` the logging acts as if the stage did not exist
859 
860 .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
861 @*/
862 PetscErrorCode PetscLogStageSetActive(PetscLogStage stage, PetscBool isActive)
863 {
864   PetscLogState state;
865 
866   PetscFunctionBegin;
867   PetscCall(PetscLogGetState(&state));
868   if (state) PetscCall(PetscLogStateStageSetActive(state, stage, isActive));
869   PetscFunctionReturn(PETSC_SUCCESS);
870 }
871 
872 /*@
873   PetscLogStageGetActive - Checks if a stage is used for `PetscLogEventBegin()` and `PetscLogEventEnd()`.
874 
875   Not Collective
876 
877   Input Parameter:
878 . stage - The stage
879 
880   Output Parameter:
881 . isActive - The activity flag, `PETSC_TRUE` for logging, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)
882 
883   Level: intermediate
884 
885 .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
886 @*/
887 PetscErrorCode PetscLogStageGetActive(PetscLogStage stage, PetscBool *isActive)
888 {
889   PetscLogState state;
890 
891   PetscFunctionBegin;
892   *isActive = PETSC_FALSE;
893   PetscCall(PetscLogGetState(&state));
894   if (state) PetscCall(PetscLogStateStageGetActive(state, stage, isActive));
895   PetscFunctionReturn(PETSC_SUCCESS);
896 }
897 
898 /*@
899   PetscLogStageSetVisible - Determines stage visibility in `PetscLogView()`
900 
901   Not Collective
902 
903   Input Parameters:
904 + stage     - The stage
905 - isVisible - The visibility flag, `PETSC_TRUE` to print, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)
906 
907   Level: intermediate
908 
909   Developer Notes:
910   Visibility only affects the default log handler in `PetscLogView()`: stages that are
911   set to invisible are suppressed from output.
912 
913 .seealso: [](ch_profiling), `PetscLogStageGetVisible()`, `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogView()`, `PetscLogGetDefaultHandler()`
914 @*/
915 PetscErrorCode PetscLogStageSetVisible(PetscLogStage stage, PetscBool isVisible)
916 
917 {
918   PetscFunctionBegin;
919   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
920     PetscLogHandler h = PetscLogHandlers[i].handler;
921 
922     if (h) PetscCall(PetscLogHandlerStageSetVisible(h, stage, isVisible));
923   }
924   PetscFunctionReturn(PETSC_SUCCESS);
925 }
926 
927 /*@
928   PetscLogStageGetVisible - Returns stage visibility in `PetscLogView()`
929 
930   Not Collective
931 
932   Input Parameter:
933 . stage - The stage
934 
935   Output Parameter:
936 . isVisible - The visibility flag, `PETSC_TRUE` to print, else `PETSC_FALSE` (defaults to `PETSC_TRUE`)
937 
938   Level: intermediate
939 
940 .seealso: [](ch_profiling), `PetscLogStageSetVisible()`, `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscLogView()`, `PetscLogGetDefaultHandler()`
941 @*/
942 PetscErrorCode PetscLogStageGetVisible(PetscLogStage stage, PetscBool *isVisible)
943 {
944   PetscLogHandler handler;
945 
946   PetscFunctionBegin;
947   *isVisible = PETSC_FALSE;
948   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
949   if (handler) PetscCall(PetscLogHandlerStageGetVisible(handler, stage, isVisible));
950   PetscFunctionReturn(PETSC_SUCCESS);
951 }
952 
953 /*@
954   PetscLogStageGetId - Returns the stage id when given the stage name.
955 
956   Not Collective
957 
958   Input Parameter:
959 . name - The stage name
960 
961   Output Parameter:
962 . stage - The stage, , or -1 if no stage with that name exists
963 
964   Level: intermediate
965 
966 .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
967 @*/
968 PetscErrorCode PetscLogStageGetId(const char name[], PetscLogStage *stage)
969 {
970   PetscLogState state;
971 
972   PetscFunctionBegin;
973   *stage = -1;
974   PetscCall(PetscLogGetState(&state));
975   if (state) PetscCall(PetscLogStateGetStageFromName(state, name, stage));
976   PetscFunctionReturn(PETSC_SUCCESS);
977 }
978 
979 /*@
980   PetscLogStageGetName - Returns the stage name when given the stage id.
981 
982   Not Collective
983 
984   Input Parameter:
985 . stage - The stage
986 
987   Output Parameter:
988 . name - The stage name
989 
990   Level: intermediate
991 
992 .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogStagePush()`, `PetscLogStagePop()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
993 @*/
994 PetscErrorCode PetscLogStageGetName(PetscLogStage stage, const char *name[])
995 {
996   PetscLogStageInfo stage_info;
997   PetscLogState     state;
998 
999   PetscFunctionBegin;
1000   *name = NULL;
1001   PetscCall(PetscLogGetState(&state));
1002   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
1003   PetscCall(PetscLogStateStageGetInfo(state, stage, &stage_info));
1004   *name = stage_info.name;
1005   PetscFunctionReturn(PETSC_SUCCESS);
1006 }
1007 
1008 /*------------------------------------------------ Event Functions --------------------------------------------------*/
1009 
1010 /*@
1011   PetscLogEventRegister - Registers an event name for logging operations
1012 
1013   Not Collective
1014 
1015   Input Parameters:
1016 + name    - The name associated with the event
1017 - classid - The classid associated to the class for this event, obtain either with
1018            `PetscClassIdRegister()` or use a predefined one such as `KSP_CLASSID`, `SNES_CLASSID`, the predefined ones
1019            are only available in C code
1020 
1021   Output Parameter:
1022 . event - The event id for use with `PetscLogEventBegin()` and `PetscLogEventEnd()`.
1023 
1024   Example Usage:
1025 .vb
1026       PetscLogEvent USER_EVENT;
1027       PetscClassId classid;
1028       PetscLogDouble user_event_flops;
1029       PetscClassIdRegister("class name",&classid);
1030       PetscLogEventRegister("User event name",classid,&USER_EVENT);
1031       PetscLogEventBegin(USER_EVENT,0,0,0,0);
1032          [code segment to monitor]
1033          PetscLogFlops(user_event_flops);
1034       PetscLogEventEnd(USER_EVENT,0,0,0,0);
1035 .ve
1036 
1037   Level: intermediate
1038 
1039   Notes:
1040   PETSc automatically logs library events if the code has been
1041   configured with --with-log (which is the default) and
1042   -log_view or -log_all is specified.  `PetscLogEventRegister()` is
1043   intended for logging user events to supplement this PETSc
1044   information.
1045 
1046   PETSc can gather data for use with the utilities Jumpshot
1047   (part of the MPICH distribution).  If PETSc has been compiled
1048   with flag -DPETSC_HAVE_MPE (MPE is an additional utility within
1049   MPICH), the user can employ another command line option, -log_mpe,
1050   to create a logfile, "mpe.log", which can be visualized
1051   Jumpshot.
1052 
1053   The classid is associated with each event so that classes of events
1054   can be disabled simultaneously, such as all matrix events. The user
1055   can either use an existing classid, such as `MAT_CLASSID`, or create
1056   their own as shown in the example.
1057 
1058   If an existing event with the same name exists, its event handle is
1059   returned instead of creating a new event.
1060 
1061 .seealso: [](ch_profiling), `PetscLogStageRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogFlops()`,
1062           `PetscLogEventActivate()`, `PetscLogEventDeactivate()`, `PetscClassIdRegister()`
1063 @*/
1064 PetscErrorCode PetscLogEventRegister(const char name[], PetscClassId classid, PetscLogEvent *event)
1065 {
1066   PetscLogState state;
1067 
1068   PetscFunctionBegin;
1069   *event = -1;
1070   PetscCall(PetscLogGetState(&state));
1071   if (state) PetscCall(PetscLogStateEventRegister(state, name, classid, event));
1072   PetscFunctionReturn(PETSC_SUCCESS);
1073 }
1074 
1075 /*@
1076   PetscLogEventSetCollective - Indicates that a particular event is collective.
1077 
1078   Logically Collective
1079 
1080   Input Parameters:
1081 + event      - The event id
1082 - collective - `PetscBool` indicating whether a particular event is collective
1083 
1084   Level: developer
1085 
1086   Notes:
1087   New events returned from `PetscLogEventRegister()` are collective by default.
1088 
1089   Collective events are handled specially if the command line option `-log_sync` is used. In that case the logging saves information about
1090   two parts of the event; the time for all the MPI ranks to synchronize and then the time for the actual computation/communication
1091   to be performed. This option is useful to debug imbalance within the computations or communications.
1092 
1093 .seealso: [](ch_profiling), `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogEventRegister()`
1094 @*/
1095 PetscErrorCode PetscLogEventSetCollective(PetscLogEvent event, PetscBool collective)
1096 {
1097   PetscLogState state;
1098 
1099   PetscFunctionBegin;
1100   PetscCall(PetscLogGetState(&state));
1101   if (state) PetscCall(PetscLogStateEventSetCollective(state, event, collective));
1102   PetscFunctionReturn(PETSC_SUCCESS);
1103 }
1104 
1105 /*
1106   PetscLogClassSetActiveAll - Activate or inactivate logging for all events associated with a PETSc object class in every stage.
1107 
1108   Not Collective
1109 
1110   Input Parameters:
1111 + classid - The object class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1112 - isActive - if `PETSC_FALSE`, events associated with this class will not be send to log handlers.
1113 
1114   Level: developer
1115 
1116 .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventActivateAll()`, `PetscLogStageSetActive()`, `PetscLogEventActivateClass()`
1117 */
1118 static PetscErrorCode PetscLogClassSetActiveAll(PetscClassId classid, PetscBool isActive)
1119 {
1120   PetscLogState state;
1121 
1122   PetscFunctionBegin;
1123   PetscCall(PetscLogGetState(&state));
1124   if (state) PetscCall(PetscLogStateClassSetActiveAll(state, classid, isActive));
1125   PetscFunctionReturn(PETSC_SUCCESS);
1126 }
1127 
1128 /*@
1129   PetscLogEventIncludeClass - Activates event logging for a PETSc object class in every stage.
1130 
1131   Not Collective
1132 
1133   Input Parameter:
1134 . classid - The object class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1135 
1136   Level: developer
1137 
1138 .seealso: [](ch_profiling), `PetscLogEventActivateClass()`, `PetscLogEventDeactivateClass()`, `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1139 @*/
1140 PetscErrorCode PetscLogEventIncludeClass(PetscClassId classid)
1141 {
1142   PetscFunctionBegin;
1143   PetscCall(PetscLogClassSetActiveAll(classid, PETSC_TRUE));
1144   PetscFunctionReturn(PETSC_SUCCESS);
1145 }
1146 
1147 /*@
1148   PetscLogEventExcludeClass - Deactivates event logging for a PETSc object class in every stage.
1149 
1150   Not Collective
1151 
1152   Input Parameter:
1153 . classid - The object class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1154 
1155   Level: developer
1156 
1157   Note:
1158   If a class is excluded then events associated with that class are not logged.
1159 
1160 .seealso: [](ch_profiling), `PetscLogEventDeactivateClass()`, `PetscLogEventActivateClass()`, `PetscLogEventDeactivate()`, `PetscLogEventActivate()`
1161 @*/
1162 PetscErrorCode PetscLogEventExcludeClass(PetscClassId classid)
1163 {
1164   PetscFunctionBegin;
1165   PetscCall(PetscLogClassSetActiveAll(classid, PETSC_FALSE));
1166   PetscFunctionReturn(PETSC_SUCCESS);
1167 }
1168 
1169 /*
1170   PetscLogEventSetActive - Activate or inactivate logging for an event in a given stage
1171 
1172   Not Collective
1173 
1174   Input Parameters:
1175 + stage - A registered `PetscLogStage` (or `PETSC_DEFAULT` for the current stage)
1176 . event - A `PetscLogEvent`
1177 - isActive - If `PETSC_FALSE`, activity from this event (`PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogEventSync()`) will not be sent to log handlers during this stage
1178 
1179   Usage:
1180 .vb
1181       PetscLogEventSetActive(VEC_SetValues, PETSC_FALSE);
1182         [code where you do not want to log VecSetValues()]
1183       PetscLogEventSetActive(VEC_SetValues, PETSC_TRUE);
1184         [code where you do want to log VecSetValues()]
1185 .ve
1186 
1187   Level: advanced
1188 
1189   Note:
1190   The event may be either a pre-defined PETSc event (found in include/petsclog.h)
1191   or an event number obtained with `PetscLogEventRegister()`.
1192 
1193 .seealso: [](ch_profiling), `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`
1194 */
1195 static PetscErrorCode PetscLogEventSetActive(PetscLogStage stage, PetscLogEvent event, PetscBool isActive)
1196 {
1197   PetscLogState state;
1198 
1199   PetscFunctionBegin;
1200   PetscCall(PetscLogGetState(&state));
1201   if (state) PetscCall(PetscLogStateEventSetActive(state, stage, event, isActive));
1202   PetscFunctionReturn(PETSC_SUCCESS);
1203 }
1204 
1205 /*@
1206   PetscLogEventActivate - Indicates that a particular event should be logged.
1207 
1208   Not Collective
1209 
1210   Input Parameter:
1211 . event - The event id
1212 
1213   Example Usage:
1214 .vb
1215       PetscLogEventDeactivate(VEC_SetValues);
1216         [code where you do not want to log VecSetValues()]
1217       PetscLogEventActivate(VEC_SetValues);
1218         [code where you do want to log VecSetValues()]
1219 .ve
1220 
1221   Level: advanced
1222 
1223   Note:
1224   The event may be either a pre-defined PETSc event (found in include/petsclog.h)
1225   or an event number obtained with `PetscLogEventRegister()`.
1226 
1227 .seealso: [](ch_profiling), `PetscLogEventDeactivate()`, `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`
1228 @*/
1229 PetscErrorCode PetscLogEventActivate(PetscLogEvent event)
1230 {
1231   PetscFunctionBegin;
1232   PetscCall(PetscLogEventSetActive(PETSC_DEFAULT, event, PETSC_TRUE));
1233   PetscFunctionReturn(PETSC_SUCCESS);
1234 }
1235 
1236 /*@
1237   PetscLogEventDeactivate - Indicates that a particular event should not be logged.
1238 
1239   Not Collective
1240 
1241   Input Parameter:
1242 . event - The event id
1243 
1244   Example Usage:
1245 .vb
1246       PetscLogEventDeactivate(VEC_SetValues);
1247         [code where you do not want to log VecSetValues()]
1248       PetscLogEventActivate(VEC_SetValues);
1249         [code where you do want to log VecSetValues()]
1250 .ve
1251 
1252   Level: advanced
1253 
1254   Note:
1255   The event may be either a pre-defined PETSc event (found in
1256   include/petsclog.h) or an event number obtained with `PetscLogEventRegister()`).
1257 
1258 .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`
1259 @*/
1260 PetscErrorCode PetscLogEventDeactivate(PetscLogEvent event)
1261 {
1262   PetscFunctionBegin;
1263   PetscCall(PetscLogEventSetActive(PETSC_DEFAULT, event, PETSC_FALSE));
1264   PetscFunctionReturn(PETSC_SUCCESS);
1265 }
1266 
1267 /*@
1268   PetscLogEventDeactivatePush - Indicates that a particular event should not be logged until `PetscLogEventDeactivatePop()` is called
1269 
1270   Not Collective
1271 
1272   Input Parameter:
1273 . event - The event id
1274 
1275   Example Usage:
1276 .vb
1277       PetscLogEventDeactivatePush(VEC_SetValues);
1278         [code where you do not want to log VecSetValues()]
1279       PetscLogEventDeactivatePop(VEC_SetValues);
1280         [code where you do want to log VecSetValues()]
1281 .ve
1282 
1283   Level: advanced
1284 
1285   Note:
1286   The event may be either a pre-defined PETSc event (found in
1287   include/petsclog.h) or an event number obtained with `PetscLogEventRegister()`).
1288 
1289   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 
1291 .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivate()`, `PetscLogEventDeactivatePop()`
1292 @*/
1293 PetscErrorCode PetscLogEventDeactivatePush(PetscLogEvent event)
1294 {
1295   PetscFunctionBegin;
1296   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1297     PetscLogHandler h = PetscLogHandlers[i].handler;
1298 
1299     if (h) PetscCall(PetscLogHandlerEventDeactivatePush(h, PETSC_DEFAULT, event));
1300   }
1301   PetscFunctionReturn(PETSC_SUCCESS);
1302 }
1303 
1304 /*@
1305   PetscLogEventDeactivatePop - Indicates that a particular event should again be logged after the logging was turned off with `PetscLogEventDeactivatePush()`
1306 
1307   Not Collective
1308 
1309   Input Parameter:
1310 . event - The event id
1311 
1312   Example Usage:
1313 .vb
1314       PetscLogEventDeactivatePush(VEC_SetValues);
1315         [code where you do not want to log VecSetValues()]
1316       PetscLogEventDeactivatePop(VEC_SetValues);
1317         [code where you do want to log VecSetValues()]
1318 .ve
1319 
1320   Level: advanced
1321 
1322   Note:
1323   The event may be either a pre-defined PETSc event (found in
1324   include/petsclog.h) or an event number obtained with `PetscLogEventRegister()`).
1325 
1326 .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivatePush()`
1327 @*/
1328 PetscErrorCode PetscLogEventDeactivatePop(PetscLogEvent event)
1329 {
1330   PetscFunctionBegin;
1331   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1332     PetscLogHandler h = PetscLogHandlers[i].handler;
1333 
1334     if (h) PetscCall(PetscLogHandlerEventDeactivatePop(h, PETSC_DEFAULT, event));
1335   }
1336   PetscFunctionReturn(PETSC_SUCCESS);
1337 }
1338 
1339 /*@
1340   PetscLogEventSetActiveAll - Turns on logging of all events
1341 
1342   Not Collective
1343 
1344   Input Parameters:
1345 + event    - The event id
1346 - isActive - The activity flag determining whether the event is logged
1347 
1348   Level: advanced
1349 
1350 .seealso: [](ch_profiling), `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1351 @*/
1352 PetscErrorCode PetscLogEventSetActiveAll(PetscLogEvent event, PetscBool isActive)
1353 {
1354   PetscLogState state;
1355 
1356   PetscFunctionBegin;
1357   PetscCall(PetscLogGetState(&state));
1358   if (state) PetscCall(PetscLogStateEventSetActiveAll(state, event, isActive));
1359   PetscFunctionReturn(PETSC_SUCCESS);
1360 }
1361 
1362 /*
1363   PetscLogClassSetActive - Activates event logging for a PETSc object class for the current stage
1364 
1365   Not Collective
1366 
1367   Input Parameters:
1368 + stage - A registered `PetscLogStage` (or `PETSC_DEFAULT` for the current stage)
1369 . classid - The event class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1370 - isActive - If `PETSC_FALSE`, events associated with this class are not sent to log handlers.
1371 
1372   Level: developer
1373 
1374 .seealso: [](ch_profiling), `PetscLogEventIncludeClass()`, `PetscLogEventActivate()`, `PetscLogEventActivateAll()`, `PetscLogStageSetActive()`
1375 */
1376 static PetscErrorCode PetscLogClassSetActive(PetscLogStage stage, PetscClassId classid, PetscBool isActive)
1377 {
1378   PetscLogState state;
1379 
1380   PetscFunctionBegin;
1381   PetscCall(PetscLogGetState(&state));
1382   if (state) PetscCall(PetscLogStateClassSetActive(state, stage, classid, isActive));
1383   PetscFunctionReturn(PETSC_SUCCESS);
1384 }
1385 
1386 /*@
1387   PetscLogEventActivateClass - Activates event logging for a PETSc object class for the current stage
1388 
1389   Not Collective
1390 
1391   Input Parameter:
1392 . classid - The event class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1393 
1394   Level: developer
1395 
1396 .seealso: [](ch_profiling), `PetscLogEventIncludeClass()`, `PetscLogEventExcludeClass()`, `PetscLogEventDeactivateClass()`, `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1397 @*/
1398 PetscErrorCode PetscLogEventActivateClass(PetscClassId classid)
1399 {
1400   PetscFunctionBegin;
1401   PetscCall(PetscLogClassSetActive(PETSC_DEFAULT, classid, PETSC_TRUE));
1402   PetscFunctionReturn(PETSC_SUCCESS);
1403 }
1404 
1405 /*@
1406   PetscLogEventDeactivateClass - Deactivates event logging for a PETSc object class for the current stage
1407 
1408   Not Collective
1409 
1410   Input Parameter:
1411 . classid - The event class, for example `MAT_CLASSID`, `SNES_CLASSID`, etc.
1412 
1413   Level: developer
1414 
1415 .seealso: [](ch_profiling), `PetscLogEventIncludeClass()`, `PetscLogEventExcludeClass()`, `PetscLogEventActivateClass()`, `PetscLogEventActivate()`, `PetscLogEventDeactivate()`
1416 @*/
1417 PetscErrorCode PetscLogEventDeactivateClass(PetscClassId classid)
1418 {
1419   PetscFunctionBegin;
1420   PetscCall(PetscLogClassSetActive(PETSC_DEFAULT, classid, PETSC_FALSE));
1421   PetscFunctionReturn(PETSC_SUCCESS);
1422 }
1423 
1424 /*MC
1425   PetscLogEventSync - Synchronizes the beginning of a user event.
1426 
1427   Synopsis:
1428   #include <petsclog.h>
1429   PetscErrorCode PetscLogEventSync(PetscLogEvent e, MPI_Comm comm)
1430 
1431   Collective
1432 
1433   Input Parameters:
1434 + e    - `PetscLogEvent` obtained from `PetscLogEventRegister()`
1435 - comm - an MPI communicator
1436 
1437   Example Usage:
1438 .vb
1439   PetscLogEvent USER_EVENT;
1440 
1441   PetscLogEventRegister("User event", 0, &USER_EVENT);
1442   PetscLogEventSync(USER_EVENT, PETSC_COMM_WORLD);
1443   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
1444   [code segment to monitor]
1445   PetscLogEventEnd(USER_EVENT, 0, 0, 0 , 0);
1446 .ve
1447 
1448   Level: developer
1449 
1450   Note:
1451   This routine should be called only if there is not a `PetscObject` available to pass to
1452   `PetscLogEventBegin()`.
1453 
1454 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`
1455 M*/
1456 
1457 /*MC
1458   PetscLogEventBegin - Logs the beginning of a user event.
1459 
1460   Synopsis:
1461   #include <petsclog.h>
1462   PetscErrorCode PetscLogEventBegin(PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)
1463 
1464   Not Collective
1465 
1466   Input Parameters:
1467 + e  - `PetscLogEvent` obtained from `PetscLogEventRegister()`
1468 . o1 - object associated with the event, or `NULL`
1469 . o2 - object associated with the event, or `NULL`
1470 . o3 - object associated with the event, or `NULL`
1471 - o4 - object associated with the event, or `NULL`
1472 
1473   Fortran Synopsis:
1474   void PetscLogEventBegin(int e, PetscErrorCode ierr)
1475 
1476   Example Usage:
1477 .vb
1478   PetscLogEvent USER_EVENT;
1479 
1480   PetscLogDouble user_event_flops;
1481   PetscLogEventRegister("User event",0, &USER_EVENT);
1482   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
1483   [code segment to monitor]
1484   PetscLogFlops(user_event_flops);
1485   PetscLogEventEnd(USER_EVENT, 0, 0, 0, 0);
1486 .ve
1487 
1488   Level: intermediate
1489 
1490   Developer Note:
1491   `PetscLogEventBegin()` and `PetscLogEventBegin()` return error codes instead of explicitly
1492   handling the errors that occur in the macro directly because other packages that use this
1493   macros have used them in their own functions or methods that do not return error codes and it
1494   would be disruptive to change the current behavior.
1495 
1496 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventEnd()`, `PetscLogFlops()`
1497 M*/
1498 
1499 /*MC
1500   PetscLogEventEnd - Log the end of a user event.
1501 
1502   Synopsis:
1503   #include <petsclog.h>
1504   PetscErrorCode PetscLogEventEnd(PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)
1505 
1506   Not Collective
1507 
1508   Input Parameters:
1509 + e  - `PetscLogEvent` obtained from `PetscLogEventRegister()`
1510 . o1 - object associated with the event, or `NULL`
1511 . o2 - object associated with the event, or `NULL`
1512 . o3 - object associated with the event, or `NULL`
1513 - o4 - object associated with the event, or `NULL`
1514 
1515   Fortran Synopsis:
1516   void PetscLogEventEnd(int e, PetscErrorCode ierr)
1517 
1518   Example Usage:
1519 .vb
1520   PetscLogEvent USER_EVENT;
1521 
1522   PetscLogDouble user_event_flops;
1523   PetscLogEventRegister("User event", 0, &USER_EVENT);
1524   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
1525   [code segment to monitor]
1526   PetscLogFlops(user_event_flops);
1527   PetscLogEventEnd(USER_EVENT, 0, 0, 0, 0);
1528 .ve
1529 
1530   Level: intermediate
1531 
1532 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogFlops()`
1533 M*/
1534 
1535 /*@C
1536   PetscLogStageGetPerfInfo - Return the performance information about the given stage
1537 
1538   No Fortran Support
1539 
1540   Input Parameters:
1541 . stage - The stage number or `PETSC_DETERMINE` for the current stage
1542 
1543   Output Parameter:
1544 . info - This structure is filled with the performance information
1545 
1546   Level: intermediate
1547 
1548   Notes:
1549   This is a low level routine used by the logging functions in PETSc.
1550 
1551   A `PETSCLOGHANDLERDEFAULT` must be running for this to work, having been started either with
1552   `PetscLogDefaultBegin()` or from the command line with `-log_view`.  If it was not started,
1553   all performance statistics in `info` will be zeroed.
1554 
1555 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogGetDefaultHandler()`
1556 @*/
1557 PetscErrorCode PetscLogStageGetPerfInfo(PetscLogStage stage, PetscEventPerfInfo *info)
1558 {
1559   PetscLogHandler     handler;
1560   PetscEventPerfInfo *event_info;
1561 
1562   PetscFunctionBegin;
1563   PetscAssertPointer(info, 2);
1564   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
1565   if (handler) {
1566     PetscCall(PetscLogHandlerGetStagePerfInfo(handler, stage, &event_info));
1567     *info = *event_info;
1568   } else {
1569     PetscCall(PetscInfo(NULL, "Default log handler is not running, PetscLogStageGetPerfInfo() returning zeros\n"));
1570     PetscCall(PetscMemzero(info, sizeof(*info)));
1571   }
1572   PetscFunctionReturn(PETSC_SUCCESS);
1573 }
1574 
1575 /*@C
1576   PetscLogEventGetPerfInfo - Return the performance information about the given event in the given stage
1577 
1578   No Fortran Support
1579 
1580   Input Parameters:
1581 + stage - The stage number or `PETSC_DETERMINE` for the current stage
1582 - event - The event number
1583 
1584   Output Parameter:
1585 . info - This structure is filled with the performance information
1586 
1587   Level: intermediate
1588 
1589   Note:
1590   This is a low level routine used by the logging functions in PETSc
1591 
1592   A `PETSCLOGHANDLERDEFAULT` must be running for this to work, having been started either with
1593   `PetscLogDefaultBegin()` or from the command line with `-log_view`.  If it was not started,
1594   all performance statistics in `info` will be zeroed.
1595 
1596 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogGetDefaultHandler()`
1597 @*/
1598 PetscErrorCode PetscLogEventGetPerfInfo(PetscLogStage stage, PetscLogEvent event, PetscEventPerfInfo *info)
1599 {
1600   PetscLogHandler     handler;
1601   PetscEventPerfInfo *event_info;
1602 
1603   PetscFunctionBegin;
1604   PetscAssertPointer(info, 3);
1605   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
1606   if (handler) {
1607     PetscCall(PetscLogHandlerGetEventPerfInfo(handler, stage, event, &event_info));
1608     *info = *event_info;
1609   } else {
1610     PetscCall(PetscInfo(NULL, "Default log handler is not running, PetscLogEventGetPerfInfo() returning zeros\n"));
1611     PetscCall(PetscMemzero(info, sizeof(*info)));
1612   }
1613   PetscFunctionReturn(PETSC_SUCCESS);
1614 }
1615 
1616 /*@
1617   PetscLogEventSetDof - Set the nth number of degrees of freedom of a numerical problem associated with this event
1618 
1619   Not Collective
1620 
1621   Input Parameters:
1622 + event - The event id to log
1623 . n     - The dof index, in [0, 8)
1624 - dof   - The number of dofs
1625 
1626   Options Database Key:
1627 . -log_view - Activates log summary
1628 
1629   Level: developer
1630 
1631   Note:
1632   This is to enable logging of convergence
1633 
1634 .seealso: `PetscLogEventSetError()`, `PetscLogEventRegister()`, `PetscLogGetDefaultHandler()`
1635 @*/
1636 PetscErrorCode PetscLogEventSetDof(PetscLogEvent event, PetscInt n, PetscLogDouble dof)
1637 {
1638   PetscFunctionBegin;
1639   PetscCheck(!(n < 0) && !(n > 7), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Error index %" PetscInt_FMT " is not in [0, 8)", n);
1640   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1641     PetscLogHandler h = PetscLogHandlers[i].handler;
1642 
1643     if (h) {
1644       PetscEventPerfInfo *event_info;
1645 
1646       PetscCall(PetscLogHandlerGetEventPerfInfo(h, PETSC_DEFAULT, event, &event_info));
1647       if (event_info) event_info->dof[n] = dof;
1648     }
1649   }
1650   PetscFunctionReturn(PETSC_SUCCESS);
1651 }
1652 
1653 /*@
1654   PetscLogEventSetError - Set the nth error associated with a numerical problem associated with this event
1655 
1656   Not Collective
1657 
1658   Input Parameters:
1659 + event - The event id to log
1660 . n     - The error index, in [0, 8)
1661 - error - The error
1662 
1663   Options Database Key:
1664 . -log_view - Activates log summary
1665 
1666   Level: developer
1667 
1668   Notes:
1669   This is to enable logging of convergence, and enable users to interpret the errors as they wish. For example,
1670   as different norms, or as errors for different fields
1671 
1672   This is a low level routine used by the logging functions in PETSc
1673 
1674 .seealso: `PetscLogEventSetDof()`, `PetscLogEventRegister()`, `PetscLogGetDefaultHandler()`
1675 @*/
1676 PetscErrorCode PetscLogEventSetError(PetscLogEvent event, PetscInt n, PetscLogDouble error)
1677 {
1678   PetscFunctionBegin;
1679   PetscCheck(!(n < 0) && !(n > 7), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Error index %" PetscInt_FMT " is not in [0, 8)", n);
1680   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1681     PetscLogHandler h = PetscLogHandlers[i].handler;
1682 
1683     if (h) {
1684       PetscEventPerfInfo *event_info;
1685 
1686       PetscCall(PetscLogHandlerGetEventPerfInfo(h, PETSC_DEFAULT, event, &event_info));
1687       if (event_info) event_info->errors[n] = error;
1688     }
1689   }
1690   PetscFunctionReturn(PETSC_SUCCESS);
1691 }
1692 
1693 /*@
1694   PetscLogEventGetId - Returns the event id when given the event name.
1695 
1696   Not Collective
1697 
1698   Input Parameter:
1699 . name - The event name
1700 
1701   Output Parameter:
1702 . event - The event, or -1 if no event with that name exists
1703 
1704   Level: intermediate
1705 
1706 .seealso: [](ch_profiling), `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogStageGetId()`
1707 @*/
1708 PetscErrorCode PetscLogEventGetId(const char name[], PetscLogEvent *event)
1709 {
1710   PetscLogState state;
1711 
1712   PetscFunctionBegin;
1713   *event = -1;
1714   PetscCall(PetscLogGetState(&state));
1715   if (state) PetscCall(PetscLogStateGetEventFromName(state, name, event));
1716   PetscFunctionReturn(PETSC_SUCCESS);
1717 }
1718 
1719 /*@
1720   PetscLogEventGetName - Returns the event name when given the event id.
1721 
1722   Not Collective
1723 
1724   Input Parameter:
1725 . event - The event
1726 
1727   Output Parameter:
1728 . name - The event name
1729 
1730   Level: intermediate
1731 
1732 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
1733 @*/
1734 PetscErrorCode PetscLogEventGetName(PetscLogEvent event, const char *name[])
1735 {
1736   PetscLogEventInfo event_info;
1737   PetscLogState     state;
1738 
1739   PetscFunctionBegin;
1740   *name = NULL;
1741   PetscCall(PetscLogGetState(&state));
1742   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
1743   PetscCall(PetscLogStateEventGetInfo(state, event, &event_info));
1744   *name = event_info.name;
1745   PetscFunctionReturn(PETSC_SUCCESS);
1746 }
1747 
1748 /*@
1749   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 
1751   Not collective
1752 
1753   Level: advanced
1754 
1755   Notes:
1756   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 
1758   Other log handlers (such as the nested handler, `PetscLogNestedBegin()`) will ignore this function.
1759 
1760 .seealso: [](ch_profiling), `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`, `PetscLogEventsResume()`, `PetscLogGetDefaultHandler()`
1761 @*/
1762 PetscErrorCode PetscLogEventsPause(void)
1763 {
1764   PetscFunctionBegin;
1765   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1766     PetscLogHandler h = PetscLogHandlers[i].handler;
1767 
1768     if (h) PetscCall(PetscLogHandlerEventsPause(h));
1769   }
1770   PetscFunctionReturn(PETSC_SUCCESS);
1771 }
1772 
1773 /*@
1774   PetscLogEventsResume - Return logging to normal behavior after it was paused with `PetscLogEventsPause()`.
1775 
1776   Not collective
1777 
1778   Level: advanced
1779 
1780 .seealso: [](ch_profiling), `PetscLogEventDeactivatePush()`, `PetscLogEventDeactivatePop()`, `PetscLogEventsPause()`, `PetscLogGetDefaultHandler()`
1781 @*/
1782 PetscErrorCode PetscLogEventsResume(void)
1783 {
1784   PetscFunctionBegin;
1785   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
1786     PetscLogHandler h = PetscLogHandlers[i].handler;
1787 
1788     if (h) PetscCall(PetscLogHandlerEventsResume(h));
1789   }
1790   PetscFunctionReturn(PETSC_SUCCESS);
1791 }
1792 
1793 /*------------------------------------------------ Class Functions --------------------------------------------------*/
1794 
1795 /*MC
1796    PetscLogObjectCreate - Log the creation of a `PetscObject`
1797 
1798    Synopsis:
1799    #include <petsclog.h>
1800    PetscErrorCode PetscLogObjectCreate(PetscObject h)
1801 
1802    Not Collective
1803 
1804    Input Parameters:
1805 .  h - A `PetscObject`
1806 
1807    Level: developer
1808 
1809    Developer Note:
1810      Called internally by PETSc when creating objects: users do not need to call this directly.
1811      Notification of the object creation is sent to each `PetscLogHandler` that is running.
1812 
1813 .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogObjectDestroy()`
1814 M*/
1815 
1816 /*MC
1817    PetscLogObjectDestroy - Logs the destruction of a `PetscObject`
1818 
1819    Synopsis:
1820    #include <petsclog.h>
1821    PetscErrorCode PetscLogObjectDestroy(PetscObject h)
1822 
1823    Not Collective
1824 
1825    Input Parameters:
1826 .  h - A `PetscObject`
1827 
1828    Level: developer
1829 
1830    Developer Note:
1831      Called internally by PETSc when destroying objects: users do not need to call this directly.
1832      Notification of the object creation is sent to each `PetscLogHandler` that is running.
1833 
1834 .seealso: [](ch_profiling), `PetscLogHandler`, `PetscLogObjectCreate()`
1835 M*/
1836 
1837 /*@
1838   PetscLogClassGetClassId - Returns the `PetscClassId` when given the class name.
1839 
1840   Not Collective
1841 
1842   Input Parameter:
1843 . name - The class name
1844 
1845   Output Parameter:
1846 . classid - The `PetscClassId` id, or -1 if no class with that name exists
1847 
1848   Level: intermediate
1849 
1850 .seealso: [](ch_profiling), `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscLogStageGetId()`
1851 @*/
1852 PetscErrorCode PetscLogClassGetClassId(const char name[], PetscClassId *classid)
1853 {
1854   PetscLogClass     log_class;
1855   PetscLogClassInfo class_info;
1856   PetscLogState     state;
1857 
1858   PetscFunctionBegin;
1859   *classid = -1;
1860   PetscCall(PetscLogGetState(&state));
1861   if (!state) PetscFunctionReturn(PETSC_SUCCESS);
1862   PetscCall(PetscLogStateGetClassFromName(state, name, &log_class));
1863   if (log_class < 0) {
1864     *classid = -1;
1865     PetscFunctionReturn(PETSC_SUCCESS);
1866   }
1867   PetscCall(PetscLogStateClassGetInfo(state, log_class, &class_info));
1868   *classid = class_info.classid;
1869   PetscFunctionReturn(PETSC_SUCCESS);
1870 }
1871 
1872 /*@C
1873   PetscLogClassIdGetName - Returns a `PetscClassId`'s name.
1874 
1875   Not Collective
1876 
1877   Input Parameter:
1878 . classid - A `PetscClassId`
1879 
1880   Output Parameter:
1881 . name - The class name
1882 
1883   Level: intermediate
1884 
1885 .seealso: [](ch_profiling), `PetscLogClassRegister()`, `PetscLogClassBegin()`, `PetscLogClassEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`, `PetscPreLoadClass()`
1886 @*/
1887 PetscErrorCode PetscLogClassIdGetName(PetscClassId classid, const char **name)
1888 {
1889   PetscLogClass     log_class;
1890   PetscLogClassInfo class_info;
1891   PetscLogState     state;
1892 
1893   PetscFunctionBegin;
1894   PetscCall(PetscLogGetState(&state));
1895   PetscCall(PetscLogStateGetClassFromClassId(state, classid, &log_class));
1896   PetscCall(PetscLogStateClassGetInfo(state, log_class, &class_info));
1897   *name = class_info.name;
1898   PetscFunctionReturn(PETSC_SUCCESS);
1899 }
1900 
1901 /*------------------------------------------------ Output Functions -------------------------------------------------*/
1902 /*@
1903   PetscLogDump - Dumps logs of objects to a file. This file is intended to
1904   be read by bin/petscview. This program no longer exists.
1905 
1906   Collective on `PETSC_COMM_WORLD`
1907 
1908   Input Parameter:
1909 . sname - an optional file name
1910 
1911   Example Usage:
1912 .vb
1913   PetscInitialize(...);
1914   PetscLogDefaultBegin();
1915   // ... code ...
1916   PetscLogDump(filename);
1917   PetscFinalize();
1918 .ve
1919 
1920   Level: advanced
1921 
1922   Note:
1923   The default file name is Log.<rank> where <rank> is the MPI process rank. If no name is specified,
1924   this file will be used.
1925 
1926 .seealso: [](ch_profiling), `PetscLogDefaultBegin()`, `PetscLogView()`, `PetscLogGetDefaultHandler()`
1927 @*/
1928 PetscErrorCode PetscLogDump(const char sname[])
1929 {
1930   PetscLogHandler handler;
1931 
1932   PetscFunctionBegin;
1933   PetscCall(PetscLogGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
1934   PetscCall(PetscLogHandlerDump(handler, sname));
1935   PetscFunctionReturn(PETSC_SUCCESS);
1936 }
1937 
1938 /*@
1939   PetscLogMPEDump - Dumps the MPE logging info to file for later use with Jumpshot.
1940 
1941   Collective on `PETSC_COMM_WORLD`
1942 
1943   Input Parameter:
1944 . sname - filename for the MPE logfile
1945 
1946   Level: advanced
1947 
1948 .seealso: [](ch_profiling), `PetscLogDump()`, `PetscLogMPEBegin()`
1949 @*/
1950 PetscErrorCode PetscLogMPEDump(const char sname[])
1951 {
1952   PetscFunctionBegin;
1953   #if defined(PETSC_HAVE_MPE)
1954   if (PetscBeganMPE) {
1955     char name[PETSC_MAX_PATH_LEN];
1956 
1957     PetscCall(PetscInfo(0, "Finalizing MPE.\n"));
1958     if (sname) {
1959       PetscCall(PetscStrncpy(name, sname, sizeof(name)));
1960     } else {
1961       PetscCall(PetscGetProgramName(name, sizeof(name)));
1962     }
1963     PetscCall(MPE_Finish_log(name));
1964   } else {
1965     PetscCall(PetscInfo(0, "Not finalizing MPE (not started by PETSc).\n"));
1966   }
1967   #else
1968   SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "PETSc was configured without MPE support, reconfigure with --with-mpe or --download-mpe");
1969   #endif
1970   PetscFunctionReturn(PETSC_SUCCESS);
1971 }
1972 
1973 /*@
1974   PetscLogView - Prints a summary of the logging.
1975 
1976   Collective
1977 
1978   Input Parameter:
1979 . viewer - an ASCII viewer
1980 
1981   Options Database Keys:
1982 + -log_view [:filename]                    - Prints summary of log information
1983 . -log_view :filename.py:ascii_info_detail - Saves logging information from each process as a Python file
1984 . -log_view :filename.xml:ascii_xml        - Saves a summary of the logging information in a nested format (see below for how to view it)
1985 . -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)
1986 . -log_view_memory                         - Also display memory usage in each event
1987 . -log_view_gpu_time                       - Also display time in each event for GPU kernels (Note this may slow the computation)
1988 . -log_view_gpu_energy                     - Also display energy (estimated with power*gtime) in Joules for GPU kernels
1989 . -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.
1990 . -log_all                                 - Saves a file Log.rank for each MPI rank with details of each step of the computation
1991 - -log_trace [filename]                    - Displays a trace of what each process is doing
1992 
1993   Level: beginner
1994 
1995   Notes:
1996   It is possible to control the logging programmatically but we recommend using the options database approach whenever possible
1997   By default the summary is printed to stdout.
1998 
1999   Before calling this routine you must have called either PetscLogDefaultBegin() or PetscLogNestedBegin()
2000 
2001   If PETSc is configured with --with-logging=0 then this functionality is not available
2002 
2003   To view the nested XML format filename.xml first copy  ${PETSC_DIR}/share/petsc/xml/performance_xml2html.xsl to the current
2004   directory then open filename.xml with your browser. Specific notes for certain browsers
2005 .vb
2006     Firefox and Internet explorer - simply open the file
2007     Google Chrome - you must start up Chrome with the option --allow-file-access-from-files
2008     Safari - see https://ccm.net/faq/36342-safari-how-to-enable-local-file-access
2009 .ve
2010   or one can use the package <http://xmlsoft.org/XSLT/xsltproc2.html> to translate the xml file to html and then open it with
2011   your browser.
2012   Alternatively, use the script ${PETSC_DIR}/lib/petsc/bin/petsc-performance-view to automatically open a new browser
2013   window and render the XML log file contents.
2014 
2015   The nested XML format was kindly donated by Koos Huijssen and Christiaan M. Klaij  MARITIME  RESEARCH  INSTITUTE  NETHERLANDS
2016 
2017   The Flame Graph output can be visualised using either the original Flame Graph script <https://github.com/brendangregg/FlameGraph>
2018   or using speedscope <https://www.speedscope.app>.
2019   Old XML profiles may be converted into this format using the script ${PETSC_DIR}/lib/petsc/bin/xml2flamegraph.py.
2020 
2021 .seealso: [](ch_profiling), `PetscLogDefaultBegin()`, `PetscLogDump()`
2022 @*/
2023 PetscErrorCode PetscLogView(PetscViewer viewer)
2024 {
2025   PetscBool         isascii;
2026   PetscViewerFormat format;
2027   int               stage;
2028   PetscLogState     state;
2029   PetscIntStack     temp_stack;
2030   PetscLogHandler   handler;
2031   PetscBool         is_empty;
2032 
2033   PetscFunctionBegin;
2034   PetscCall(PetscLogGetState(&state));
2035   /* Pop off any stages the user forgot to remove */
2036   PetscCall(PetscIntStackCreate(&temp_stack));
2037   PetscCall(PetscLogStateGetCurrentStage(state, &stage));
2038   while (stage >= 0) {
2039     PetscCall(PetscLogStagePop());
2040     PetscCall(PetscIntStackPush(temp_stack, stage));
2041     PetscCall(PetscLogStateGetCurrentStage(state, &stage));
2042   }
2043   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
2044   PetscCheck(isascii, PetscObjectComm((PetscObject)viewer), PETSC_ERR_SUP, "Currently can only view logging to ASCII");
2045   PetscCall(PetscViewerGetFormat(viewer, &format));
2046   if (format == PETSC_VIEWER_ASCII_XML || format == PETSC_VIEWER_ASCII_FLAMEGRAPH) {
2047     PetscCall(PetscLogGetHandler(PETSCLOGHANDLERNESTED, &handler));
2048     PetscCall(PetscLogHandlerView(handler, viewer));
2049   } else {
2050     PetscCall(PetscLogGetHandler(PETSCLOGHANDLERDEFAULT, &handler));
2051     PetscCall(PetscLogHandlerView(handler, viewer));
2052   }
2053   PetscCall(PetscIntStackEmpty(temp_stack, &is_empty));
2054   while (!is_empty) {
2055     PetscCall(PetscIntStackPop(temp_stack, &stage));
2056     PetscCall(PetscLogStagePush(stage));
2057     PetscCall(PetscIntStackEmpty(temp_stack, &is_empty));
2058   }
2059   PetscCall(PetscIntStackDestroy(temp_stack));
2060   PetscFunctionReturn(PETSC_SUCCESS);
2061 }
2062 
2063 /*@C
2064   PetscLogViewFromOptions - Processes command line options to determine if/how a `PetscLog` is to be viewed.
2065 
2066   Collective on `PETSC_COMM_WORLD`
2067 
2068   Level: developer
2069 
2070 .seealso: [](ch_profiling), `PetscLogView()`
2071 @*/
2072 PetscErrorCode PetscLogViewFromOptions(void)
2073 {
2074   PetscInt          n_max = PETSC_LOG_VIEW_FROM_OPTIONS_MAX;
2075   PetscViewer       viewers[PETSC_LOG_VIEW_FROM_OPTIONS_MAX];
2076   PetscViewerFormat formats[PETSC_LOG_VIEW_FROM_OPTIONS_MAX];
2077   PetscBool         flg;
2078 
2079   PetscFunctionBegin;
2080   PetscCall(PetscOptionsCreateViewers(PETSC_COMM_WORLD, NULL, NULL, "-log_view", &n_max, viewers, formats, &flg));
2081   /*
2082      PetscLogHandlerView_Default_Info() wants to be sure that the only objects still around are these viewers, so keep track of how many there are
2083    */
2084   PetscLogNumViewersCreated = n_max;
2085   for (PetscInt i = 0; i < n_max; i++) {
2086     PetscInt refct;
2087 
2088     PetscCall(PetscViewerPushFormat(viewers[i], formats[i]));
2089     PetscCall(PetscLogView(viewers[i]));
2090     PetscCall(PetscViewerPopFormat(viewers[i]));
2091     PetscCall(PetscObjectGetReference((PetscObject)viewers[i], &refct));
2092     PetscCall(PetscViewerDestroy(&viewers[i]));
2093     if (refct == 1) PetscLogNumViewersDestroyed++;
2094   }
2095   PetscLogNumViewersDestroyed = 0;
2096   PetscFunctionReturn(PETSC_SUCCESS);
2097 }
2098 
2099 PETSC_INTERN PetscErrorCode PetscLogHandlerNestedSetThreshold(PetscLogHandler, PetscLogDouble, PetscLogDouble *);
2100 
2101 /*@
2102   PetscLogSetThreshold - Set the threshold time for logging the events; this is a percentage out of 100, so 1. means any event
2103   that takes 1 or more percent of the time.
2104 
2105   Logically Collective on `PETSC_COMM_WORLD`
2106 
2107   Input Parameter:
2108 . newThresh - the threshold to use
2109 
2110   Output Parameter:
2111 . oldThresh - the previously set threshold value
2112 
2113   Options Database Keys:
2114 . -log_view :filename.xml:ascii_xml - Prints an XML summary of flop and timing information to the file
2115 
2116   Example Usage:
2117 .vb
2118   PetscInitialize(...);
2119   PetscLogNestedBegin();
2120   PetscLogSetThreshold(0.1,&oldthresh);
2121   // ... code ...
2122   PetscLogView(viewer);
2123   PetscFinalize();
2124 .ve
2125 
2126   Level: advanced
2127 
2128   Note:
2129   This threshold is only used by the nested log handler
2130 
2131 .seealso: `PetscLogDump()`, `PetscLogView()`, `PetscLogTraceBegin()`, `PetscLogDefaultBegin()`,
2132           `PetscLogNestedBegin()`
2133 @*/
2134 PetscErrorCode PetscLogSetThreshold(PetscLogDouble newThresh, PetscLogDouble *oldThresh)
2135 {
2136   PetscLogHandler handler;
2137 
2138   PetscFunctionBegin;
2139   PetscCall(PetscLogTryGetHandler(PETSCLOGHANDLERNESTED, &handler));
2140   PetscCall(PetscLogHandlerNestedSetThreshold(handler, newThresh, oldThresh));
2141   PetscFunctionReturn(PETSC_SUCCESS);
2142 }
2143 
2144 /*----------------------------------------------- Counter Functions -------------------------------------------------*/
2145 /*@
2146   PetscGetFlops - Returns the number of flops used on this processor
2147   since the program began.
2148 
2149   Not Collective
2150 
2151   Output Parameter:
2152 . flops - number of floating point operations
2153 
2154   Level: intermediate
2155 
2156   Notes:
2157   A global counter logs all PETSc flop counts.  The user can use
2158   `PetscLogFlops()` to increment this counter to include flops for the
2159   application code.
2160 
2161   A separate counter `PetscLogGpuFlops()` logs the flops that occur on any GPU associated with this MPI rank
2162 
2163 .seealso: [](ch_profiling), `PetscLogGpuFlops()`, `PetscTime()`, `PetscLogFlops()`
2164 @*/
2165 PetscErrorCode PetscGetFlops(PetscLogDouble *flops)
2166 {
2167   PetscFunctionBegin;
2168   *flops = petsc_TotalFlops;
2169   PetscFunctionReturn(PETSC_SUCCESS);
2170 }
2171 
2172 /*@C
2173   PetscLogObjectState - Record information about an object with the default log handler
2174 
2175   Not Collective
2176 
2177   Input Parameters:
2178 + obj    - the `PetscObject`
2179 . format - a printf-style format string
2180 - ...    - printf arguments to format
2181 
2182   Level: developer
2183 
2184 .seealso: [](ch_profiling), `PetscLogObjectCreate()`, `PetscLogObjectDestroy()`, `PetscLogGetDefaultHandler()`
2185 @*/
2186 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...)
2187 {
2188   PetscFunctionBegin;
2189   for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
2190     PetscLogHandler h = PetscLogHandlers[i].handler;
2191 
2192     if (h) {
2193       va_list Argp;
2194       va_start(Argp, format);
2195       PetscCall(PetscLogHandlerLogObjectState_Internal(h, obj, format, Argp));
2196       va_end(Argp);
2197     }
2198   }
2199   PetscFunctionReturn(PETSC_SUCCESS);
2200 }
2201 
2202 /*MC
2203   PetscLogFlops - Adds floating point operations to the global counter.
2204 
2205   Synopsis:
2206   #include <petsclog.h>
2207   PetscErrorCode PetscLogFlops(PetscLogDouble f)
2208 
2209   Not Collective
2210 
2211   Input Parameter:
2212 . f - flop counter
2213 
2214   Example Usage:
2215 .vb
2216   PetscLogEvent USER_EVENT;
2217 
2218   PetscLogEventRegister("User event", 0, &USER_EVENT);
2219   PetscLogEventBegin(USER_EVENT, 0, 0, 0, 0);
2220   [code segment to monitor]
2221   PetscLogFlops(user_flops)
2222   PetscLogEventEnd(USER_EVENT, 0, 0, 0, 0);
2223 .ve
2224 
2225   Level: intermediate
2226 
2227   Note:
2228    A global counter logs all PETSc flop counts. The user can use PetscLogFlops() to increment
2229    this counter to include flops for the application code.
2230 
2231 .seealso: [](ch_profiling), `PetscLogGpuFlops()`, `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscGetFlops()`
2232 M*/
2233 
2234 /*MC
2235   PetscPreLoadBegin - Begin a segment of code that may be preloaded (run twice) to get accurate
2236   timings
2237 
2238   Synopsis:
2239   #include <petsclog.h>
2240   void PetscPreLoadBegin(PetscBool flag, char *name);
2241 
2242   Not Collective
2243 
2244   Input Parameters:
2245 + flag - `PETSC_TRUE` to run twice, `PETSC_FALSE` to run once, may be overridden with command
2246          line option `-preload true|false`
2247 - name - name of first stage (lines of code timed separately with `-log_view`) to be preloaded
2248 
2249   Example Usage:
2250 .vb
2251   PetscPreLoadBegin(PETSC_TRUE, "first stage");
2252   // lines of code
2253   PetscPreLoadStage("second stage");
2254   // lines of code
2255   PetscPreLoadEnd();
2256 .ve
2257 
2258   Level: intermediate
2259 
2260   Note:
2261   Only works in C/C++, not Fortran
2262 
2263   Flags available within the macro\:
2264 + PetscPreLoadingUsed - `PETSC_TRUE` if we are or have done preloading
2265 . PetscPreLoadingOn   - `PETSC_TRUE` if it is CURRENTLY doing preload
2266 . PetscPreLoadIt      - `0` for the first computation (with preloading turned off it is only
2267                         `0`) `1`  for the second
2268 - PetscPreLoadMax     - number of times it will do the computation, only one when preloading is
2269                         turned on
2270 
2271   The first two variables are available throughout the program, the second two only between the
2272   `PetscPreLoadBegin()` and `PetscPreLoadEnd()`
2273 
2274 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadEnd()`, `PetscPreLoadStage()`
2275 M*/
2276 
2277 /*MC
2278   PetscPreLoadEnd - End a segment of code that may be preloaded (run twice) to get accurate
2279   timings
2280 
2281   Synopsis:
2282   #include <petsclog.h>
2283   void PetscPreLoadEnd(void);
2284 
2285   Not Collective
2286 
2287   Example Usage:
2288 .vb
2289   PetscPreLoadBegin(PETSC_TRUE, "first stage");
2290   // lines of code
2291   PetscPreLoadStage("second stage");
2292   // lines of code
2293   PetscPreLoadEnd();
2294 .ve
2295 
2296   Level: intermediate
2297 
2298   Note:
2299   Only works in C/C++ not Fortran
2300 
2301 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadStage()`
2302 M*/
2303 
2304 /*MC
2305   PetscPreLoadStage - Start a new segment of code to be timed separately to get accurate timings
2306 
2307   Synopsis:
2308   #include <petsclog.h>
2309   void PetscPreLoadStage(char *name);
2310 
2311   Not Collective
2312 
2313   Example Usage:
2314 .vb
2315   PetscPreLoadBegin(PETSC_TRUE,"first stage");
2316   // lines of code
2317   PetscPreLoadStage("second stage");
2318   // lines of code
2319   PetscPreLoadEnd();
2320 .ve
2321 
2322   Level: intermediate
2323 
2324   Note:
2325   Only works in C/C++ not Fortran
2326 
2327 .seealso: [](ch_profiling), `PetscLogEventRegister()`, `PetscLogEventBegin()`, `PetscLogEventEnd()`, `PetscPreLoadBegin()`, `PetscPreLoadEnd()`
2328 M*/
2329 
2330   #if PetscDefined(HAVE_DEVICE)
2331     #include <petsc/private/deviceimpl.h>
2332 
2333 /*@
2334   PetscLogGpuTime - turn on the logging of GPU time for GPU kernels
2335 
2336   Options Database Key:
2337 . -log_view_gpu_time - provide the GPU times for all events in the `-log_view` output
2338 
2339   Level: advanced
2340 
2341   Notes:
2342   Turning on the timing of the GPU kernels can slow down the entire computation and should only
2343   be used when studying the performance of individual operations on GPU such as vector operations and
2344   matrix-vector operations.
2345 
2346   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 
2348   This routine should only be called once near the beginning of the program. Once it is started
2349   it cannot be turned off.
2350 
2351 .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuFlops()`, `PetscLogGpuTimeEnd()`, `PetscLogGpuTimeBegin()`
2352 @*/
2353 PetscErrorCode PetscLogGpuTime(void)
2354 {
2355   PetscFunctionBegin;
2356   PetscCheck(petsc_gtime == 0.0, PETSC_COMM_SELF, PETSC_ERR_SUP, "GPU logging has already been turned on");
2357   PetscLogGpuTimeFlag = PETSC_TRUE;
2358   PetscFunctionReturn(PETSC_SUCCESS);
2359 }
2360 
2361 /*@
2362   PetscLogGpuTimeBegin - Start timer for device
2363 
2364   Level: intermediate
2365 
2366   Notes:
2367   When GPU is enabled, the timer is run on the GPU, it is a separate logging of time
2368   devoted to GPU computations (excluding kernel launch times).
2369 
2370   When GPU is not available, the timer is run on the CPU, it is a separate logging of
2371   time devoted to GPU computations (including kernel launch times).
2372 
2373   There is no need to call WaitForCUDA() or WaitForHIP() between `PetscLogGpuTimeBegin()` and
2374   `PetscLogGpuTimeEnd()`
2375 
2376   This timer should NOT include times for data transfers between the GPU and CPU, nor setup
2377   actions such as allocating space.
2378 
2379   The regular logging captures the time for data transfers and any CPU activities during the
2380   event. It is used to compute the flop rate on the GPU as it is actively engaged in running a
2381   kernel.
2382 
2383   Developer Notes:
2384   The GPU event timer captures the execution time of all the kernels launched in the default
2385   stream by the CPU between `PetscLogGpuTimeBegin()` and `PetscLogGpuTimeEnd()`.
2386 
2387   `PetscLogGpuTimeBegin()` and `PetscLogGpuTimeEnd()` insert the begin and end events into the
2388   default stream (stream 0). The device will record a time stamp for the event when it reaches
2389   that event in the stream. The function xxxEventSynchronize() is called in
2390   `PetscLogGpuTimeEnd()` to block CPU execution, but not continued GPU execution, until the
2391   timer event is recorded.
2392 
2393 .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuFlops()`, `PetscLogGpuTimeEnd()`, `PetscLogGpuTime()`
2394 @*/
2395 PetscErrorCode PetscLogGpuTimeBegin(void)
2396 {
2397   PetscBool isActive;
2398 
2399   PetscFunctionBegin;
2400   PetscCall(PetscLogEventBeginIsActive(&isActive));
2401   if (!isActive || !PetscLogGpuTimeFlag) PetscFunctionReturn(PETSC_SUCCESS);
2402   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2403     PetscDeviceContext dctx;
2404 
2405     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2406     PetscCall(PetscDeviceContextBeginTimer_Internal(dctx));
2407   } else {
2408     PetscCall(PetscTimeSubtract(&petsc_gtime));
2409   }
2410   PetscFunctionReturn(PETSC_SUCCESS);
2411 }
2412 
2413 /*@
2414   PetscLogGpuTimeEnd - Stop timer for device
2415 
2416   Level: intermediate
2417 
2418 .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuFlops()`, `PetscLogGpuTimeBegin()`
2419 @*/
2420 PetscErrorCode PetscLogGpuTimeEnd(void)
2421 {
2422   PetscBool isActive;
2423 
2424   PetscFunctionBegin;
2425   PetscCall(PetscLogEventEndIsActive(&isActive));
2426   if (!isActive || !PetscLogGpuTimeFlag) PetscFunctionReturn(PETSC_SUCCESS);
2427   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2428     PetscDeviceContext dctx;
2429     PetscLogDouble     elapsed;
2430 
2431     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2432     PetscCall(PetscDeviceContextEndTimer_Internal(dctx, &elapsed));
2433     petsc_gtime += (elapsed / 1000.0);
2434     #if PetscDefined(HAVE_CUDA_VERSION_12_2PLUS)
2435     if (PetscLogGpuEnergyFlag) {
2436       PetscLogDouble power;
2437       PetscCall(PetscDeviceContextGetPower_Internal(dctx, &power));
2438       petsc_genergy += (power * elapsed / 1000000.0); // convert to Joules
2439     }
2440     #endif
2441   } else {
2442     PetscCall(PetscTimeAdd(&petsc_gtime));
2443   }
2444   PetscFunctionReturn(PETSC_SUCCESS);
2445 }
2446 
2447 /*@
2448   PetscLogGpuEnergy - turn on the logging of GPU energy (estimated with power*gtime) for GPU kernels
2449 
2450   Options Database Key:
2451 . -log_view_gpu_energy - provide the GPU energy consumption (estimated with power*gtime) for all events in the `-log_view` output
2452 
2453   Level: advanced
2454 
2455   Note:
2456   This option is mutually exclusive to `-log_view_gpu_energy_meter`.
2457 
2458   Developer Note:
2459   This option turns on energy monitoring of GPU kernels and requires CUDA version >= 12.2. The energy consumption is estimated as
2460   instant_power * gpu_kernel_time. Due to the delay in NVML power sampling, we read the instantaneous power draw at the end of each
2461   event using `nvmlDeviceGetFieldValues()` with the field ID `NVML_FI_DEV_POWER_INSTANT`.
2462 
2463 .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuEnergyMeter()`
2464 @*/
2465 PetscErrorCode PetscLogGpuEnergy(void)
2466 {
2467   PetscFunctionBegin;
2468   PetscCheck(PetscDefined(HAVE_CUDA_VERSION_12_2PLUS), PETSC_COMM_WORLD, PETSC_ERR_SUP_SYS, "-log_view_gpu_energy requires CUDA version >= 12.2");
2469   PetscCheck(petsc_genergy == 0.0, PETSC_COMM_SELF, PETSC_ERR_SUP, "GPU energy logging has already been turned on");
2470   PetscLogGpuEnergyFlag      = PETSC_TRUE;
2471   PetscLogGpuEnergyMeterFlag = PETSC_FALSE;
2472   PetscFunctionReturn(PETSC_SUCCESS);
2473 }
2474 
2475 /*@
2476   PetscLogGpuEnergyMeter - turn on the logging of GPU energy (readings from energy meters) for GPU kernels
2477 
2478   Options Database Key:
2479 . -log_view_gpu_energy_meter - provide the GPU energy (readings from energy meters) consumption for all events in the `-log_view` output
2480 
2481   Level: advanced
2482 
2483   Note:
2484   This option is mutually exclusive to `-log_view_gpu_energy`.
2485 
2486   Developer Note:
2487   This option turns on energy monitoring of GPU kernels. The energy consumption is measured directly using the NVML API
2488   `nvmlDeviceGetTotalEnergyConsumption()`, which returns the total energy used by the GPU since the driver was last initialized.
2489   For newer GPUs, energy readings are updated every 20-100ms, so this approach may be inaccurate for short-duration GPU events.
2490 
2491 @*/
2492 PetscErrorCode PetscLogGpuEnergyMeter(void)
2493 {
2494   PetscFunctionBegin;
2495   PetscCheck(petsc_genergy == 0.0, PETSC_COMM_SELF, PETSC_ERR_SUP, "GPU energy logging has already been turned on");
2496   PetscLogGpuEnergyMeterFlag = PETSC_TRUE;
2497   PetscLogGpuEnergyFlag      = PETSC_FALSE;
2498   PetscFunctionReturn(PETSC_SUCCESS);
2499 }
2500 
2501 /*@
2502   PetscLogGpuEnergyMeterBegin - Start energy meter for device
2503 
2504   Level: intermediate
2505 
2506   Notes:
2507   The GPU event energy meter captures the energy used by the GPU between `PetscLogGpuEnergyMeterBegin()` and `PetscLogGpuEnergyMeterEnd()`.
2508 
2509   `PetscLogGpuEnergyMeterBegin()` and `PetscLogGpuEnergyMeterEnd()` collect the energy readings using `nvmlDeviceGetTotalEnergyConsumption()`.
2510   The function `cupmStreamSynchronize()` is called before the energy query to ensure completion.
2511 
2512 .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuEnergyMeterEnd()`, `PetscLogGpuEnergyMeter()`
2513 @*/
2514 PetscErrorCode PetscLogGpuEnergyMeterBegin(void)
2515 {
2516   PetscBool isActive;
2517 
2518   PetscFunctionBegin;
2519   PetscCall(PetscLogEventBeginIsActive(&isActive));
2520   if (!isActive || !PetscLogGpuEnergyMeterFlag) PetscFunctionReturn(PETSC_SUCCESS);
2521   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2522     PetscDeviceContext dctx;
2523 
2524     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2525     PetscCall(PetscDeviceContextBeginEnergyMeter_Internal(dctx));
2526   }
2527   PetscFunctionReturn(PETSC_SUCCESS);
2528 }
2529 
2530 /*@
2531   PetscLogGpuEnergyMeterEnd - Stop energy meter for device
2532 
2533   Level: intermediate
2534 
2535 .seealso: [](ch_profiling), `PetscLogView()`, `PetscLogGpuEnergyMeterBegin()`
2536 @*/
2537 PetscErrorCode PetscLogGpuEnergyMeterEnd(void)
2538 {
2539   PetscBool isActive;
2540 
2541   PetscFunctionBegin;
2542   PetscCall(PetscLogEventEndIsActive(&isActive));
2543   if (!isActive || !PetscLogGpuEnergyMeterFlag) PetscFunctionReturn(PETSC_SUCCESS);
2544   if (!PetscDefined(HAVE_KOKKOS_WITHOUT_GPU)) {
2545     PetscDeviceContext dctx;
2546     PetscLogDouble     energy;
2547 
2548     PetscCall(PetscDeviceContextGetCurrentContext(&dctx));
2549     PetscCall(PetscDeviceContextEndEnergyMeter_Internal(dctx, &energy));
2550     petsc_genergy_meter += (energy / 1000.0); // convert to Joules
2551   }
2552   PetscFunctionReturn(PETSC_SUCCESS);
2553 }
2554   #endif /* end of PETSC_HAVE_DEVICE */
2555 
2556 #endif /* PETSC_USE_LOG*/
2557 
2558 /* -- Utility functions for logging from Fortran -- */
2559 
2560 PETSC_EXTERN PetscErrorCode PetscASend(int count, int datatype)
2561 {
2562   PetscFunctionBegin;
2563 #if PetscDefined(USE_LOG)
2564   PetscCall(PetscAddLogDouble(&petsc_send_ct, &petsc_send_ct_th, 1));
2565   #if !defined(MPIUNI_H) && !defined(PETSC_HAVE_BROKEN_RECURSIVE_MACRO)
2566   PetscCall(PetscMPITypeSize(count, MPI_Type_f2c((MPI_Fint)datatype), &petsc_send_len, &petsc_send_len_th));
2567   #endif
2568 #endif
2569   PetscFunctionReturn(PETSC_SUCCESS);
2570 }
2571 
2572 PETSC_EXTERN PetscErrorCode PetscARecv(int count, int datatype)
2573 {
2574   PetscFunctionBegin;
2575 #if PetscDefined(USE_LOG)
2576   PetscCall(PetscAddLogDouble(&petsc_recv_ct, &petsc_recv_ct_th, 1));
2577   #if !defined(MPIUNI_H) && !defined(PETSC_HAVE_BROKEN_RECURSIVE_MACRO)
2578   PetscCall(PetscMPITypeSize(count, MPI_Type_f2c((MPI_Fint)datatype), &petsc_recv_len, &petsc_recv_len_th));
2579   #endif
2580 #endif
2581   PetscFunctionReturn(PETSC_SUCCESS);
2582 }
2583 
2584 PETSC_EXTERN PetscErrorCode PetscAReduce(void)
2585 {
2586   PetscFunctionBegin;
2587   if (PetscDefined(USE_LOG)) PetscCall(PetscAddLogDouble(&petsc_allreduce_ct, &petsc_allreduce_ct_th, 1));
2588   PetscFunctionReturn(PETSC_SUCCESS);
2589 }
2590 
2591 PetscClassId PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID;
2592 PetscClassId PETSC_OBJECT_CLASSID  = 0;
2593 
2594 static PetscBool PetscLogInitializeCalled = PETSC_FALSE;
2595 
2596 PETSC_INTERN PetscErrorCode PetscLogInitialize(void)
2597 {
2598   int stage;
2599 
2600   PetscFunctionBegin;
2601   if (PetscLogInitializeCalled) PetscFunctionReturn(PETSC_SUCCESS);
2602   PetscLogInitializeCalled = PETSC_TRUE;
2603   if (PetscDefined(USE_LOG)) {
2604     /* Setup default logging structures */
2605     PetscCall(PetscLogStateCreate(&petsc_log_state));
2606     for (PetscInt i = 0; i < PETSC_LOG_HANDLER_MAX; i++) {
2607       if (PetscLogHandlers[i].handler) PetscCall(PetscLogHandlerSetState(PetscLogHandlers[i].handler, petsc_log_state));
2608     }
2609     PetscCall(PetscLogStateStageRegister(petsc_log_state, "Main Stage", &stage));
2610     PetscCall(PetscSpinlockCreate(&PetscLogSpinLock));
2611 #if defined(PETSC_HAVE_THREADSAFETY)
2612     petsc_log_tid = 0;
2613     petsc_log_gid = 0;
2614 #endif
2615 
2616     /* All processors sync here for more consistent logging */
2617     PetscCallMPI(MPI_Barrier(PETSC_COMM_WORLD));
2618     PetscCall(PetscTime(&petsc_BaseTime));
2619     PetscCall(PetscLogStagePush(stage));
2620   }
2621   PetscFunctionReturn(PETSC_SUCCESS);
2622 }
2623 
2624 PETSC_INTERN PetscErrorCode PetscLogFinalize(void)
2625 {
2626   PetscFunctionBegin;
2627   if (PetscDefined(USE_LOG)) {
2628     /* Resetting phase */
2629     // pop remaining stages
2630     if (petsc_log_state) {
2631       while (petsc_log_state->current_stage >= 0) PetscCall(PetscLogStagePop());
2632     }
2633     for (int i = 0; i < PETSC_LOG_HANDLER_MAX; i++) PetscCall(PetscLogHandlerDestroy(&PetscLogHandlers[i].handler));
2634     PetscCall(PetscArrayzero(PetscLogHandlers, PETSC_LOG_HANDLER_MAX));
2635     PetscCall(PetscLogStateDestroy(&petsc_log_state));
2636 
2637     petsc_TotalFlops         = 0.0;
2638     petsc_BaseTime           = 0.0;
2639     petsc_TotalFlops         = 0.0;
2640     petsc_send_ct            = 0.0;
2641     petsc_recv_ct            = 0.0;
2642     petsc_send_len           = 0.0;
2643     petsc_recv_len           = 0.0;
2644     petsc_isend_ct           = 0.0;
2645     petsc_irecv_ct           = 0.0;
2646     petsc_isend_len          = 0.0;
2647     petsc_irecv_len          = 0.0;
2648     petsc_wait_ct            = 0.0;
2649     petsc_wait_any_ct        = 0.0;
2650     petsc_wait_all_ct        = 0.0;
2651     petsc_sum_of_waits_ct    = 0.0;
2652     petsc_allreduce_ct       = 0.0;
2653     petsc_gather_ct          = 0.0;
2654     petsc_scatter_ct         = 0.0;
2655     petsc_TotalFlops_th      = 0.0;
2656     petsc_send_ct_th         = 0.0;
2657     petsc_recv_ct_th         = 0.0;
2658     petsc_send_len_th        = 0.0;
2659     petsc_recv_len_th        = 0.0;
2660     petsc_isend_ct_th        = 0.0;
2661     petsc_irecv_ct_th        = 0.0;
2662     petsc_isend_len_th       = 0.0;
2663     petsc_irecv_len_th       = 0.0;
2664     petsc_wait_ct_th         = 0.0;
2665     petsc_wait_any_ct_th     = 0.0;
2666     petsc_wait_all_ct_th     = 0.0;
2667     petsc_sum_of_waits_ct_th = 0.0;
2668     petsc_allreduce_ct_th    = 0.0;
2669     petsc_gather_ct_th       = 0.0;
2670     petsc_scatter_ct_th      = 0.0;
2671 
2672     petsc_ctog_ct       = 0.0;
2673     petsc_gtoc_ct       = 0.0;
2674     petsc_ctog_sz       = 0.0;
2675     petsc_gtoc_sz       = 0.0;
2676     petsc_gflops        = 0.0;
2677     petsc_gtime         = 0.0;
2678     petsc_genergy       = 0.0;
2679     petsc_genergy_meter = 0.0;
2680     petsc_ctog_ct_th    = 0.0;
2681     petsc_gtoc_ct_th    = 0.0;
2682     petsc_ctog_sz_th    = 0.0;
2683     petsc_gtoc_sz_th    = 0.0;
2684     petsc_gflops_th     = 0.0;
2685     petsc_gtime_th      = 0.0;
2686   }
2687   PETSC_LARGEST_CLASSID    = PETSC_SMALLEST_CLASSID;
2688   PETSC_OBJECT_CLASSID     = 0;
2689   PetscLogInitializeCalled = PETSC_FALSE;
2690   PetscFunctionReturn(PETSC_SUCCESS);
2691 }
2692 
2693 /*@
2694   PetscClassIdRegister - Registers a new class name for objects and logging operations in an application code.
2695 
2696   Not Collective
2697 
2698   Input Parameter:
2699 . name - The class name
2700 
2701   Output Parameter:
2702 . oclass - The class id or classid
2703 
2704   Level: developer
2705 
2706 .seealso: [](ch_profiling), `PetscLogEventRegister()`
2707 @*/
2708 PetscErrorCode PetscClassIdRegister(const char name[], PetscClassId *oclass)
2709 {
2710   PetscFunctionBegin;
2711   *oclass = ++PETSC_LARGEST_CLASSID;
2712 #if defined(PETSC_USE_LOG)
2713   {
2714     PetscLogState state;
2715     PetscLogClass logclass;
2716 
2717     PetscCall(PetscLogGetState(&state));
2718     if (state) PetscCall(PetscLogStateClassRegister(state, name, *oclass, &logclass));
2719   }
2720 #endif
2721   PetscFunctionReturn(PETSC_SUCCESS);
2722 }
2723