Actual source code: lognested.c
1: #include <petscviewer.h>
2: #include "lognested.h"
3: #include "xmlviewer.h"
5: PETSC_INTERN PetscErrorCode PetscLogHandlerNestedSetThreshold(PetscLogHandler h, PetscLogDouble newThresh, PetscLogDouble *oldThresh)
6: {
7: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
9: PetscFunctionBegin;
10: if (oldThresh) *oldThresh = nested->threshold;
11: if (newThresh == (PetscLogDouble)PETSC_DECIDE) newThresh = 0.01;
12: if (newThresh == (PetscLogDouble)PETSC_DEFAULT) newThresh = 0.01;
13: nested->threshold = PetscMax(newThresh, 0.0);
14: PetscFunctionReturn(PETSC_SUCCESS);
15: }
17: static PetscErrorCode PetscLogEventGetNestedEvent(PetscLogHandler h, PetscLogEvent e, PetscLogEvent *nested_event)
18: {
19: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
20: NestedIdPair key;
21: PetscHashIter iter;
22: PetscBool missing;
23: PetscLogState state;
25: PetscFunctionBegin;
26: PetscCall(PetscLogHandlerGetState(h, &state));
27: PetscCall(PetscIntStackTop(nested->nested_stack, &key.root));
28: key.leaf = NestedIdFromEvent(e);
29: PetscCall(PetscNestedHashPut(nested->pair_map, key, &iter, &missing));
30: if (missing) {
31: // register a new nested event
32: char name[BUFSIZ];
33: PetscLogEventInfo event_info;
34: PetscLogEventInfo nested_event_info;
36: PetscCall(PetscLogStateEventGetInfo(state, e, &event_info));
37: PetscCall(PetscLogStateEventGetInfo(nested->state, key.root, &nested_event_info));
38: PetscCall(PetscSNPrintf(name, sizeof(name) - 1, "%s;%s", nested_event_info.name, event_info.name));
39: PetscCall(PetscLogStateEventRegister(nested->state, name, event_info.classid, nested_event));
40: PetscCall(PetscLogStateEventSetCollective(nested->state, *nested_event, event_info.collective));
41: PetscCall(PetscNestedHashIterSet(nested->pair_map, iter, *nested_event));
42: } else {
43: PetscCall(PetscNestedHashIterGet(nested->pair_map, iter, nested_event));
44: }
45: PetscFunctionReturn(PETSC_SUCCESS);
46: }
48: static PetscErrorCode PetscLogStageGetNestedEvent(PetscLogHandler h, PetscLogStage stage, PetscLogEvent *nested_event)
49: {
50: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
51: NestedIdPair key;
52: PetscHashIter iter;
53: PetscBool missing;
54: PetscLogState state;
56: PetscFunctionBegin;
57: PetscCall(PetscLogHandlerGetState(h, &state));
58: PetscCall(PetscIntStackTop(nested->nested_stack, &key.root));
59: key.leaf = NestedIdFromStage(stage);
60: PetscCall(PetscNestedHashPut(nested->pair_map, key, &iter, &missing));
61: if (missing) {
62: PetscLogStageInfo stage_info;
63: char name[BUFSIZ];
64: PetscBool collective = PETSC_TRUE;
66: PetscCall(PetscLogStateStageGetInfo(state, stage, &stage_info));
67: if (key.root >= 0) {
68: PetscLogEventInfo nested_event_info;
70: PetscCall(PetscLogStateEventGetInfo(nested->state, key.root, &nested_event_info));
71: PetscCall(PetscSNPrintf(name, sizeof(name) - 1, "%s;%s", nested_event_info.name, stage_info.name));
72: collective = nested_event_info.collective;
73: } else {
74: PetscCall(PetscSNPrintf(name, sizeof(name) - 1, "%s", stage_info.name));
75: }
76: PetscCall(PetscLogStateEventRegister(nested->state, name, nested->nested_stage_id, nested_event));
77: PetscCall(PetscLogStateEventSetCollective(nested->state, *nested_event, collective));
78: PetscCall(PetscNestedHashIterSet(nested->pair_map, iter, *nested_event));
79: } else {
80: PetscCall(PetscNestedHashIterGet(nested->pair_map, iter, nested_event));
81: }
82: PetscFunctionReturn(PETSC_SUCCESS);
83: }
85: static PetscErrorCode PetscLogNestedFindNestedId(PetscLogHandler h, NestedId orig_id, PetscInt *pop_count)
86: {
87: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
88: PetscInt count, i;
90: PetscFunctionBegin;
91: // stop before zero cause there is a null event at the bottom of the stack
92: for (i = nested->orig_stack->top, count = 0; i > 0; i--) {
93: count++;
94: if (nested->orig_stack->stack[i] == orig_id) break;
95: }
96: *pop_count = count;
97: if (count == 1) PetscFunctionReturn(PETSC_SUCCESS); // Normal function, just the top of the stack is being popped.
98: if (orig_id > 0) {
99: PetscLogEvent event_id = NestedIdToEvent(orig_id);
100: PetscLogState state;
101: PetscLogEventInfo event_info;
103: PetscCall(PetscLogHandlerGetState(h, &state));
104: PetscCall(PetscLogStateEventGetInfo(state, event_id, &event_info));
105: PetscCheck(i > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Tried to end event %s, but it is not in the event stack", event_info.name);
106: } else {
107: PetscLogStage stage_id = NestedIdToStage(orig_id);
108: PetscLogState state;
109: PetscLogStageInfo stage_info;
111: PetscCall(PetscLogHandlerGetState(h, &state));
112: PetscCall(PetscLogStateStageGetInfo(state, stage_id, &stage_info));
113: PetscCheck(i > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Tried to pop stage %s, but it is not in the stage stack", stage_info.name);
114: }
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: static PetscErrorCode PetscLogNestedCheckNested(PetscLogHandler h, NestedId leaf, PetscLogEvent nested_event)
119: {
120: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
121: NestedIdPair key;
122: NestedId val;
124: PetscFunctionBegin;
125: PetscCall(PetscIntStackTop(nested->nested_stack, &key.root));
126: key.leaf = leaf;
127: PetscCall(PetscNestedHashGet(nested->pair_map, key, &val));
128: PetscCheck(val == nested_event, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Logging events and stages are not nested, nested logging cannot be used");
129: PetscFunctionReturn(PETSC_SUCCESS);
130: }
132: static PetscErrorCode PetscLogHandlerEventBegin_Nested(PetscLogHandler h, PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)
133: {
134: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
135: PetscLogEvent nested_event;
137: PetscFunctionBegin;
138: PetscCall(PetscLogEventGetNestedEvent(h, e, &nested_event));
139: PetscCall(PetscLogHandlerEventBegin(nested->handler, nested_event, o1, o2, o3, o4));
140: PetscCall(PetscIntStackPush(nested->nested_stack, nested_event));
141: PetscCall(PetscIntStackPush(nested->orig_stack, NestedIdFromEvent(e)));
142: PetscFunctionReturn(PETSC_SUCCESS);
143: }
145: static PetscErrorCode PetscLogHandlerNestedEventEnd(PetscLogHandler h, NestedId id, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)
146: {
147: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
148: PetscInt pop_count;
150: PetscFunctionBegin;
151: PetscCall(PetscLogNestedFindNestedId(h, id, &pop_count));
152: for (PetscInt c = 0; c < pop_count; c++) {
153: PetscLogEvent nested_event;
154: PetscLogEvent nested_id;
156: PetscCall(PetscIntStackPop(nested->nested_stack, &nested_event));
157: PetscCall(PetscIntStackPop(nested->orig_stack, &nested_id));
158: if (PetscDefined(USE_DEBUG)) PetscCall(PetscLogNestedCheckNested(h, nested_id, nested_event));
159: if ((pop_count > 1) && (c + 1 < pop_count)) {
160: if (nested_id > 0) {
161: PetscLogEvent event_id = NestedIdToEvent(nested_id);
162: PetscLogState state;
163: PetscLogEventInfo event_info;
165: PetscCall(PetscLogHandlerGetState(h, &state));
166: PetscCall(PetscLogStateEventGetInfo(state, event_id, &event_info));
167: PetscCall(PetscInfo(h, "Log event %s wasn't ended, ending it to maintain stack property for nested log handler\n", event_info.name));
168: }
169: }
170: PetscCall(PetscLogHandlerEventEnd(nested->handler, nested_event, o1, o2, o3, o4));
171: }
172: PetscFunctionReturn(PETSC_SUCCESS);
173: }
175: static PetscErrorCode PetscLogHandlerEventEnd_Nested(PetscLogHandler h, PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4)
176: {
177: PetscFunctionBegin;
178: PetscCall(PetscLogHandlerNestedEventEnd(h, NestedIdFromEvent(e), o1, o2, o3, o4));
179: PetscFunctionReturn(PETSC_SUCCESS);
180: }
182: static PetscErrorCode PetscLogHandlerEventSync_Nested(PetscLogHandler h, PetscLogEvent e, MPI_Comm comm)
183: {
184: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
185: PetscLogEvent nested_event;
187: PetscFunctionBegin;
188: PetscCall(PetscLogEventGetNestedEvent(h, e, &nested_event));
189: PetscCall(PetscLogHandlerEventSync(nested->handler, nested_event, comm));
190: PetscFunctionReturn(PETSC_SUCCESS);
191: }
193: static PetscErrorCode PetscLogHandlerStagePush_Nested(PetscLogHandler h, PetscLogStage stage)
194: {
195: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
196: PetscLogEvent nested_event;
198: PetscFunctionBegin;
199: if (nested->nested_stage_id == -1) PetscCall(PetscClassIdRegister("LogNestedStage", &nested->nested_stage_id));
200: PetscCall(PetscLogStageGetNestedEvent(h, stage, &nested_event));
201: PetscCall(PetscLogHandlerEventBegin(nested->handler, nested_event, NULL, NULL, NULL, NULL));
202: PetscCall(PetscIntStackPush(nested->nested_stack, nested_event));
203: PetscCall(PetscIntStackPush(nested->orig_stack, NestedIdFromStage(stage)));
204: PetscFunctionReturn(PETSC_SUCCESS);
205: }
207: static PetscErrorCode PetscLogHandlerStagePop_Nested(PetscLogHandler h, PetscLogStage stage)
208: {
209: PetscFunctionBegin;
210: PetscCall(PetscLogHandlerNestedEventEnd(h, NestedIdFromStage(stage), NULL, NULL, NULL, NULL));
211: PetscFunctionReturn(PETSC_SUCCESS);
212: }
214: static PetscErrorCode PetscLogHandlerContextCreate_Nested(MPI_Comm comm, PetscLogHandler_Nested *nested_p)
215: {
216: PetscLogStage root_stage;
217: PetscLogHandler_Nested nested;
219: PetscFunctionBegin;
220: PetscCall(PetscNew(nested_p));
221: nested = *nested_p;
222: PetscCall(PetscLogStateCreate(&nested->state));
223: PetscCall(PetscIntStackCreate(&nested->nested_stack));
224: PetscCall(PetscIntStackCreate(&nested->orig_stack));
225: nested->nested_stage_id = -1;
226: nested->threshold = 0.01;
227: PetscCall(PetscNestedHashCreate(&nested->pair_map));
228: PetscCall(PetscLogHandlerCreate(comm, &nested->handler));
229: PetscCall(PetscLogHandlerSetType(nested->handler, PETSCLOGHANDLERDEFAULT));
230: PetscCall(PetscLogHandlerSetState(nested->handler, nested->state));
231: PetscCall(PetscLogStateStageRegister(nested->state, "", &root_stage));
232: PetscAssert(root_stage == 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "root stage not zero");
233: PetscCall(PetscLogHandlerStagePush(nested->handler, root_stage));
234: PetscCall(PetscLogStateStagePush(nested->state, root_stage));
235: PetscCall(PetscIntStackPush(nested->nested_stack, -1));
236: PetscCall(PetscIntStackPush(nested->orig_stack, -1));
237: PetscFunctionReturn(PETSC_SUCCESS);
238: }
240: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wconversion")
241: static PetscErrorCode PetscLogHandlerObjectCreate_Nested(PetscLogHandler h, PetscObject obj)
242: {
243: PetscClassId classid;
244: PetscInt num_registered, num_nested_registered;
245: PetscLogState state;
246: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
248: PetscFunctionBegin;
249: // register missing objects
250: PetscCall(PetscObjectGetClassId(obj, &classid));
251: PetscCall(PetscLogHandlerGetState(h, &state));
252: PetscCall(PetscLogStateGetNumClasses(nested->state, &num_nested_registered));
253: PetscCall(PetscLogStateGetNumClasses(state, &num_registered));
254: for (PetscLogClass c = num_nested_registered; c < num_registered; c++) {
255: PetscLogClassInfo class_info;
256: PetscLogClass nested_c;
258: PetscCall(PetscLogStateClassGetInfo(state, c, &class_info));
259: PetscCall(PetscLogStateClassRegister(nested->state, class_info.name, class_info.classid, &nested_c));
260: }
261: PetscCall(PetscLogHandlerObjectCreate(nested->handler, obj));
262: PetscFunctionReturn(PETSC_SUCCESS);
263: }
265: static PetscErrorCode PetscLogHandlerObjectDestroy_Nested(PetscLogHandler h, PetscObject obj)
266: {
267: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
269: PetscFunctionBegin;
270: PetscCall(PetscLogHandlerObjectDestroy(nested->handler, obj));
271: PetscFunctionReturn(PETSC_SUCCESS);
272: }
274: static PetscErrorCode PetscLogHandlerDestroy_Nested(PetscLogHandler h)
275: {
276: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)h->data;
278: PetscFunctionBegin;
279: PetscCall(PetscLogStateStagePop(nested->state));
280: PetscCall(PetscLogHandlerStagePop(nested->handler, 0));
281: PetscCall(PetscLogStateDestroy(&nested->state));
282: PetscCall(PetscIntStackDestroy(nested->nested_stack));
283: PetscCall(PetscIntStackDestroy(nested->orig_stack));
284: PetscCall(PetscNestedHashDestroy(&nested->pair_map));
285: PetscCall(PetscLogHandlerDestroy(&nested->handler));
286: PetscCall(PetscFree(nested));
287: PetscFunctionReturn(PETSC_SUCCESS);
288: }
290: static PetscErrorCode PetscLogNestedEventNodesOrderDepthFirst(PetscInt num_nodes, PetscInt parent, PetscNestedEventNode tree[], PetscInt *num_descendants)
291: {
292: PetscInt node, start_loc;
294: PetscFunctionBegin;
295: node = 0;
296: start_loc = 0;
297: while (node < num_nodes) {
298: if (tree[node].parent == parent) {
299: PetscInt num_this_descendants = 0;
300: PetscNestedEventNode tmp = tree[start_loc];
301: tree[start_loc] = tree[node];
302: tree[node] = tmp;
303: PetscCall(PetscLogNestedEventNodesOrderDepthFirst(num_nodes - start_loc - 1, tree[start_loc].id, &tree[start_loc + 1], &num_this_descendants));
304: tree[start_loc].num_descendants = num_this_descendants;
305: *num_descendants += 1 + num_this_descendants;
306: start_loc += 1 + num_this_descendants;
307: node = start_loc;
308: } else {
309: node++;
310: }
311: }
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
315: static PetscErrorCode PetscLogNestedCreatePerfNodes(MPI_Comm comm, PetscLogHandler_Nested nested, PetscLogGlobalNames global_events, PetscNestedEventNode **tree_p, PetscEventPerfInfo **perf_p)
316: {
317: PetscMPIInt size;
318: PetscInt num_nodes;
319: PetscInt num_map_entries;
320: PetscEventPerfInfo *perf;
321: NestedIdPair *keys;
322: NestedId *vals;
323: PetscInt offset;
324: PetscInt num_descendants;
325: PetscNestedEventNode *tree;
327: PetscFunctionBegin;
328: PetscCall(PetscLogGlobalNamesGetSize(global_events, NULL, &num_nodes));
329: PetscCall(PetscCalloc1(num_nodes, &tree));
330: for (PetscInt node = 0; node < num_nodes; node++) {
331: tree[node].id = node;
332: PetscCall(PetscLogGlobalNamesGlobalGetName(global_events, node, &tree[node].name));
333: tree[node].parent = -1;
334: }
335: PetscCall(PetscNestedHashGetSize(nested->pair_map, &num_map_entries));
336: PetscCall(PetscMalloc2(num_map_entries, &keys, num_map_entries, &vals));
337: offset = 0;
338: PetscCall(PetscNestedHashGetPairs(nested->pair_map, &offset, keys, vals));
339: for (PetscInt k = 0; k < num_map_entries; k++) {
340: NestedId root_local = keys[k].root;
341: NestedId leaf_local = vals[k];
342: PetscInt root_global;
343: PetscInt leaf_global;
345: PetscCall(PetscLogGlobalNamesLocalGetGlobal(global_events, leaf_local, &leaf_global));
346: if (root_local >= 0) {
347: PetscCall(PetscLogGlobalNamesLocalGetGlobal(global_events, root_local, &root_global));
348: tree[leaf_global].parent = root_global;
349: }
350: }
351: PetscCall(PetscFree2(keys, vals));
352: PetscCallMPI(MPI_Comm_size(comm, &size));
353: if (size > 1) { // get missing parents from other processes
354: PetscInt *parents;
356: PetscCall(PetscMalloc1(num_nodes, &parents));
357: for (PetscInt node = 0; node < num_nodes; node++) parents[node] = tree[node].parent;
358: PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, parents, num_nodes, MPIU_INT, MPI_MAX, comm));
359: for (PetscInt node = 0; node < num_nodes; node++) tree[node].parent = parents[node];
360: PetscCall(PetscFree(parents));
361: }
363: num_descendants = 0;
364: PetscCall(PetscLogNestedEventNodesOrderDepthFirst(num_nodes, -1, tree, &num_descendants));
365: PetscAssert(num_descendants == num_nodes, comm, PETSC_ERR_PLIB, "Failed tree ordering invariant");
367: PetscCall(PetscCalloc1(num_nodes, &perf));
368: for (PetscInt tree_node = 0; tree_node < num_nodes; tree_node++) {
369: PetscInt global_id = tree[tree_node].id;
370: PetscInt event_id;
372: PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_events, global_id, &event_id));
373: if (event_id >= 0) {
374: PetscEventPerfInfo *event_info;
376: PetscCall(PetscLogHandlerGetEventPerfInfo(nested->handler, 0, event_id, &event_info));
377: perf[tree_node] = *event_info;
378: } else {
379: PetscCall(PetscArrayzero(&perf[tree_node], 1));
380: }
381: }
382: *tree_p = tree;
383: *perf_p = perf;
384: PetscFunctionReturn(PETSC_SUCCESS);
385: }
386: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
388: static PetscErrorCode PetscLogHandlerView_Nested(PetscLogHandler handler, PetscViewer viewer)
389: {
390: PetscLogHandler_Nested nested = (PetscLogHandler_Nested)handler->data;
391: PetscNestedEventNode *nodes;
392: PetscEventPerfInfo *perf;
393: PetscLogGlobalNames global_events;
394: PetscNestedEventTree tree;
395: PetscViewerFormat format;
396: MPI_Comm comm = PetscObjectComm((PetscObject)viewer);
398: PetscFunctionBegin;
399: PetscCall(PetscLogRegistryCreateGlobalEventNames(comm, nested->state->registry, &global_events));
400: PetscCall(PetscLogNestedCreatePerfNodes(comm, nested, global_events, &nodes, &perf));
401: tree.comm = comm;
402: tree.global_events = global_events;
403: tree.perf = perf;
404: tree.nodes = nodes;
405: PetscCall(PetscViewerGetFormat(viewer, &format));
406: if (format == PETSC_VIEWER_ASCII_XML) {
407: PetscCall(PetscLogHandlerView_Nested_XML(nested, &tree, viewer));
408: } else if (format == PETSC_VIEWER_ASCII_FLAMEGRAPH) {
409: PetscCall(PetscLogHandlerView_Nested_Flamegraph(nested, &tree, viewer));
410: } else SETERRQ(comm, PETSC_ERR_ARG_INCOMP, "No nested viewer for this format");
411: PetscCall(PetscLogGlobalNamesDestroy(&global_events));
412: PetscCall(PetscFree(tree.nodes));
413: PetscCall(PetscFree(tree.perf));
414: PetscFunctionReturn(PETSC_SUCCESS);
415: }
417: /*MC
418: PETSCLOGHANDLERNESTED - PETSCLOGHANDLERNESTED = "nested" - A `PetscLogHandler` that collects data for PETSc's
419: XML and flamegraph profiling log viewers. A log handler of this type is created and started by
420: by `PetscLogNestedBegin()`.
422: Level: developer
424: .seealso: [](ch_profiling), `PetscLogHandler`
425: M*/
427: PETSC_INTERN PetscErrorCode PetscLogHandlerCreate_Nested(PetscLogHandler handler)
428: {
429: PetscFunctionBegin;
430: PetscCall(PetscLogHandlerContextCreate_Nested(PetscObjectComm((PetscObject)handler), (PetscLogHandler_Nested *)&handler->data));
431: handler->ops->destroy = PetscLogHandlerDestroy_Nested;
432: handler->ops->stagepush = PetscLogHandlerStagePush_Nested;
433: handler->ops->stagepop = PetscLogHandlerStagePop_Nested;
434: handler->ops->eventbegin = PetscLogHandlerEventBegin_Nested;
435: handler->ops->eventend = PetscLogHandlerEventEnd_Nested;
436: handler->ops->eventsync = PetscLogHandlerEventSync_Nested;
437: handler->ops->objectcreate = PetscLogHandlerObjectCreate_Nested;
438: handler->ops->objectdestroy = PetscLogHandlerObjectDestroy_Nested;
439: handler->ops->view = PetscLogHandlerView_Nested;
440: PetscFunctionReturn(PETSC_SUCCESS);
441: }