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