1 2 /* 3 PETSc code to log object creation and destruction and PETSc events. 4 5 This provides the public API used by the rest of PETSc and by users. 6 7 These routines use a private API that is not used elsewhere in PETSc and is not 8 accessible to users. The private API is defined in logimpl.h and the utils directory. 9 10 */ 11 #include <petsc/private/logimpl.h> /*I "petscsys.h" I*/ 12 #include <petsctime.h> 13 #include <petscviewer.h> 14 15 PetscErrorCode PetscLogObjectParent(PetscObject p,PetscObject c) 16 { 17 if (!c || !p) return 0; 18 c->parent = p; 19 c->parentid = p->id; 20 return 0; 21 } 22 23 /*@C 24 PetscLogObjectMemory - Adds to an object a count of additional amount of memory that is used by the object. 25 26 Not collective. 27 28 Input Parameters: 29 + obj - the PETSc object 30 - mem - the amount of memory that is being added to the object 31 32 Level: developer 33 34 Developer Notes: 35 Currently we do not always do a good job of associating all memory allocations with an object. 36 37 .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscInitializeNoArguments() 38 39 @*/ 40 PetscErrorCode PetscLogObjectMemory(PetscObject p,PetscLogDouble m) 41 { 42 if (!p) return 0; 43 p->mem += m; 44 return 0; 45 } 46 47 PetscLogEvent PETSC_LARGEST_EVENT = PETSC_EVENT; 48 49 #if defined(PETSC_USE_LOG) 50 #include <petscmachineinfo.h> 51 #include <petscconfiginfo.h> 52 53 /* used in the MPI_XXX() count macros in petsclog.h */ 54 55 /* Action and object logging variables */ 56 Action *petsc_actions = NULL; 57 Object *petsc_objects = NULL; 58 PetscBool petsc_logActions = PETSC_FALSE; 59 PetscBool petsc_logObjects = PETSC_FALSE; 60 int petsc_numActions = 0, petsc_maxActions = 100; 61 int petsc_numObjects = 0, petsc_maxObjects = 100; 62 int petsc_numObjectsDestroyed = 0; 63 64 /* Global counters */ 65 PetscLogDouble petsc_BaseTime = 0.0; 66 PetscLogDouble petsc_TotalFlops = 0.0; /* The number of flops */ 67 PetscLogDouble petsc_tmp_flops = 0.0; /* The incremental number of flops */ 68 PetscLogDouble petsc_send_ct = 0.0; /* The number of sends */ 69 PetscLogDouble petsc_recv_ct = 0.0; /* The number of receives */ 70 PetscLogDouble petsc_send_len = 0.0; /* The total length of all sent messages */ 71 PetscLogDouble petsc_recv_len = 0.0; /* The total length of all received messages */ 72 PetscLogDouble petsc_isend_ct = 0.0; /* The number of immediate sends */ 73 PetscLogDouble petsc_irecv_ct = 0.0; /* The number of immediate receives */ 74 PetscLogDouble petsc_isend_len = 0.0; /* The total length of all immediate send messages */ 75 PetscLogDouble petsc_irecv_len = 0.0; /* The total length of all immediate receive messages */ 76 PetscLogDouble petsc_wait_ct = 0.0; /* The number of waits */ 77 PetscLogDouble petsc_wait_any_ct = 0.0; /* The number of anywaits */ 78 PetscLogDouble petsc_wait_all_ct = 0.0; /* The number of waitalls */ 79 PetscLogDouble petsc_sum_of_waits_ct = 0.0; /* The total number of waits */ 80 PetscLogDouble petsc_allreduce_ct = 0.0; /* The number of reductions */ 81 PetscLogDouble petsc_gather_ct = 0.0; /* The number of gathers and gathervs */ 82 PetscLogDouble petsc_scatter_ct = 0.0; /* The number of scatters and scattervs */ 83 #if defined(PETSC_HAVE_DEVICE) 84 PetscLogDouble petsc_ctog_ct = 0.0; /* The total number of CPU to GPU copies */ 85 PetscLogDouble petsc_gtoc_ct = 0.0; /* The total number of GPU to CPU copies */ 86 PetscLogDouble petsc_ctog_sz = 0.0; /* The total size of CPU to GPU copies */ 87 PetscLogDouble petsc_gtoc_sz = 0.0; /* The total size of GPU to CPU copies */ 88 PetscLogDouble petsc_ctog_ct_scalar = 0.0; /* The total number of CPU to GPU copies */ 89 PetscLogDouble petsc_gtoc_ct_scalar = 0.0; /* The total number of GPU to CPU copies */ 90 PetscLogDouble petsc_ctog_sz_scalar = 0.0; /* The total size of CPU to GPU copies */ 91 PetscLogDouble petsc_gtoc_sz_scalar = 0.0; /* The total size of GPU to CPU copies */ 92 PetscLogDouble petsc_gflops = 0.0; /* The flops done on a GPU */ 93 PetscLogDouble petsc_gtime = 0.0; /* The time spent on a GPU */ 94 #endif 95 96 /* Logging functions */ 97 PetscErrorCode (*PetscLogPHC)(PetscObject) = NULL; 98 PetscErrorCode (*PetscLogPHD)(PetscObject) = NULL; 99 PetscErrorCode (*PetscLogPLB)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL; 100 PetscErrorCode (*PetscLogPLE)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL; 101 102 /* Tracing event logging variables */ 103 FILE *petsc_tracefile = NULL; 104 int petsc_tracelevel = 0; 105 const char *petsc_traceblanks = " "; 106 char petsc_tracespace[128] = " "; 107 PetscLogDouble petsc_tracetime = 0.0; 108 static PetscBool PetscLogInitializeCalled = PETSC_FALSE; 109 110 PETSC_INTERN PetscErrorCode PetscLogInitialize(void) 111 { 112 int stage; 113 PetscBool opt; 114 115 PetscFunctionBegin; 116 if (PetscLogInitializeCalled) PetscFunctionReturn(0); 117 PetscLogInitializeCalled = PETSC_TRUE; 118 119 PetscCall(PetscOptionsHasName(NULL,NULL, "-log_exclude_actions", &opt)); 120 if (opt) petsc_logActions = PETSC_FALSE; 121 PetscCall(PetscOptionsHasName(NULL,NULL, "-log_exclude_objects", &opt)); 122 if (opt) petsc_logObjects = PETSC_FALSE; 123 if (petsc_logActions) { 124 PetscCall(PetscMalloc1(petsc_maxActions, &petsc_actions)); 125 } 126 if (petsc_logObjects) { 127 PetscCall(PetscMalloc1(petsc_maxObjects, &petsc_objects)); 128 } 129 PetscLogPHC = PetscLogObjCreateDefault; 130 PetscLogPHD = PetscLogObjDestroyDefault; 131 /* Setup default logging structures */ 132 PetscCall(PetscStageLogCreate(&petsc_stageLog)); 133 PetscCall(PetscStageLogRegister(petsc_stageLog, "Main Stage", &stage)); 134 135 /* All processors sync here for more consistent logging */ 136 PetscCallMPI(MPI_Barrier(PETSC_COMM_WORLD)); 137 PetscTime(&petsc_BaseTime); 138 PetscCall(PetscLogStagePush(stage)); 139 PetscFunctionReturn(0); 140 } 141 142 PETSC_INTERN PetscErrorCode PetscLogFinalize(void) 143 { 144 PetscStageLog stageLog; 145 146 PetscFunctionBegin; 147 PetscCall(PetscFree(petsc_actions)); 148 PetscCall(PetscFree(petsc_objects)); 149 PetscCall(PetscLogNestedEnd()); 150 PetscCall(PetscLogSet(NULL, NULL)); 151 152 /* Resetting phase */ 153 PetscCall(PetscLogGetStageLog(&stageLog)); 154 PetscCall(PetscStageLogDestroy(stageLog)); 155 156 petsc_TotalFlops = 0.0; 157 petsc_numActions = 0; 158 petsc_numObjects = 0; 159 petsc_numObjectsDestroyed = 0; 160 petsc_maxActions = 100; 161 petsc_maxObjects = 100; 162 petsc_actions = NULL; 163 petsc_objects = NULL; 164 petsc_logActions = PETSC_FALSE; 165 petsc_logObjects = PETSC_FALSE; 166 petsc_BaseTime = 0.0; 167 petsc_TotalFlops = 0.0; 168 petsc_tmp_flops = 0.0; 169 petsc_send_ct = 0.0; 170 petsc_recv_ct = 0.0; 171 petsc_send_len = 0.0; 172 petsc_recv_len = 0.0; 173 petsc_isend_ct = 0.0; 174 petsc_irecv_ct = 0.0; 175 petsc_isend_len = 0.0; 176 petsc_irecv_len = 0.0; 177 petsc_wait_ct = 0.0; 178 petsc_wait_any_ct = 0.0; 179 petsc_wait_all_ct = 0.0; 180 petsc_sum_of_waits_ct = 0.0; 181 petsc_allreduce_ct = 0.0; 182 petsc_gather_ct = 0.0; 183 petsc_scatter_ct = 0.0; 184 #if defined(PETSC_HAVE_DEVICE) 185 petsc_ctog_ct = 0.0; 186 petsc_gtoc_ct = 0.0; 187 petsc_ctog_sz = 0.0; 188 petsc_gtoc_sz = 0.0; 189 petsc_gflops = 0.0; 190 petsc_gtime = 0.0; 191 #endif 192 PETSC_LARGEST_EVENT = PETSC_EVENT; 193 PetscLogPHC = NULL; 194 PetscLogPHD = NULL; 195 petsc_tracefile = NULL; 196 petsc_tracelevel = 0; 197 petsc_traceblanks = " "; 198 petsc_tracespace[0] = ' '; petsc_tracespace[1] = 0; 199 petsc_tracetime = 0.0; 200 PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID; 201 PETSC_OBJECT_CLASSID = 0; 202 petsc_stageLog = NULL; 203 PetscLogInitializeCalled = PETSC_FALSE; 204 PetscFunctionReturn(0); 205 } 206 207 /*@C 208 PetscLogSet - Sets the logging functions called at the beginning and ending of every event. 209 210 Not Collective 211 212 Input Parameters: 213 + b - The function called at beginning of event 214 - e - The function called at end of event 215 216 Level: developer 217 218 .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogTraceBegin() 219 @*/ 220 PetscErrorCode PetscLogSet(PetscErrorCode (*b)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject), 221 PetscErrorCode (*e)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject)) 222 { 223 PetscFunctionBegin; 224 PetscLogPLB = b; 225 PetscLogPLE = e; 226 PetscFunctionReturn(0); 227 } 228 229 /*@C 230 PetscLogIsActive - Check if logging is currently in progress. 231 232 Not Collective 233 234 Output Parameter: 235 . isActive - PETSC_TRUE if logging is in progress, PETSC_FALSE otherwise 236 237 Level: beginner 238 239 .seealso: PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogSet() 240 @*/ 241 PetscErrorCode PetscLogIsActive(PetscBool *isActive) 242 { 243 PetscFunctionBegin; 244 *isActive = (PetscLogPLB && PetscLogPLE) ? PETSC_TRUE : PETSC_FALSE; 245 PetscFunctionReturn(0); 246 } 247 248 /*@C 249 PetscLogDefaultBegin - Turns on logging of objects and events. This logs flop 250 rates and object creation and should not slow programs down too much. 251 This routine may be called more than once. 252 253 Logically Collective over PETSC_COMM_WORLD 254 255 Options Database Keys: 256 . -log_view [viewertype:filename:viewerformat] - Prints summary of flop and timing information to the 257 screen (for code configured with --with-log=1 (which is the default)) 258 259 Usage: 260 .vb 261 PetscInitialize(...); 262 PetscLogDefaultBegin(); 263 ... code ... 264 PetscLogView(viewer); or PetscLogDump(); 265 PetscFinalize(); 266 .ve 267 268 Notes: 269 PetscLogView(viewer) or PetscLogDump() actually cause the printing of 270 the logging information. 271 272 Level: advanced 273 274 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogTraceBegin() 275 @*/ 276 PetscErrorCode PetscLogDefaultBegin(void) 277 { 278 PetscFunctionBegin; 279 PetscCall(PetscLogSet(PetscLogEventBeginDefault, PetscLogEventEndDefault)); 280 PetscFunctionReturn(0); 281 } 282 283 /*@C 284 PetscLogAllBegin - Turns on extensive logging of objects and events. Logs 285 all events. This creates large log files and slows the program down. 286 287 Logically Collective on PETSC_COMM_WORLD 288 289 Options Database Keys: 290 . -log_all - Prints extensive log information 291 292 Usage: 293 .vb 294 PetscInitialize(...); 295 PetscLogAllBegin(); 296 ... code ... 297 PetscLogDump(filename); 298 PetscFinalize(); 299 .ve 300 301 Notes: 302 A related routine is PetscLogDefaultBegin() (with the options key -log), which is 303 intended for production runs since it logs only flop rates and object 304 creation (and shouldn't significantly slow the programs). 305 306 Level: advanced 307 308 .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogTraceBegin() 309 @*/ 310 PetscErrorCode PetscLogAllBegin(void) 311 { 312 PetscFunctionBegin; 313 PetscCall(PetscLogSet(PetscLogEventBeginComplete, PetscLogEventEndComplete)); 314 PetscFunctionReturn(0); 315 } 316 317 /*@C 318 PetscLogTraceBegin - Activates trace logging. Every time a PETSc event 319 begins or ends, the event name is printed. 320 321 Logically Collective on PETSC_COMM_WORLD 322 323 Input Parameter: 324 . file - The file to print trace in (e.g. stdout) 325 326 Options Database Key: 327 . -log_trace [filename] - Activates PetscLogTraceBegin() 328 329 Notes: 330 PetscLogTraceBegin() prints the processor number, the execution time (sec), 331 then "Event begin:" or "Event end:" followed by the event name. 332 333 PetscLogTraceBegin() allows tracing of all PETSc calls, which is useful 334 to determine where a program is hanging without running in the 335 debugger. Can be used in conjunction with the -info option. 336 337 Level: intermediate 338 339 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogDefaultBegin() 340 @*/ 341 PetscErrorCode PetscLogTraceBegin(FILE *file) 342 { 343 PetscFunctionBegin; 344 petsc_tracefile = file; 345 346 PetscCall(PetscLogSet(PetscLogEventBeginTrace, PetscLogEventEndTrace)); 347 PetscFunctionReturn(0); 348 } 349 350 /*@ 351 PetscLogActions - Determines whether actions are logged for the graphical viewer. 352 353 Not Collective 354 355 Input Parameter: 356 . flag - PETSC_TRUE if actions are to be logged 357 358 Level: intermediate 359 360 Note: Logging of actions continues to consume more memory as the program 361 runs. Long running programs should consider turning this feature off. 362 363 Options Database Keys: 364 . -log_exclude_actions - Turns off actions logging 365 366 .seealso: PetscLogStagePush(), PetscLogStagePop() 367 @*/ 368 PetscErrorCode PetscLogActions(PetscBool flag) 369 { 370 PetscFunctionBegin; 371 petsc_logActions = flag; 372 PetscFunctionReturn(0); 373 } 374 375 /*@ 376 PetscLogObjects - Determines whether objects are logged for the graphical viewer. 377 378 Not Collective 379 380 Input Parameter: 381 . flag - PETSC_TRUE if objects are to be logged 382 383 Level: intermediate 384 385 Note: Logging of objects continues to consume more memory as the program 386 runs. Long running programs should consider turning this feature off. 387 388 Options Database Keys: 389 . -log_exclude_objects - Turns off objects logging 390 391 .seealso: PetscLogStagePush(), PetscLogStagePop() 392 @*/ 393 PetscErrorCode PetscLogObjects(PetscBool flag) 394 { 395 PetscFunctionBegin; 396 petsc_logObjects = flag; 397 PetscFunctionReturn(0); 398 } 399 400 /*------------------------------------------------ Stage Functions --------------------------------------------------*/ 401 /*@C 402 PetscLogStageRegister - Attaches a character string name to a logging stage. 403 404 Not Collective 405 406 Input Parameter: 407 . sname - The name to associate with that stage 408 409 Output Parameter: 410 . stage - The stage number 411 412 Level: intermediate 413 414 .seealso: PetscLogStagePush(), PetscLogStagePop() 415 @*/ 416 PetscErrorCode PetscLogStageRegister(const char sname[],PetscLogStage *stage) 417 { 418 PetscStageLog stageLog; 419 PetscLogEvent event; 420 421 PetscFunctionBegin; 422 PetscCall(PetscLogGetStageLog(&stageLog)); 423 PetscCall(PetscStageLogRegister(stageLog, sname, stage)); 424 /* Copy events already changed in the main stage, this sucks */ 425 PetscCall(PetscEventPerfLogEnsureSize(stageLog->stageInfo[*stage].eventLog, stageLog->eventLog->numEvents)); 426 for (event = 0; event < stageLog->eventLog->numEvents; event++) { 427 PetscCall(PetscEventPerfInfoCopy(&stageLog->stageInfo[0].eventLog->eventInfo[event],&stageLog->stageInfo[*stage].eventLog->eventInfo[event])); 428 } 429 PetscCall(PetscClassPerfLogEnsureSize(stageLog->stageInfo[*stage].classLog, stageLog->classLog->numClasses)); 430 PetscFunctionReturn(0); 431 } 432 433 /*@C 434 PetscLogStagePush - This function pushes a stage on the stack. 435 436 Not Collective 437 438 Input Parameter: 439 . stage - The stage on which to log 440 441 Usage: 442 If the option -log_sumary is used to run the program containing the 443 following code, then 2 sets of summary data will be printed during 444 PetscFinalize(). 445 .vb 446 PetscInitialize(int *argc,char ***args,0,0); 447 [stage 0 of code] 448 PetscLogStagePush(1); 449 [stage 1 of code] 450 PetscLogStagePop(); 451 PetscBarrier(...); 452 [more stage 0 of code] 453 PetscFinalize(); 454 .ve 455 456 Notes: 457 Use PetscLogStageRegister() to register a stage. 458 459 Level: intermediate 460 461 .seealso: PetscLogStagePop(), PetscLogStageRegister(), PetscBarrier() 462 @*/ 463 PetscErrorCode PetscLogStagePush(PetscLogStage stage) 464 { 465 PetscStageLog stageLog; 466 467 PetscFunctionBegin; 468 PetscCall(PetscLogGetStageLog(&stageLog)); 469 PetscCall(PetscStageLogPush(stageLog, stage)); 470 PetscFunctionReturn(0); 471 } 472 473 /*@C 474 PetscLogStagePop - This function pops a stage from the stack. 475 476 Not Collective 477 478 Usage: 479 If the option -log_sumary is used to run the program containing the 480 following code, then 2 sets of summary data will be printed during 481 PetscFinalize(). 482 .vb 483 PetscInitialize(int *argc,char ***args,0,0); 484 [stage 0 of code] 485 PetscLogStagePush(1); 486 [stage 1 of code] 487 PetscLogStagePop(); 488 PetscBarrier(...); 489 [more stage 0 of code] 490 PetscFinalize(); 491 .ve 492 493 Notes: 494 Use PetscLogStageRegister() to register a stage. 495 496 Level: intermediate 497 498 .seealso: PetscLogStagePush(), PetscLogStageRegister(), PetscBarrier() 499 @*/ 500 PetscErrorCode PetscLogStagePop(void) 501 { 502 PetscStageLog stageLog; 503 504 PetscFunctionBegin; 505 PetscCall(PetscLogGetStageLog(&stageLog)); 506 PetscCall(PetscStageLogPop(stageLog)); 507 PetscFunctionReturn(0); 508 } 509 510 /*@ 511 PetscLogStageSetActive - Determines stage activity for PetscLogEventBegin() and PetscLogEventEnd(). 512 513 Not Collective 514 515 Input Parameters: 516 + stage - The stage 517 - isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE) 518 519 Level: intermediate 520 521 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 522 @*/ 523 PetscErrorCode PetscLogStageSetActive(PetscLogStage stage, PetscBool isActive) 524 { 525 PetscStageLog stageLog; 526 527 PetscFunctionBegin; 528 PetscCall(PetscLogGetStageLog(&stageLog)); 529 PetscCall(PetscStageLogSetActive(stageLog, stage, isActive)); 530 PetscFunctionReturn(0); 531 } 532 533 /*@ 534 PetscLogStageGetActive - Returns stage activity for PetscLogEventBegin() and PetscLogEventEnd(). 535 536 Not Collective 537 538 Input Parameter: 539 . stage - The stage 540 541 Output Parameter: 542 . isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE) 543 544 Level: intermediate 545 546 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 547 @*/ 548 PetscErrorCode PetscLogStageGetActive(PetscLogStage stage, PetscBool *isActive) 549 { 550 PetscStageLog stageLog; 551 552 PetscFunctionBegin; 553 PetscCall(PetscLogGetStageLog(&stageLog)); 554 PetscCall(PetscStageLogGetActive(stageLog, stage, isActive)); 555 PetscFunctionReturn(0); 556 } 557 558 /*@ 559 PetscLogStageSetVisible - Determines stage visibility in PetscLogView() 560 561 Not Collective 562 563 Input Parameters: 564 + stage - The stage 565 - isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE) 566 567 Level: intermediate 568 569 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView() 570 @*/ 571 PetscErrorCode PetscLogStageSetVisible(PetscLogStage stage, PetscBool isVisible) 572 { 573 PetscStageLog stageLog; 574 575 PetscFunctionBegin; 576 PetscCall(PetscLogGetStageLog(&stageLog)); 577 PetscCall(PetscStageLogSetVisible(stageLog, stage, isVisible)); 578 PetscFunctionReturn(0); 579 } 580 581 /*@ 582 PetscLogStageGetVisible - Returns stage visibility in PetscLogView() 583 584 Not Collective 585 586 Input Parameter: 587 . stage - The stage 588 589 Output Parameter: 590 . isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE) 591 592 Level: intermediate 593 594 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView() 595 @*/ 596 PetscErrorCode PetscLogStageGetVisible(PetscLogStage stage, PetscBool *isVisible) 597 { 598 PetscStageLog stageLog; 599 600 PetscFunctionBegin; 601 PetscCall(PetscLogGetStageLog(&stageLog)); 602 PetscCall(PetscStageLogGetVisible(stageLog, stage, isVisible)); 603 PetscFunctionReturn(0); 604 } 605 606 /*@C 607 PetscLogStageGetId - Returns the stage id when given the stage name. 608 609 Not Collective 610 611 Input Parameter: 612 . name - The stage name 613 614 Output Parameter: 615 . stage - The stage, , or -1 if no stage with that name exists 616 617 Level: intermediate 618 619 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 620 @*/ 621 PetscErrorCode PetscLogStageGetId(const char name[], PetscLogStage *stage) 622 { 623 PetscStageLog stageLog; 624 625 PetscFunctionBegin; 626 PetscCall(PetscLogGetStageLog(&stageLog)); 627 PetscCall(PetscStageLogGetStage(stageLog, name, stage)); 628 PetscFunctionReturn(0); 629 } 630 631 /*------------------------------------------------ Event Functions --------------------------------------------------*/ 632 /*@C 633 PetscLogEventRegister - Registers an event name for logging operations in an application code. 634 635 Not Collective 636 637 Input Parameters: 638 + name - The name associated with the event 639 - classid - The classid associated to the class for this event, obtain either with 640 PetscClassIdRegister() or use a predefined one such as KSP_CLASSID, SNES_CLASSID, the predefined ones 641 are only available in C code 642 643 Output Parameter: 644 . event - The event id for use with PetscLogEventBegin() and PetscLogEventEnd(). 645 646 Example of Usage: 647 .vb 648 PetscLogEvent USER_EVENT; 649 PetscClassId classid; 650 PetscLogDouble user_event_flops; 651 PetscClassIdRegister("class name",&classid); 652 PetscLogEventRegister("User event name",classid,&USER_EVENT); 653 PetscLogEventBegin(USER_EVENT,0,0,0,0); 654 [code segment to monitor] 655 PetscLogFlops(user_event_flops); 656 PetscLogEventEnd(USER_EVENT,0,0,0,0); 657 .ve 658 659 Notes: 660 PETSc automatically logs library events if the code has been 661 configured with --with-log (which is the default) and 662 -log_view or -log_all is specified. PetscLogEventRegister() is 663 intended for logging user events to supplement this PETSc 664 information. 665 666 PETSc can gather data for use with the utilities Jumpshot 667 (part of the MPICH distribution). If PETSc has been compiled 668 with flag -DPETSC_HAVE_MPE (MPE is an additional utility within 669 MPICH), the user can employ another command line option, -log_mpe, 670 to create a logfile, "mpe.log", which can be visualized 671 Jumpshot. 672 673 The classid is associated with each event so that classes of events 674 can be disabled simultaneously, such as all matrix events. The user 675 can either use an existing classid, such as MAT_CLASSID, or create 676 their own as shown in the example. 677 678 If an existing event with the same name exists, its event handle is 679 returned instead of creating a new event. 680 681 Level: intermediate 682 683 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogFlops(), 684 PetscLogEventActivate(), PetscLogEventDeactivate(), PetscClassIdRegister() 685 @*/ 686 PetscErrorCode PetscLogEventRegister(const char name[],PetscClassId classid,PetscLogEvent *event) 687 { 688 PetscStageLog stageLog; 689 int stage; 690 691 PetscFunctionBegin; 692 *event = PETSC_DECIDE; 693 PetscCall(PetscLogGetStageLog(&stageLog)); 694 PetscCall(PetscEventRegLogGetEvent(stageLog->eventLog, name, event)); 695 if (*event > 0) PetscFunctionReturn(0); 696 PetscCall(PetscEventRegLogRegister(stageLog->eventLog, name, classid, event)); 697 for (stage = 0; stage < stageLog->numStages; stage++) { 698 PetscCall(PetscEventPerfLogEnsureSize(stageLog->stageInfo[stage].eventLog, stageLog->eventLog->numEvents)); 699 PetscCall(PetscClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses)); 700 } 701 PetscFunctionReturn(0); 702 } 703 704 /*@ 705 PetscLogEventSetCollective - Indicates that a particular event is collective. 706 707 Not Collective 708 709 Input Parameters: 710 + event - The event id 711 - collective - Bolean flag indicating whether a particular event is collective 712 713 Note: 714 New events returned from PetscLogEventRegister() are collective by default. 715 716 Level: developer 717 718 .seealso: PetscLogEventRegister() 719 @*/ 720 PetscErrorCode PetscLogEventSetCollective(PetscLogEvent event,PetscBool collective) 721 { 722 PetscStageLog stageLog; 723 PetscEventRegLog eventRegLog; 724 725 PetscFunctionBegin; 726 PetscCall(PetscLogGetStageLog(&stageLog)); 727 PetscCall(PetscStageLogGetEventRegLog(stageLog,&eventRegLog)); 728 PetscCheckFalse(event < 0 || event > eventRegLog->numEvents,PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid event id"); 729 eventRegLog->eventInfo[event].collective = collective; 730 PetscFunctionReturn(0); 731 } 732 733 /*@ 734 PetscLogEventIncludeClass - Activates event logging for a PETSc object class in every stage. 735 736 Not Collective 737 738 Input Parameter: 739 . classid - The object class, for example MAT_CLASSID, SNES_CLASSID, etc. 740 741 Level: developer 742 743 .seealso: PetscLogEventActivateClass(),PetscLogEventDeactivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate() 744 @*/ 745 PetscErrorCode PetscLogEventIncludeClass(PetscClassId classid) 746 { 747 PetscStageLog stageLog; 748 int stage; 749 750 PetscFunctionBegin; 751 PetscCall(PetscLogGetStageLog(&stageLog)); 752 for (stage = 0; stage < stageLog->numStages; stage++) { 753 PetscCall(PetscEventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid)); 754 } 755 PetscFunctionReturn(0); 756 } 757 758 /*@ 759 PetscLogEventExcludeClass - Deactivates event logging for a PETSc object class in every stage. 760 761 Not Collective 762 763 Input Parameter: 764 . classid - The object class, for example MAT_CLASSID, SNES_CLASSID, etc. 765 766 Level: developer 767 768 .seealso: PetscLogEventDeactivateClass(),PetscLogEventActivateClass(),PetscLogEventDeactivate(),PetscLogEventActivate() 769 @*/ 770 PetscErrorCode PetscLogEventExcludeClass(PetscClassId classid) 771 { 772 PetscStageLog stageLog; 773 int stage; 774 775 PetscFunctionBegin; 776 PetscCall(PetscLogGetStageLog(&stageLog)); 777 for (stage = 0; stage < stageLog->numStages; stage++) { 778 PetscCall(PetscEventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid)); 779 } 780 PetscFunctionReturn(0); 781 } 782 783 /*@ 784 PetscLogEventActivate - Indicates that a particular event should be logged. 785 786 Not Collective 787 788 Input Parameter: 789 . event - The event id 790 791 Usage: 792 .vb 793 PetscLogEventDeactivate(VEC_SetValues); 794 [code where you do not want to log VecSetValues()] 795 PetscLogEventActivate(VEC_SetValues); 796 [code where you do want to log VecSetValues()] 797 .ve 798 799 Note: 800 The event may be either a pre-defined PETSc event (found in include/petsclog.h) 801 or an event number obtained with PetscLogEventRegister(). 802 803 Level: advanced 804 805 .seealso: PlogEventDeactivate(), PlogEventDeactivatePush(), PetscLogEventDeactivatePop() 806 @*/ 807 PetscErrorCode PetscLogEventActivate(PetscLogEvent event) 808 { 809 PetscStageLog stageLog; 810 int stage; 811 812 PetscFunctionBegin; 813 PetscCall(PetscLogGetStageLog(&stageLog)); 814 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 815 PetscCall(PetscEventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event)); 816 PetscFunctionReturn(0); 817 } 818 819 /*@ 820 PetscLogEventDeactivate - Indicates that a particular event should not be logged. 821 822 Not Collective 823 824 Input Parameter: 825 . event - The event id 826 827 Usage: 828 .vb 829 PetscLogEventDeactivate(VEC_SetValues); 830 [code where you do not want to log VecSetValues()] 831 PetscLogEventActivate(VEC_SetValues); 832 [code where you do want to log VecSetValues()] 833 .ve 834 835 Note: 836 The event may be either a pre-defined PETSc event (found in 837 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 838 839 Level: advanced 840 841 .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePush(), PetscLogEventDeactivatePop() 842 @*/ 843 PetscErrorCode PetscLogEventDeactivate(PetscLogEvent event) 844 { 845 PetscStageLog stageLog; 846 int stage; 847 848 PetscFunctionBegin; 849 PetscCall(PetscLogGetStageLog(&stageLog)); 850 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 851 PetscCall(PetscEventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event)); 852 PetscFunctionReturn(0); 853 } 854 855 /*@ 856 PetscLogEventDeactivatePush - Indicates that a particular event should not be logged. 857 858 Not Collective 859 860 Input Parameter: 861 . event - The event id 862 863 Usage: 864 .vb 865 PetscLogEventDeactivatePush(VEC_SetValues); 866 [code where you do not want to log VecSetValues()] 867 PetscLogEventDeactivatePop(VEC_SetValues); 868 [code where you do want to log VecSetValues()] 869 .ve 870 871 Note: 872 The event may be either a pre-defined PETSc event (found in 873 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 874 875 Level: advanced 876 877 .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePop() 878 @*/ 879 PetscErrorCode PetscLogEventDeactivatePush(PetscLogEvent event) 880 { 881 PetscStageLog stageLog; 882 int stage; 883 884 PetscFunctionBegin; 885 PetscCall(PetscLogGetStageLog(&stageLog)); 886 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 887 PetscCall(PetscEventPerfLogDeactivatePush(stageLog->stageInfo[stage].eventLog, event)); 888 PetscFunctionReturn(0); 889 } 890 891 /*@ 892 PetscLogEventDeactivatePop - Indicates that a particular event shouldbe logged. 893 894 Not Collective 895 896 Input Parameter: 897 . event - The event id 898 899 Usage: 900 .vb 901 PetscLogEventDeactivatePush(VEC_SetValues); 902 [code where you do not want to log VecSetValues()] 903 PetscLogEventDeactivatePop(VEC_SetValues); 904 [code where you do want to log VecSetValues()] 905 .ve 906 907 Note: 908 The event may be either a pre-defined PETSc event (found in 909 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 910 911 Level: advanced 912 913 .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePush() 914 @*/ 915 PetscErrorCode PetscLogEventDeactivatePop(PetscLogEvent event) 916 { 917 PetscStageLog stageLog; 918 int stage; 919 920 PetscFunctionBegin; 921 PetscCall(PetscLogGetStageLog(&stageLog)); 922 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 923 PetscCall(PetscEventPerfLogDeactivatePop(stageLog->stageInfo[stage].eventLog, event)); 924 PetscFunctionReturn(0); 925 } 926 927 /*@ 928 PetscLogEventSetActiveAll - Sets the event activity in every stage. 929 930 Not Collective 931 932 Input Parameters: 933 + event - The event id 934 - isActive - The activity flag determining whether the event is logged 935 936 Level: advanced 937 938 .seealso: PlogEventActivate(),PlogEventDeactivate() 939 @*/ 940 PetscErrorCode PetscLogEventSetActiveAll(PetscLogEvent event, PetscBool isActive) 941 { 942 PetscStageLog stageLog; 943 int stage; 944 945 PetscFunctionBegin; 946 PetscCall(PetscLogGetStageLog(&stageLog)); 947 for (stage = 0; stage < stageLog->numStages; stage++) { 948 if (isActive) { 949 PetscCall(PetscEventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event)); 950 } else { 951 PetscCall(PetscEventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event)); 952 } 953 } 954 PetscFunctionReturn(0); 955 } 956 957 /*@ 958 PetscLogEventActivateClass - Activates event logging for a PETSc object class. 959 960 Not Collective 961 962 Input Parameter: 963 . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc. 964 965 Level: developer 966 967 .seealso: PetscLogEventDeactivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate() 968 @*/ 969 PetscErrorCode PetscLogEventActivateClass(PetscClassId classid) 970 { 971 PetscStageLog stageLog; 972 int stage; 973 974 PetscFunctionBegin; 975 PetscCall(PetscLogGetStageLog(&stageLog)); 976 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 977 PetscCall(PetscEventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid)); 978 PetscFunctionReturn(0); 979 } 980 981 /*@ 982 PetscLogEventDeactivateClass - Deactivates event logging for a PETSc object class. 983 984 Not Collective 985 986 Input Parameter: 987 . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc. 988 989 Level: developer 990 991 .seealso: PetscLogEventActivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate() 992 @*/ 993 PetscErrorCode PetscLogEventDeactivateClass(PetscClassId classid) 994 { 995 PetscStageLog stageLog; 996 int stage; 997 998 PetscFunctionBegin; 999 PetscCall(PetscLogGetStageLog(&stageLog)); 1000 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 1001 PetscCall(PetscEventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid)); 1002 PetscFunctionReturn(0); 1003 } 1004 1005 /*MC 1006 PetscLogEventSync - Synchronizes the beginning of a user event. 1007 1008 Synopsis: 1009 #include <petsclog.h> 1010 PetscErrorCode PetscLogEventSync(int e,MPI_Comm comm) 1011 1012 Collective 1013 1014 Input Parameters: 1015 + e - integer associated with the event obtained from PetscLogEventRegister() 1016 - comm - an MPI communicator 1017 1018 Usage: 1019 .vb 1020 PetscLogEvent USER_EVENT; 1021 PetscLogEventRegister("User event",0,&USER_EVENT); 1022 PetscLogEventSync(USER_EVENT,PETSC_COMM_WORLD); 1023 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1024 [code segment to monitor] 1025 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1026 .ve 1027 1028 Notes: 1029 This routine should be called only if there is not a 1030 PetscObject available to pass to PetscLogEventBegin(). 1031 1032 Level: developer 1033 1034 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd() 1035 1036 M*/ 1037 1038 /*MC 1039 PetscLogEventBegin - Logs the beginning of a user event. 1040 1041 Synopsis: 1042 #include <petsclog.h> 1043 PetscErrorCode PetscLogEventBegin(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 1044 1045 Not Collective 1046 1047 Input Parameters: 1048 + e - integer associated with the event obtained from PetscLogEventRegister() 1049 - o1,o2,o3,o4 - objects associated with the event, or 0 1050 1051 Fortran Synopsis: 1052 void PetscLogEventBegin(int e,PetscErrorCode ierr) 1053 1054 Usage: 1055 .vb 1056 PetscLogEvent USER_EVENT; 1057 PetscLogDouble user_event_flops; 1058 PetscLogEventRegister("User event",0,&USER_EVENT); 1059 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1060 [code segment to monitor] 1061 PetscLogFlops(user_event_flops); 1062 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1063 .ve 1064 1065 Notes: 1066 You need to register each integer event with the command 1067 PetscLogEventRegister(). 1068 1069 Developer Notes: 1070 PetscLogEventBegin() and PetscLogEventBegin() return error codes instead of explicitly handling the 1071 errors that occur in the macro directly because other packages that use this macros have used them in their 1072 own functions or methods that do not return error codes and it would be disruptive to change the current 1073 behavior. 1074 1075 Level: intermediate 1076 1077 .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops() 1078 1079 M*/ 1080 1081 /*MC 1082 PetscLogEventEnd - Log the end of a user event. 1083 1084 Synopsis: 1085 #include <petsclog.h> 1086 PetscErrorCode PetscLogEventEnd(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 1087 1088 Not Collective 1089 1090 Input Parameters: 1091 + e - integer associated with the event obtained with PetscLogEventRegister() 1092 - o1,o2,o3,o4 - objects associated with the event, or 0 1093 1094 Fortran Synopsis: 1095 void PetscLogEventEnd(int e,PetscErrorCode ierr) 1096 1097 Usage: 1098 .vb 1099 PetscLogEvent USER_EVENT; 1100 PetscLogDouble user_event_flops; 1101 PetscLogEventRegister("User event",0,&USER_EVENT,); 1102 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1103 [code segment to monitor] 1104 PetscLogFlops(user_event_flops); 1105 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1106 .ve 1107 1108 Notes: 1109 You should also register each additional integer event with the command 1110 PetscLogEventRegister(). 1111 1112 Level: intermediate 1113 1114 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogFlops() 1115 1116 M*/ 1117 1118 /*@C 1119 PetscLogEventGetId - Returns the event id when given the event name. 1120 1121 Not Collective 1122 1123 Input Parameter: 1124 . name - The event name 1125 1126 Output Parameter: 1127 . event - The event, or -1 if no event with that name exists 1128 1129 Level: intermediate 1130 1131 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogStageGetId() 1132 @*/ 1133 PetscErrorCode PetscLogEventGetId(const char name[], PetscLogEvent *event) 1134 { 1135 PetscStageLog stageLog; 1136 1137 PetscFunctionBegin; 1138 PetscCall(PetscLogGetStageLog(&stageLog)); 1139 PetscCall(PetscEventRegLogGetEvent(stageLog->eventLog, name, event)); 1140 PetscFunctionReturn(0); 1141 } 1142 1143 /*------------------------------------------------ Output Functions -------------------------------------------------*/ 1144 /*@C 1145 PetscLogDump - Dumps logs of objects to a file. This file is intended to 1146 be read by bin/petscview. This program no longer exists. 1147 1148 Collective on PETSC_COMM_WORLD 1149 1150 Input Parameter: 1151 . name - an optional file name 1152 1153 Usage: 1154 .vb 1155 PetscInitialize(...); 1156 PetscLogDefaultBegin(); or PetscLogAllBegin(); 1157 ... code ... 1158 PetscLogDump(filename); 1159 PetscFinalize(); 1160 .ve 1161 1162 Notes: 1163 The default file name is 1164 $ Log.<rank> 1165 where <rank> is the processor number. If no name is specified, 1166 this file will be used. 1167 1168 Level: advanced 1169 1170 .seealso: PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogView() 1171 @*/ 1172 PetscErrorCode PetscLogDump(const char sname[]) 1173 { 1174 PetscStageLog stageLog; 1175 PetscEventPerfInfo *eventInfo; 1176 FILE *fd; 1177 char file[PETSC_MAX_PATH_LEN], fname[PETSC_MAX_PATH_LEN]; 1178 PetscLogDouble flops, _TotalTime; 1179 PetscMPIInt rank; 1180 int action, object, curStage; 1181 PetscLogEvent event; 1182 1183 PetscFunctionBegin; 1184 /* Calculate the total elapsed time */ 1185 PetscTime(&_TotalTime); 1186 _TotalTime -= petsc_BaseTime; 1187 /* Open log file */ 1188 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 1189 if (sname && sname[0]) sprintf(file, "%s.%d", sname, rank); 1190 else sprintf(file, "Log.%d", rank); 1191 PetscCall(PetscFixFilename(file, fname)); 1192 PetscCall(PetscFOpen(PETSC_COMM_WORLD, fname, "w", &fd)); 1193 PetscCheck(!(rank == 0) || !(!fd),PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN, "Cannot open file: %s", fname); 1194 /* Output totals */ 1195 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Total Flop %14e %16.8e\n", petsc_TotalFlops, _TotalTime)); 1196 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Clock Resolution %g\n", 0.0)); 1197 /* Output actions */ 1198 if (petsc_logActions) { 1199 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Actions accomplished %d\n", petsc_numActions)); 1200 for (action = 0; action < petsc_numActions; action++) { 1201 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "%g %d %d %d %d %d %d %g %g %g\n", 1202 petsc_actions[action].time, petsc_actions[action].action, (int)petsc_actions[action].event, (int)petsc_actions[action].classid, petsc_actions[action].id1, 1203 petsc_actions[action].id2, petsc_actions[action].id3, petsc_actions[action].flops, petsc_actions[action].mem, petsc_actions[action].maxmem)); 1204 } 1205 } 1206 /* Output objects */ 1207 if (petsc_logObjects) { 1208 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Objects created %d destroyed %d\n", petsc_numObjects, petsc_numObjectsDestroyed)); 1209 for (object = 0; object < petsc_numObjects; object++) { 1210 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Parent ID: %d Memory: %d\n", petsc_objects[object].parent, (int) petsc_objects[object].mem)); 1211 if (!petsc_objects[object].name[0]) { 1212 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd,"No Name\n")); 1213 } else { 1214 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Name: %s\n", petsc_objects[object].name)); 1215 } 1216 if (petsc_objects[object].info[0] != 0) { 1217 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "No Info\n")); 1218 } else { 1219 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Info: %s\n", petsc_objects[object].info)); 1220 } 1221 } 1222 } 1223 /* Output events */ 1224 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Event log:\n")); 1225 PetscCall(PetscLogGetStageLog(&stageLog)); 1226 PetscCall(PetscIntStackTop(stageLog->stack, &curStage)); 1227 eventInfo = stageLog->stageInfo[curStage].eventLog->eventInfo; 1228 for (event = 0; event < stageLog->stageInfo[curStage].eventLog->numEvents; event++) { 1229 if (eventInfo[event].time != 0.0) flops = eventInfo[event].flops/eventInfo[event].time; 1230 else flops = 0.0; 1231 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "%d %16d %16g %16g %16g\n", event, eventInfo[event].count,eventInfo[event].flops, eventInfo[event].time, flops)); 1232 } 1233 PetscCall(PetscFClose(PETSC_COMM_WORLD, fd)); 1234 PetscFunctionReturn(0); 1235 } 1236 1237 /* 1238 PetscLogView_Detailed - Each process prints the times for its own events 1239 1240 */ 1241 PetscErrorCode PetscLogView_Detailed(PetscViewer viewer) 1242 { 1243 PetscStageLog stageLog; 1244 PetscEventPerfInfo *eventInfo = NULL, *stageInfo = NULL; 1245 PetscLogDouble locTotalTime, numRed, maxMem; 1246 int numStages,numEvents,stage,event; 1247 MPI_Comm comm = PetscObjectComm((PetscObject) viewer); 1248 PetscMPIInt rank,size; 1249 1250 PetscFunctionBegin; 1251 PetscCallMPI(MPI_Comm_size(comm, &size)); 1252 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1253 /* Must preserve reduction count before we go on */ 1254 numRed = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1255 /* Get the total elapsed time */ 1256 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1257 PetscCall(PetscViewerASCIIPrintf(viewer,"size = %d\n",size)); 1258 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalTimes = {}\n")); 1259 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalMessages = {}\n")); 1260 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalMessageLens = {}\n")); 1261 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalReductions = {}\n")); 1262 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalFlop = {}\n")); 1263 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalObjects = {}\n")); 1264 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalMemory = {}\n")); 1265 PetscCall(PetscLogGetStageLog(&stageLog)); 1266 PetscCallMPI(MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm)); 1267 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages = {}\n")); 1268 for (stage=0; stage<numStages; stage++) { 1269 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"] = {}\n",stageLog->stageInfo[stage].name)); 1270 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"summary\"] = {}\n",stageLog->stageInfo[stage].name)); 1271 PetscCallMPI(MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1272 for (event = 0; event < numEvents; event++) { 1273 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"%s\"] = {}\n",stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name)); 1274 } 1275 } 1276 PetscCall(PetscMallocGetMaximumUsage(&maxMem)); 1277 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 1278 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalTimes[%d] = %g\n",rank,locTotalTime)); 1279 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessages[%d] = %g\n",rank,(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct))); 1280 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessageLens[%d] = %g\n",rank,(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len))); 1281 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalReductions[%d] = %g\n",rank,numRed)); 1282 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalFlop[%d] = %g\n",rank,petsc_TotalFlops)); 1283 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalObjects[%d] = %d\n",rank,petsc_numObjects)); 1284 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalMemory[%d] = %g\n",rank,maxMem)); 1285 PetscCall(PetscViewerFlush(viewer)); 1286 for (stage=0; stage<numStages; stage++) { 1287 stageInfo = &stageLog->stageInfo[stage].perfInfo; 1288 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"summary\"][%d] = {\"time\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g}\n", 1289 stageLog->stageInfo[stage].name,rank, 1290 stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops)); 1291 PetscCallMPI(MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1292 for (event = 0; event < numEvents; event++) { 1293 eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event]; 1294 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"%s\"][%d] = {\"count\" : %d, \"time\" : %g, \"syncTime\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g", 1295 stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name,rank, 1296 eventInfo->count,eventInfo->time,eventInfo->syncTime,eventInfo->numMessages,eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops)); 1297 if (eventInfo->dof[0] >= 0.) { 1298 PetscInt d, e; 1299 1300 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", \"dof\" : [")); 1301 for (d = 0; d < 8; ++d) { 1302 if (d > 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", ")); 1303 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->dof[d])); 1304 } 1305 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "]")); 1306 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", \"error\" : [")); 1307 for (e = 0; e < 8; ++e) { 1308 if (e > 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", ")); 1309 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->errors[e])); 1310 } 1311 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "]")); 1312 } 1313 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"}\n")); 1314 } 1315 } 1316 PetscCall(PetscViewerFlush(viewer)); 1317 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 1318 PetscFunctionReturn(0); 1319 } 1320 1321 /* 1322 PetscLogView_CSV - Each process prints the times for its own events in Comma-Separated Value Format 1323 */ 1324 PetscErrorCode PetscLogView_CSV(PetscViewer viewer) 1325 { 1326 PetscStageLog stageLog; 1327 PetscEventPerfInfo *eventInfo = NULL; 1328 PetscLogDouble locTotalTime, maxMem; 1329 int numStages,numEvents,stage,event; 1330 MPI_Comm comm = PetscObjectComm((PetscObject) viewer); 1331 PetscMPIInt rank,size; 1332 1333 PetscFunctionBegin; 1334 PetscCallMPI(MPI_Comm_size(comm, &size)); 1335 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1336 /* Must preserve reduction count before we go on */ 1337 /* Get the total elapsed time */ 1338 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1339 PetscCall(PetscLogGetStageLog(&stageLog)); 1340 PetscCallMPI(MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm)); 1341 PetscCall(PetscMallocGetMaximumUsage(&maxMem)); 1342 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 1343 PetscCall(PetscViewerASCIIPrintf(viewer,"Stage Name,Event Name,Rank,Count,Time,Num Messages,Message Length,Num Reductions,FLOP,dof0,dof1,dof2,dof3,dof4,dof5,dof6,dof7,e0,e1,e2,e3,e4,e5,e6,e7,%d\n", size)); 1344 PetscCall(PetscViewerFlush(viewer)); 1345 for (stage=0; stage<numStages; stage++) { 1346 PetscEventPerfInfo *stageInfo = &stageLog->stageInfo[stage].perfInfo; 1347 1348 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"%s,summary,%d,1,%g,%g,%g,%g,%g\n", 1349 stageLog->stageInfo[stage].name,rank,stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops)); 1350 PetscCallMPI(MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1351 for (event = 0; event < numEvents; event++) { 1352 eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event]; 1353 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"%s,%s,%d,%d,%g,%g,%g,%g,%g",stageLog->stageInfo[stage].name, 1354 stageLog->eventLog->eventInfo[event].name,rank,eventInfo->count,eventInfo->time,eventInfo->numMessages, 1355 eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops)); 1356 if (eventInfo->dof[0] >= 0.) { 1357 PetscInt d, e; 1358 1359 for (d = 0; d < 8; ++d) { 1360 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->dof[d])); 1361 } 1362 for (e = 0; e < 8; ++e) { 1363 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->errors[e])); 1364 } 1365 } 1366 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"\n")); 1367 } 1368 } 1369 PetscCall(PetscViewerFlush(viewer)); 1370 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 1371 PetscFunctionReturn(0); 1372 } 1373 1374 static PetscErrorCode PetscLogViewWarnSync(MPI_Comm comm,FILE *fd) 1375 { 1376 PetscFunctionBegin; 1377 if (!PetscLogSyncOn) PetscFunctionReturn(0); 1378 PetscCall(PetscFPrintf(comm, fd, "\n\n")); 1379 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n")); 1380 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1381 PetscCall(PetscFPrintf(comm, fd, " # WARNING!!! #\n")); 1382 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1383 PetscCall(PetscFPrintf(comm, fd, " # This program was run with logging synchronization. #\n")); 1384 PetscCall(PetscFPrintf(comm, fd, " # This option provides more meaningful imbalance #\n")); 1385 PetscCall(PetscFPrintf(comm, fd, " # figures at the expense of slowing things down and #\n")); 1386 PetscCall(PetscFPrintf(comm, fd, " # providing a distorted view of the overall runtime. #\n")); 1387 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1388 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n\n\n")); 1389 PetscFunctionReturn(0); 1390 } 1391 1392 static PetscErrorCode PetscLogViewWarnDebugging(MPI_Comm comm,FILE *fd) 1393 { 1394 PetscFunctionBegin; 1395 if (PetscDefined(USE_DEBUG)) { 1396 PetscCall(PetscFPrintf(comm, fd, "\n\n")); 1397 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n")); 1398 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1399 PetscCall(PetscFPrintf(comm, fd, " # WARNING!!! #\n")); 1400 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1401 PetscCall(PetscFPrintf(comm, fd, " # This code was compiled with a debugging option. #\n")); 1402 PetscCall(PetscFPrintf(comm, fd, " # To get timing results run ./configure #\n")); 1403 PetscCall(PetscFPrintf(comm, fd, " # using --with-debugging=no, the performance will #\n")); 1404 PetscCall(PetscFPrintf(comm, fd, " # be generally two or three times faster. #\n")); 1405 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1406 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n\n\n")); 1407 } 1408 PetscFunctionReturn(0); 1409 } 1410 1411 static PetscErrorCode PetscLogViewWarnNoGpuAwareMpi(MPI_Comm comm,FILE *fd) 1412 { 1413 #if defined(PETSC_HAVE_DEVICE) 1414 PetscMPIInt size; 1415 1416 PetscFunctionBegin; 1417 PetscCallMPI(MPI_Comm_size(comm, &size)); 1418 if (use_gpu_aware_mpi || size == 1) PetscFunctionReturn(0); 1419 PetscCall(PetscFPrintf(comm, fd, "\n\n")); 1420 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n")); 1421 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1422 PetscCall(PetscFPrintf(comm, fd, " # WARNING!!! #\n")); 1423 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1424 PetscCall(PetscFPrintf(comm, fd, " # This code was compiled with GPU support and you've #\n")); 1425 PetscCall(PetscFPrintf(comm, fd, " # created PETSc/GPU objects, but you intentionally #\n")); 1426 PetscCall(PetscFPrintf(comm, fd, " # used -use_gpu_aware_mpi 0, requiring PETSc to copy #\n")); 1427 PetscCall(PetscFPrintf(comm, fd, " # additional data between the GPU and CPU. To obtain #\n")); 1428 PetscCall(PetscFPrintf(comm, fd, " # meaningful timing results on multi-rank runs, use #\n")); 1429 PetscCall(PetscFPrintf(comm, fd, " # GPU-aware MPI instead. #\n")); 1430 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1431 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n\n\n")); 1432 PetscFunctionReturn(0); 1433 #else 1434 return 0; 1435 #endif 1436 } 1437 1438 PetscErrorCode PetscLogView_Default(PetscViewer viewer) 1439 { 1440 FILE *fd; 1441 PetscLogDouble zero = 0.0; 1442 PetscStageLog stageLog; 1443 PetscStageInfo *stageInfo = NULL; 1444 PetscEventPerfInfo *eventInfo = NULL; 1445 PetscClassPerfInfo *classInfo; 1446 char arch[128],hostname[128],username[128],pname[PETSC_MAX_PATH_LEN],date[128]; 1447 const char *name; 1448 PetscLogDouble locTotalTime, TotalTime, TotalFlops; 1449 PetscLogDouble numMessages, messageLength, avgMessLen, numReductions; 1450 PetscLogDouble stageTime, flops, flopr, mem, mess, messLen, red; 1451 PetscLogDouble fracTime, fracFlops, fracMessages, fracLength, fracReductions, fracMess, fracMessLen, fracRed; 1452 PetscLogDouble fracStageTime, fracStageFlops, fracStageMess, fracStageMessLen, fracStageRed; 1453 PetscLogDouble min, max, tot, ratio, avg, x, y; 1454 PetscLogDouble minf, maxf, totf, ratf, mint, maxt, tott, ratt, ratC, totm, totml, totr, mal, malmax, emalmax; 1455 #if defined(PETSC_HAVE_DEVICE) 1456 PetscLogDouble cct, gct, csz, gsz, gmaxt, gflops, gflopr, fracgflops; 1457 #endif 1458 PetscMPIInt minC, maxC; 1459 PetscMPIInt size, rank; 1460 PetscBool *localStageUsed, *stageUsed; 1461 PetscBool *localStageVisible, *stageVisible; 1462 int numStages, localNumEvents, numEvents; 1463 int stage, oclass; 1464 PetscLogEvent event; 1465 PetscErrorCode ierr = 0; 1466 char version[256]; 1467 MPI_Comm comm; 1468 1469 PetscFunctionBegin; 1470 PetscCall(PetscObjectGetComm((PetscObject)viewer,&comm)); 1471 PetscCall(PetscViewerASCIIGetPointer(viewer,&fd)); 1472 PetscCallMPI(MPI_Comm_size(comm, &size)); 1473 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1474 /* Get the total elapsed time */ 1475 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1476 1477 PetscCall(PetscFPrintf(comm, fd, "**************************************** ***********************************************************************************************************************\n")); 1478 PetscCall(PetscFPrintf(comm, fd, "*** WIDEN YOUR WINDOW TO 160 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***\n")); 1479 PetscCall(PetscFPrintf(comm, fd, "****************************************************************************************************************************************************************\n")); 1480 PetscCall(PetscFPrintf(comm, fd, "\n------------------------------------------------------------------ PETSc Performance Summary: -------------------------------------------------------------------\n\n")); 1481 PetscCall(PetscLogViewWarnSync(comm,fd)); 1482 PetscCall(PetscLogViewWarnDebugging(comm,fd)); 1483 PetscCall(PetscLogViewWarnNoGpuAwareMpi(comm,fd)); 1484 PetscCall(PetscGetArchType(arch,sizeof(arch))); 1485 PetscCall(PetscGetHostName(hostname,sizeof(hostname))); 1486 PetscCall(PetscGetUserName(username,sizeof(username))); 1487 PetscCall(PetscGetProgramName(pname,sizeof(pname))); 1488 PetscCall(PetscGetDate(date,sizeof(date))); 1489 PetscCall(PetscGetVersion(version,sizeof(version))); 1490 if (size == 1) { 1491 PetscCall(PetscFPrintf(comm,fd,"%s on a %s named %s with %d processor, by %s %s\n", pname, arch, hostname, size, username, date)); 1492 } else { 1493 PetscCall(PetscFPrintf(comm,fd,"%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date)); 1494 } 1495 #if defined(PETSC_HAVE_OPENMP) 1496 PetscCall(PetscFPrintf(comm,fd,"Using %" PetscInt_FMT " OpenMP threads\n", PetscNumOMPThreads)); 1497 #endif 1498 PetscCall(PetscFPrintf(comm, fd, "Using %s\n", version)); 1499 1500 /* Must preserve reduction count before we go on */ 1501 red = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1502 1503 /* Calculate summary information */ 1504 PetscCall(PetscFPrintf(comm, fd, "\n Max Max/Min Avg Total\n")); 1505 /* Time */ 1506 PetscCallMPI(MPI_Allreduce(&locTotalTime, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1507 PetscCallMPI(MPI_Allreduce(&locTotalTime, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1508 PetscCallMPI(MPI_Allreduce(&locTotalTime, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1509 avg = tot/((PetscLogDouble) size); 1510 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1511 PetscCall(PetscFPrintf(comm, fd, "Time (sec): %5.3e %7.3f %5.3e\n", max, ratio, avg)); 1512 TotalTime = tot; 1513 /* Objects */ 1514 avg = (PetscLogDouble) petsc_numObjects; 1515 PetscCallMPI(MPI_Allreduce(&avg, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1516 PetscCallMPI(MPI_Allreduce(&avg, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1517 PetscCallMPI(MPI_Allreduce(&avg, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1518 avg = tot/((PetscLogDouble) size); 1519 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1520 PetscCall(PetscFPrintf(comm, fd, "Objects: %5.3e %7.3f %5.3e\n", max, ratio, avg)); 1521 /* Flops */ 1522 PetscCallMPI(MPI_Allreduce(&petsc_TotalFlops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1523 PetscCallMPI(MPI_Allreduce(&petsc_TotalFlops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1524 PetscCallMPI(MPI_Allreduce(&petsc_TotalFlops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1525 avg = tot/((PetscLogDouble) size); 1526 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1527 PetscCall(PetscFPrintf(comm, fd, "Flops: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1528 TotalFlops = tot; 1529 /* Flops/sec -- Must talk to Barry here */ 1530 if (locTotalTime != 0.0) flops = petsc_TotalFlops/locTotalTime; else flops = 0.0; 1531 PetscCallMPI(MPI_Allreduce(&flops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1532 PetscCallMPI(MPI_Allreduce(&flops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1533 PetscCallMPI(MPI_Allreduce(&flops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1534 avg = tot/((PetscLogDouble) size); 1535 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1536 PetscCall(PetscFPrintf(comm, fd, "Flops/sec: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1537 /* Memory */ 1538 PetscCall(PetscMallocGetMaximumUsage(&mem)); 1539 if (mem > 0.0) { 1540 PetscCallMPI(MPI_Allreduce(&mem, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1541 PetscCallMPI(MPI_Allreduce(&mem, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1542 PetscCallMPI(MPI_Allreduce(&mem, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1543 avg = tot/((PetscLogDouble) size); 1544 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1545 PetscCall(PetscFPrintf(comm, fd, "Memory (bytes): %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1546 } 1547 /* Messages */ 1548 mess = 0.5*(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct); 1549 PetscCallMPI(MPI_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1550 PetscCallMPI(MPI_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1551 PetscCallMPI(MPI_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1552 avg = tot/((PetscLogDouble) size); 1553 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1554 PetscCall(PetscFPrintf(comm, fd, "MPI Msg Count: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1555 numMessages = tot; 1556 /* Message Lengths */ 1557 mess = 0.5*(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len); 1558 PetscCallMPI(MPI_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1559 PetscCallMPI(MPI_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1560 PetscCallMPI(MPI_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1561 if (numMessages != 0) avg = tot/numMessages; else avg = 0.0; 1562 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1563 PetscCall(PetscFPrintf(comm, fd, "MPI Msg Len (bytes): %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1564 messageLength = tot; 1565 /* Reductions */ 1566 PetscCallMPI(MPI_Allreduce(&red, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1567 PetscCallMPI(MPI_Allreduce(&red, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1568 PetscCallMPI(MPI_Allreduce(&red, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1569 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1570 PetscCall(PetscFPrintf(comm, fd, "MPI Reductions: %5.3e %7.3f\n", max, ratio)); 1571 numReductions = red; /* wrong because uses count from process zero */ 1572 PetscCall(PetscFPrintf(comm, fd, "\nFlop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)\n")); 1573 PetscCall(PetscFPrintf(comm, fd, " e.g., VecAXPY() for real vectors of length N --> 2N flops\n")); 1574 PetscCall(PetscFPrintf(comm, fd, " and VecAXPY() for complex vectors of length N --> 8N flops\n")); 1575 1576 /* Get total number of stages -- 1577 Currently, a single processor can register more stages than another, but stages must all be registered in order. 1578 We can removed this requirement if necessary by having a global stage numbering and indirection on the stage ID. 1579 This seems best accomplished by assoicating a communicator with each stage. 1580 */ 1581 PetscCall(PetscLogGetStageLog(&stageLog)); 1582 PetscCallMPI(MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm)); 1583 PetscCall(PetscMalloc1(numStages, &localStageUsed)); 1584 PetscCall(PetscMalloc1(numStages, &stageUsed)); 1585 PetscCall(PetscMalloc1(numStages, &localStageVisible)); 1586 PetscCall(PetscMalloc1(numStages, &stageVisible)); 1587 if (numStages > 0) { 1588 stageInfo = stageLog->stageInfo; 1589 for (stage = 0; stage < numStages; stage++) { 1590 if (stage < stageLog->numStages) { 1591 localStageUsed[stage] = stageInfo[stage].used; 1592 localStageVisible[stage] = stageInfo[stage].perfInfo.visible; 1593 } else { 1594 localStageUsed[stage] = PETSC_FALSE; 1595 localStageVisible[stage] = PETSC_TRUE; 1596 } 1597 } 1598 PetscCallMPI(MPI_Allreduce(localStageUsed, stageUsed, numStages, MPIU_BOOL, MPI_LOR, comm)); 1599 PetscCallMPI(MPI_Allreduce(localStageVisible, stageVisible, numStages, MPIU_BOOL, MPI_LAND, comm)); 1600 for (stage = 0; stage < numStages; stage++) { 1601 if (stageUsed[stage]) { 1602 PetscCall(PetscFPrintf(comm, fd, "\nSummary of Stages: ----- Time ------ ----- Flop ------ --- Messages --- -- Message Lengths -- -- Reductions --\n")); 1603 PetscCall(PetscFPrintf(comm, fd, " Avg %%Total Avg %%Total Count %%Total Avg %%Total Count %%Total\n")); 1604 break; 1605 } 1606 } 1607 for (stage = 0; stage < numStages; stage++) { 1608 if (!stageUsed[stage]) continue; 1609 /* CANNOT use MPI_Allreduce() since it might fail the line number check */ 1610 if (localStageUsed[stage]) { 1611 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1612 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1613 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1614 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1615 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1616 name = stageInfo[stage].name; 1617 } else { 1618 PetscCallMPI(MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1619 PetscCallMPI(MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1620 PetscCallMPI(MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1621 PetscCallMPI(MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1622 PetscCallMPI(MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1623 name = ""; 1624 } 1625 mess *= 0.5; messLen *= 0.5; red /= size; 1626 if (TotalTime != 0.0) fracTime = stageTime/TotalTime; else fracTime = 0.0; 1627 if (TotalFlops != 0.0) fracFlops = flops/TotalFlops; else fracFlops = 0.0; 1628 /* Talk to Barry if (stageTime != 0.0) flops = (size*flops)/stageTime; else flops = 0.0; */ 1629 if (numMessages != 0.0) fracMessages = mess/numMessages; else fracMessages = 0.0; 1630 if (mess != 0.0) avgMessLen = messLen/mess; else avgMessLen = 0.0; 1631 if (messageLength != 0.0) fracLength = messLen/messageLength; else fracLength = 0.0; 1632 if (numReductions != 0.0) fracReductions = red/numReductions; else fracReductions = 0.0; 1633 PetscCall(PetscFPrintf(comm, fd, "%2d: %15s: %6.4e %5.1f%% %6.4e %5.1f%% %5.3e %5.1f%% %5.3e %5.1f%% %5.3e %5.1f%%\n", 1634 stage, name, stageTime/size, 100.0*fracTime, flops, 100.0*fracFlops, 1635 mess, 100.0*fracMessages, avgMessLen, 100.0*fracLength, red, 100.0*fracReductions)); 1636 } 1637 } 1638 1639 PetscCall(PetscFPrintf(comm, fd,"\n------------------------------------------------------------------------------------------------------------------------\n")); 1640 PetscCall(PetscFPrintf(comm, fd, "See the 'Profiling' chapter of the users' manual for details on interpreting output.\n")); 1641 PetscCall(PetscFPrintf(comm, fd, "Phase summary info:\n")); 1642 PetscCall(PetscFPrintf(comm, fd, " Count: number of times phase was executed\n")); 1643 PetscCall(PetscFPrintf(comm, fd, " Time and Flop: Max - maximum over all processors\n")); 1644 PetscCall(PetscFPrintf(comm, fd, " Ratio - ratio of maximum to minimum over all processors\n")); 1645 PetscCall(PetscFPrintf(comm, fd, " Mess: number of messages sent\n")); 1646 PetscCall(PetscFPrintf(comm, fd, " AvgLen: average message length (bytes)\n")); 1647 PetscCall(PetscFPrintf(comm, fd, " Reduct: number of global reductions\n")); 1648 PetscCall(PetscFPrintf(comm, fd, " Global: entire computation\n")); 1649 PetscCall(PetscFPrintf(comm, fd, " Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().\n")); 1650 PetscCall(PetscFPrintf(comm, fd, " %%T - percent time in this phase %%F - percent flop in this phase\n")); 1651 PetscCall(PetscFPrintf(comm, fd, " %%M - percent messages in this phase %%L - percent message lengths in this phase\n")); 1652 PetscCall(PetscFPrintf(comm, fd, " %%R - percent reductions in this phase\n")); 1653 PetscCall(PetscFPrintf(comm, fd, " Total Mflop/s: 10e-6 * (sum of flop over all processors)/(max time over all processors)\n")); 1654 if (PetscLogMemory) { 1655 PetscCall(PetscFPrintf(comm, fd, " Malloc Mbytes: Memory allocated and kept during event (sum over all calls to event)\n")); 1656 PetscCall(PetscFPrintf(comm, fd, " EMalloc Mbytes: extra memory allocated during event and then freed (maximum over all calls to events)\n")); 1657 PetscCall(PetscFPrintf(comm, fd, " MMalloc Mbytes: Increase in high water mark of allocated memory (sum over all calls to event)\n")); 1658 PetscCall(PetscFPrintf(comm, fd, " RMI Mbytes: Increase in resident memory (sum over all calls to event)\n")); 1659 } 1660 #if defined(PETSC_HAVE_DEVICE) 1661 PetscCall(PetscFPrintf(comm, fd, " GPU Mflop/s: 10e-6 * (sum of flop on GPU over all processors)/(max GPU time over all processors)\n")); 1662 PetscCall(PetscFPrintf(comm, fd, " CpuToGpu Count: total number of CPU to GPU copies per processor\n")); 1663 PetscCall(PetscFPrintf(comm, fd, " CpuToGpu Size (Mbytes): 10e-6 * (total size of CPU to GPU copies per processor)\n")); 1664 PetscCall(PetscFPrintf(comm, fd, " GpuToCpu Count: total number of GPU to CPU copies per processor\n")); 1665 PetscCall(PetscFPrintf(comm, fd, " GpuToCpu Size (Mbytes): 10e-6 * (total size of GPU to CPU copies per processor)\n")); 1666 PetscCall(PetscFPrintf(comm, fd, " GPU %%F: percent flops on GPU in this event\n")); 1667 #endif 1668 PetscCall(PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------\n")); 1669 1670 PetscCall(PetscLogViewWarnDebugging(comm,fd)); 1671 1672 /* Report events */ 1673 PetscCall(PetscFPrintf(comm, fd,"Event Count Time (sec) Flop --- Global --- --- Stage ---- Total")); 1674 if (PetscLogMemory) { 1675 PetscCall(PetscFPrintf(comm, fd," Malloc EMalloc MMalloc RMI")); 1676 } 1677 #if defined(PETSC_HAVE_DEVICE) 1678 PetscCall(PetscFPrintf(comm, fd," GPU - CpuToGpu - - GpuToCpu - GPU")); 1679 #endif 1680 PetscCall(PetscFPrintf(comm, fd,"\n")); 1681 PetscCall(PetscFPrintf(comm, fd," Max Ratio Max Ratio Max Ratio Mess AvgLen Reduct %%T %%F %%M %%L %%R %%T %%F %%M %%L %%R Mflop/s")); 1682 if (PetscLogMemory) { 1683 PetscCall(PetscFPrintf(comm, fd," Mbytes Mbytes Mbytes Mbytes")); 1684 } 1685 #if defined(PETSC_HAVE_DEVICE) 1686 PetscCall(PetscFPrintf(comm, fd," Mflop/s Count Size Count Size %%F")); 1687 #endif 1688 PetscCall(PetscFPrintf(comm, fd,"\n")); 1689 PetscCall(PetscFPrintf(comm, fd,"------------------------------------------------------------------------------------------------------------------------")); 1690 if (PetscLogMemory) { 1691 PetscCall(PetscFPrintf(comm, fd,"-----------------------------")); 1692 } 1693 #if defined(PETSC_HAVE_DEVICE) 1694 PetscCall(PetscFPrintf(comm, fd,"---------------------------------------")); 1695 #endif 1696 PetscCall(PetscFPrintf(comm, fd,"\n")); 1697 1698 /* Problem: The stage name will not show up unless the stage executed on proc 1 */ 1699 for (stage = 0; stage < numStages; stage++) { 1700 if (!stageVisible[stage]) continue; 1701 /* CANNOT use MPI_Allreduce() since it might fail the line number check */ 1702 if (localStageUsed[stage]) { 1703 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name)); 1704 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1705 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1706 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1707 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1708 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1709 } else { 1710 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage)); 1711 PetscCallMPI(MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1712 PetscCallMPI(MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1713 PetscCallMPI(MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1714 PetscCallMPI(MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1715 PetscCallMPI(MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1716 } 1717 mess *= 0.5; messLen *= 0.5; red /= size; 1718 1719 /* Get total number of events in this stage -- 1720 Currently, a single processor can register more events than another, but events must all be registered in order, 1721 just like stages. We can removed this requirement if necessary by having a global event numbering and indirection 1722 on the event ID. This seems best accomplished by associating a communicator with each stage. 1723 1724 Problem: If the event did not happen on proc 1, its name will not be available. 1725 Problem: Event visibility is not implemented 1726 */ 1727 if (localStageUsed[stage]) { 1728 eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo; 1729 localNumEvents = stageLog->stageInfo[stage].eventLog->numEvents; 1730 } else localNumEvents = 0; 1731 PetscCallMPI(MPI_Allreduce(&localNumEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1732 for (event = 0; event < numEvents; event++) { 1733 /* CANNOT use MPI_Allreduce() since it might fail the line number check */ 1734 if (localStageUsed[stage] && (event < stageLog->stageInfo[stage].eventLog->numEvents) && (eventInfo[event].depth == 0)) { 1735 if ((eventInfo[event].count > 0) && (eventInfo[event].time > 0.0)) flopr = eventInfo[event].flops; else flopr = 0.0; 1736 PetscCallMPI(MPI_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1737 PetscCallMPI(MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1738 PetscCallMPI(MPI_Allreduce(&eventInfo[event].flops, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1739 PetscCallMPI(MPI_Allreduce(&eventInfo[event].time, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1740 PetscCallMPI(MPI_Allreduce(&eventInfo[event].time, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1741 PetscCallMPI(MPI_Allreduce(&eventInfo[event].time, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1742 PetscCallMPI(MPI_Allreduce(&eventInfo[event].numMessages, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1743 PetscCallMPI(MPI_Allreduce(&eventInfo[event].messageLength, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1744 PetscCallMPI(MPI_Allreduce(&eventInfo[event].numReductions, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1745 PetscCallMPI(MPI_Allreduce(&eventInfo[event].count, &minC, 1, MPI_INT, MPI_MIN, comm)); 1746 PetscCallMPI(MPI_Allreduce(&eventInfo[event].count, &maxC, 1, MPI_INT, MPI_MAX, comm)); 1747 if (PetscLogMemory) { 1748 PetscCallMPI(MPI_Allreduce(&eventInfo[event].memIncrease, &mem, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1749 PetscCallMPI(MPI_Allreduce(&eventInfo[event].mallocSpace, &mal, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1750 PetscCallMPI(MPI_Allreduce(&eventInfo[event].mallocIncrease, &malmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1751 PetscCallMPI(MPI_Allreduce(&eventInfo[event].mallocIncreaseEvent, &emalmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1752 } 1753 #if defined(PETSC_HAVE_DEVICE) 1754 PetscCallMPI(MPI_Allreduce(&eventInfo[event].CpuToGpuCount, &cct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1755 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuToCpuCount, &gct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1756 PetscCallMPI(MPI_Allreduce(&eventInfo[event].CpuToGpuSize, &csz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1757 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuToCpuSize, &gsz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1758 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuFlops, &gflops,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1759 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuTime, &gmaxt ,1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1760 #endif 1761 name = stageLog->eventLog->eventInfo[event].name; 1762 } else { 1763 flopr = 0.0; 1764 PetscCallMPI(MPI_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1765 PetscCallMPI(MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1766 PetscCallMPI(MPI_Allreduce(&zero, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1767 PetscCallMPI(MPI_Allreduce(&zero, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1768 PetscCallMPI(MPI_Allreduce(&zero, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1769 PetscCallMPI(MPI_Allreduce(&zero, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1770 PetscCallMPI(MPI_Allreduce(&zero, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1771 PetscCallMPI(MPI_Allreduce(&zero, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1772 PetscCallMPI(MPI_Allreduce(&zero, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1773 PetscCallMPI(MPI_Allreduce(&ierr, &minC, 1, MPI_INT, MPI_MIN, comm)); 1774 PetscCallMPI(MPI_Allreduce(&ierr, &maxC, 1, MPI_INT, MPI_MAX, comm)); 1775 if (PetscLogMemory) { 1776 PetscCallMPI(MPI_Allreduce(&zero, &mem, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1777 PetscCallMPI(MPI_Allreduce(&zero, &mal, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1778 PetscCallMPI(MPI_Allreduce(&zero, &malmax, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1779 PetscCallMPI(MPI_Allreduce(&zero, &emalmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1780 } 1781 #if defined(PETSC_HAVE_DEVICE) 1782 PetscCallMPI(MPI_Allreduce(&zero, &cct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1783 PetscCallMPI(MPI_Allreduce(&zero, &gct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1784 PetscCallMPI(MPI_Allreduce(&zero, &csz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1785 PetscCallMPI(MPI_Allreduce(&zero, &gsz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1786 PetscCallMPI(MPI_Allreduce(&zero, &gflops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1787 PetscCallMPI(MPI_Allreduce(&zero, &gmaxt , 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1788 #endif 1789 name = ""; 1790 } 1791 if (mint < 0.0) { 1792 PetscCall(PetscFPrintf(comm, fd, "WARNING!!! Minimum time %g over all processors for %s is negative! This happens\n on some machines whose times cannot handle too rapid calls.!\n artificially changing minimum to zero.\n",mint,name)); 1793 mint = 0; 1794 } 1795 PetscCheck(minf >= 0.0,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Minimum flop %g over all processors for %s is negative! Not possible!",minf,name); 1796 totm *= 0.5; totml *= 0.5; totr /= size; 1797 1798 if (maxC != 0) { 1799 if (minC != 0) ratC = ((PetscLogDouble)maxC)/minC;else ratC = 0.0; 1800 if (mint != 0.0) ratt = maxt/mint; else ratt = 0.0; 1801 if (minf != 0.0) ratf = maxf/minf; else ratf = 0.0; 1802 if (TotalTime != 0.0) fracTime = tott/TotalTime; else fracTime = 0.0; 1803 if (TotalFlops != 0.0) fracFlops = totf/TotalFlops; else fracFlops = 0.0; 1804 if (stageTime != 0.0) fracStageTime = tott/stageTime; else fracStageTime = 0.0; 1805 if (flops != 0.0) fracStageFlops = totf/flops; else fracStageFlops = 0.0; 1806 if (numMessages != 0.0) fracMess = totm/numMessages; else fracMess = 0.0; 1807 if (messageLength != 0.0) fracMessLen = totml/messageLength; else fracMessLen = 0.0; 1808 if (numReductions != 0.0) fracRed = totr/numReductions; else fracRed = 0.0; 1809 if (mess != 0.0) fracStageMess = totm/mess; else fracStageMess = 0.0; 1810 if (messLen != 0.0) fracStageMessLen = totml/messLen; else fracStageMessLen = 0.0; 1811 if (red != 0.0) fracStageRed = totr/red; else fracStageRed = 0.0; 1812 if (totm != 0.0) totml /= totm; else totml = 0.0; 1813 if (maxt != 0.0) flopr = totf/maxt; else flopr = 0.0; 1814 if (fracStageTime > 1.00) PetscCall(PetscFPrintf(comm, fd,"Warning -- total time of event greater than time of entire stage -- something is wrong with the timer\n")); 1815 PetscCall(PetscFPrintf(comm, fd, 1816 "%-16s %7d%4.1f %5.4e%4.1f %3.2e%4.1f %2.1e %2.1e %2.1e%3.0f%3.0f%3.0f%3.0f%3.0f %3.0f%3.0f%3.0f%3.0f%3.0f %5.0f", 1817 name, maxC, ratC, maxt, ratt, maxf, ratf, totm, totml, totr, 1818 100.0*fracTime, 100.0*fracFlops, 100.0*fracMess, 100.0*fracMessLen, 100.0*fracRed, 1819 100.0*fracStageTime, 100.0*fracStageFlops, 100.0*fracStageMess, 100.0*fracStageMessLen, 100.0*fracStageRed, 1820 PetscAbs(flopr)/1.0e6)); 1821 if (PetscLogMemory) { 1822 PetscCall(PetscFPrintf(comm, fd," %5.0f %5.0f %5.0f %5.0f",mal/1.0e6,emalmax/1.0e6,malmax/1.0e6,mem/1.0e6)); 1823 } 1824 #if defined(PETSC_HAVE_DEVICE) 1825 if (totf != 0.0) fracgflops = gflops/totf; else fracgflops = 0.0; 1826 if (gmaxt != 0.0) gflopr = gflops/gmaxt; else gflopr = 0.0; 1827 PetscCall(PetscFPrintf(comm, fd," %5.0f %4.0f %3.2e %4.0f %3.2e% 3.0f",PetscAbs(gflopr)/1.0e6,cct/size,csz/(1.0e6*size),gct/size,gsz/(1.0e6*size),100.0*fracgflops)); 1828 #endif 1829 PetscCall(PetscFPrintf(comm, fd,"\n")); 1830 } 1831 } 1832 } 1833 1834 /* Memory usage and object creation */ 1835 PetscCall(PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------")); 1836 if (PetscLogMemory) { 1837 PetscCall(PetscFPrintf(comm, fd, "-----------------------------")); 1838 } 1839 #if defined(PETSC_HAVE_DEVICE) 1840 PetscCall(PetscFPrintf(comm, fd, "---------------------------------------")); 1841 #endif 1842 PetscCall(PetscFPrintf(comm, fd, "\n")); 1843 PetscCall(PetscFPrintf(comm, fd, "\n")); 1844 PetscCall(PetscFPrintf(comm, fd, "Memory usage is given in bytes:\n\n")); 1845 1846 /* Right now, only stages on the first processor are reported here, meaning only objects associated with 1847 the global communicator, or MPI_COMM_SELF for proc 1. We really should report global stats and then 1848 stats for stages local to processor sets. 1849 */ 1850 /* We should figure out the longest object name here (now 20 characters) */ 1851 PetscCall(PetscFPrintf(comm, fd, "Object Type Creations Destructions Memory Descendants' Mem.\n")); 1852 PetscCall(PetscFPrintf(comm, fd, "Reports information only for process 0.\n")); 1853 for (stage = 0; stage < numStages; stage++) { 1854 if (localStageUsed[stage]) { 1855 classInfo = stageLog->stageInfo[stage].classLog->classInfo; 1856 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name)); 1857 for (oclass = 0; oclass < stageLog->stageInfo[stage].classLog->numClasses; oclass++) { 1858 if ((classInfo[oclass].creations > 0) || (classInfo[oclass].destructions > 0)) { 1859 PetscCall(PetscFPrintf(comm, fd, "%20s %5d %5d %11.0f %g\n", stageLog->classLog->classInfo[oclass].name, 1860 classInfo[oclass].creations, classInfo[oclass].destructions, classInfo[oclass].mem, 1861 classInfo[oclass].descMem)); 1862 } 1863 } 1864 } else { 1865 if (!localStageVisible[stage]) continue; 1866 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage)); 1867 } 1868 } 1869 1870 PetscCall(PetscFree(localStageUsed)); 1871 PetscCall(PetscFree(stageUsed)); 1872 PetscCall(PetscFree(localStageVisible)); 1873 PetscCall(PetscFree(stageVisible)); 1874 1875 /* Information unrelated to this particular run */ 1876 PetscCall(PetscFPrintf(comm, fd, "========================================================================================================================\n")); 1877 PetscTime(&y); 1878 PetscTime(&x); 1879 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 1880 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 1881 PetscCall(PetscFPrintf(comm,fd,"Average time to get PetscTime(): %g\n", (y-x)/10.0)); 1882 /* MPI information */ 1883 if (size > 1) { 1884 MPI_Status status; 1885 PetscMPIInt tag; 1886 MPI_Comm newcomm; 1887 1888 PetscCallMPI(MPI_Barrier(comm)); 1889 PetscTime(&x); 1890 PetscCallMPI(MPI_Barrier(comm)); 1891 PetscCallMPI(MPI_Barrier(comm)); 1892 PetscCallMPI(MPI_Barrier(comm)); 1893 PetscCallMPI(MPI_Barrier(comm)); 1894 PetscCallMPI(MPI_Barrier(comm)); 1895 PetscTime(&y); 1896 PetscCall(PetscFPrintf(comm, fd, "Average time for MPI_Barrier(): %g\n", (y-x)/5.0)); 1897 PetscCall(PetscCommDuplicate(comm,&newcomm, &tag)); 1898 PetscCallMPI(MPI_Barrier(comm)); 1899 if (rank) { 1900 PetscCallMPI(MPI_Recv(NULL, 0, MPI_INT, rank-1, tag, newcomm, &status)); 1901 PetscCallMPI(MPI_Send(NULL, 0, MPI_INT, (rank+1)%size, tag, newcomm)); 1902 } else { 1903 PetscTime(&x); 1904 PetscCallMPI(MPI_Send(NULL, 0, MPI_INT, 1, tag, newcomm)); 1905 PetscCallMPI(MPI_Recv(NULL, 0, MPI_INT, size-1, tag, newcomm, &status)); 1906 PetscTime(&y); 1907 PetscCall(PetscFPrintf(comm,fd,"Average time for zero size MPI_Send(): %g\n", (y-x)/size)); 1908 } 1909 PetscCall(PetscCommDestroy(&newcomm)); 1910 } 1911 PetscCall(PetscOptionsView(NULL,viewer)); 1912 1913 /* Machine and compile information */ 1914 #if defined(PETSC_USE_FORTRAN_KERNELS) 1915 PetscCall(PetscFPrintf(comm, fd, "Compiled with FORTRAN kernels\n")); 1916 #else 1917 PetscCall(PetscFPrintf(comm, fd, "Compiled without FORTRAN kernels\n")); 1918 #endif 1919 #if defined(PETSC_USE_64BIT_INDICES) 1920 PetscCall(PetscFPrintf(comm, fd, "Compiled with 64 bit PetscInt\n")); 1921 #elif defined(PETSC_USE___FLOAT128) 1922 PetscCall(PetscFPrintf(comm, fd, "Compiled with 32 bit PetscInt\n")); 1923 #endif 1924 #if defined(PETSC_USE_REAL_SINGLE) 1925 PetscCall(PetscFPrintf(comm, fd, "Compiled with single precision PetscScalar and PetscReal\n")); 1926 #elif defined(PETSC_USE___FLOAT128) 1927 PetscCall(PetscFPrintf(comm, fd, "Compiled with 128 bit precision PetscScalar and PetscReal\n")); 1928 #endif 1929 #if defined(PETSC_USE_REAL_MAT_SINGLE) 1930 PetscCall(PetscFPrintf(comm, fd, "Compiled with single precision matrices\n")); 1931 #else 1932 PetscCall(PetscFPrintf(comm, fd, "Compiled with full precision matrices (default)\n")); 1933 #endif 1934 PetscCall(PetscFPrintf(comm, fd, "sizeof(short) %d sizeof(int) %d sizeof(long) %d sizeof(void*) %d sizeof(PetscScalar) %d sizeof(PetscInt) %d\n", 1935 (int) sizeof(short), (int) sizeof(int), (int) sizeof(long), (int) sizeof(void*),(int) sizeof(PetscScalar),(int) sizeof(PetscInt))); 1936 1937 PetscCall(PetscFPrintf(comm, fd, "Configure options: %s",petscconfigureoptions)); 1938 PetscCall(PetscFPrintf(comm, fd, "%s", petscmachineinfo)); 1939 PetscCall(PetscFPrintf(comm, fd, "%s", petsccompilerinfo)); 1940 PetscCall(PetscFPrintf(comm, fd, "%s", petsccompilerflagsinfo)); 1941 PetscCall(PetscFPrintf(comm, fd, "%s", petsclinkerinfo)); 1942 1943 /* Cleanup */ 1944 PetscCall(PetscFPrintf(comm, fd, "\n")); 1945 PetscCall(PetscLogViewWarnNoGpuAwareMpi(comm,fd)); 1946 PetscCall(PetscLogViewWarnDebugging(comm,fd)); 1947 PetscFunctionReturn(0); 1948 } 1949 1950 /*@C 1951 PetscLogView - Prints a summary of the logging. 1952 1953 Collective over MPI_Comm 1954 1955 Input Parameter: 1956 . viewer - an ASCII viewer 1957 1958 Options Database Keys: 1959 + -log_view [:filename] - Prints summary of log information 1960 . -log_view :filename.py:ascii_info_detail - Saves logging information from each process as a Python file 1961 . -log_view :filename.xml:ascii_xml - Saves a summary of the logging information in a nested format (see below for how to view it) 1962 . -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) 1963 . -log_all - Saves a file Log.rank for each MPI process with details of each step of the computation 1964 - -log_trace [filename] - Displays a trace of what each process is doing 1965 1966 Notes: 1967 It is possible to control the logging programatically but we recommend using the options database approach whenever possible 1968 By default the summary is printed to stdout. 1969 1970 Before calling this routine you must have called either PetscLogDefaultBegin() or PetscLogNestedBegin() 1971 1972 If PETSc is configured with --with-logging=0 then this functionality is not available 1973 1974 To view the nested XML format filename.xml first copy ${PETSC_DIR}/share/petsc/xml/performance_xml2html.xsl to the current 1975 directory then open filename.xml with your browser. Specific notes for certain browsers 1976 $ Firefox and Internet explorer - simply open the file 1977 $ Google Chrome - you must start up Chrome with the option --allow-file-access-from-files 1978 $ Safari - see https://ccm.net/faq/36342-safari-how-to-enable-local-file-access 1979 or one can use the package http://xmlsoft.org/XSLT/xsltproc2.html to translate the xml file to html and then open it with 1980 your browser. 1981 Alternatively, use the script ${PETSC_DIR}/lib/petsc/bin/petsc-performance-view to automatically open a new browser 1982 window and render the XML log file contents. 1983 1984 The nested XML format was kindly donated by Koos Huijssen and Christiaan M. Klaij MARITIME RESEARCH INSTITUTE NETHERLANDS 1985 1986 The Flame Graph output can be visualised using either the original Flame Graph script (https://github.com/brendangregg/FlameGraph) 1987 or using speedscope (https://www.speedscope.app). 1988 Old XML profiles may be converted into this format using the script ${PETSC_DIR}/lib/petsc/bin/xml2flamegraph.py. 1989 1990 Level: beginner 1991 1992 .seealso: PetscLogDefaultBegin(), PetscLogDump() 1993 @*/ 1994 PetscErrorCode PetscLogView(PetscViewer viewer) 1995 { 1996 PetscBool isascii; 1997 PetscViewerFormat format; 1998 int stage, lastStage; 1999 PetscStageLog stageLog; 2000 2001 PetscFunctionBegin; 2002 PetscCheck(PetscLogPLB,PETSC_COMM_SELF,PETSC_ERR_SUP,"Must use -log_view or PetscLogDefaultBegin() before calling this routine"); 2003 /* Pop off any stages the user forgot to remove */ 2004 lastStage = 0; 2005 PetscCall(PetscLogGetStageLog(&stageLog)); 2006 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 2007 while (stage >= 0) { 2008 lastStage = stage; 2009 PetscCall(PetscStageLogPop(stageLog)); 2010 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 2011 } 2012 PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii)); 2013 PetscCheck(isascii,PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Currently can only view logging to ASCII"); 2014 PetscCall(PetscViewerGetFormat(viewer,&format)); 2015 if (format == PETSC_VIEWER_DEFAULT || format == PETSC_VIEWER_ASCII_INFO) { 2016 PetscCall(PetscLogView_Default(viewer)); 2017 } else if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) { 2018 PetscCall(PetscLogView_Detailed(viewer)); 2019 } else if (format == PETSC_VIEWER_ASCII_CSV) { 2020 PetscCall(PetscLogView_CSV(viewer)); 2021 } else if (format == PETSC_VIEWER_ASCII_XML) { 2022 PetscCall(PetscLogView_Nested(viewer)); 2023 } else if (format == PETSC_VIEWER_ASCII_FLAMEGRAPH) { 2024 PetscCall(PetscLogView_Flamegraph(viewer)); 2025 } 2026 PetscCall(PetscStageLogPush(stageLog, lastStage)); 2027 PetscFunctionReturn(0); 2028 } 2029 2030 /*@C 2031 PetscLogViewFromOptions - Processes command line options to determine if/how a PetscLog is to be viewed. 2032 2033 Collective on PETSC_COMM_WORLD 2034 2035 Not normally called by user 2036 2037 Level: intermediate 2038 2039 @*/ 2040 PetscErrorCode PetscLogViewFromOptions(void) 2041 { 2042 PetscViewer viewer; 2043 PetscBool flg; 2044 PetscViewerFormat format; 2045 2046 PetscFunctionBegin; 2047 PetscCall(PetscOptionsGetViewer(PETSC_COMM_WORLD,NULL,NULL,"-log_view",&viewer,&format,&flg)); 2048 if (flg) { 2049 PetscCall(PetscViewerPushFormat(viewer,format)); 2050 PetscCall(PetscLogView(viewer)); 2051 PetscCall(PetscViewerPopFormat(viewer)); 2052 PetscCall(PetscViewerDestroy(&viewer)); 2053 } 2054 PetscFunctionReturn(0); 2055 } 2056 2057 /*----------------------------------------------- Counter Functions -------------------------------------------------*/ 2058 /*@C 2059 PetscGetFlops - Returns the number of flops used on this processor 2060 since the program began. 2061 2062 Not Collective 2063 2064 Output Parameter: 2065 flops - number of floating point operations 2066 2067 Notes: 2068 A global counter logs all PETSc flop counts. The user can use 2069 PetscLogFlops() to increment this counter to include flops for the 2070 application code. 2071 2072 Level: intermediate 2073 2074 .seealso: PetscTime(), PetscLogFlops() 2075 @*/ 2076 PetscErrorCode PetscGetFlops(PetscLogDouble *flops) 2077 { 2078 PetscFunctionBegin; 2079 *flops = petsc_TotalFlops; 2080 PetscFunctionReturn(0); 2081 } 2082 2083 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...) 2084 { 2085 size_t fullLength; 2086 va_list Argp; 2087 2088 PetscFunctionBegin; 2089 if (!petsc_logObjects) PetscFunctionReturn(0); 2090 va_start(Argp, format); 2091 PetscCall(PetscVSNPrintf(petsc_objects[obj->id].info, 64,format,&fullLength, Argp)); 2092 va_end(Argp); 2093 PetscFunctionReturn(0); 2094 } 2095 2096 /*MC 2097 PetscLogFlops - Adds floating point operations to the global counter. 2098 2099 Synopsis: 2100 #include <petsclog.h> 2101 PetscErrorCode PetscLogFlops(PetscLogDouble f) 2102 2103 Not Collective 2104 2105 Input Parameter: 2106 . f - flop counter 2107 2108 Usage: 2109 .vb 2110 PetscLogEvent USER_EVENT; 2111 PetscLogEventRegister("User event",0,&USER_EVENT); 2112 PetscLogEventBegin(USER_EVENT,0,0,0,0); 2113 [code segment to monitor] 2114 PetscLogFlops(user_flops) 2115 PetscLogEventEnd(USER_EVENT,0,0,0,0); 2116 .ve 2117 2118 Notes: 2119 A global counter logs all PETSc flop counts. The user can use 2120 PetscLogFlops() to increment this counter to include flops for the 2121 application code. 2122 2123 Level: intermediate 2124 2125 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscGetFlops() 2126 2127 M*/ 2128 2129 /*MC 2130 PetscPreLoadBegin - Begin a segment of code that may be preloaded (run twice) 2131 to get accurate timings 2132 2133 Synopsis: 2134 #include <petsclog.h> 2135 void PetscPreLoadBegin(PetscBool flag,char *name); 2136 2137 Not Collective 2138 2139 Input Parameters: 2140 + flag - PETSC_TRUE to run twice, PETSC_FALSE to run once, may be overridden 2141 with command line option -preload true or -preload false 2142 - name - name of first stage (lines of code timed separately with -log_view) to 2143 be preloaded 2144 2145 Usage: 2146 .vb 2147 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2148 lines of code 2149 PetscPreLoadStage("second stage"); 2150 lines of code 2151 PetscPreLoadEnd(); 2152 .ve 2153 2154 Notes: 2155 Only works in C/C++, not Fortran 2156 2157 Flags available within the macro. 2158 + PetscPreLoadingUsed - true if we are or have done preloading 2159 . PetscPreLoadingOn - true if it is CURRENTLY doing preload 2160 . PetscPreLoadIt - 0 for the first computation (with preloading turned off it is only 0) 1 for the second 2161 - PetscPreLoadMax - number of times it will do the computation, only one when preloading is turned on 2162 The first two variables are available throughout the program, the second two only between the PetscPreLoadBegin() 2163 and PetscPreLoadEnd() 2164 2165 Level: intermediate 2166 2167 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadEnd(), PetscPreLoadStage() 2168 2169 M*/ 2170 2171 /*MC 2172 PetscPreLoadEnd - End a segment of code that may be preloaded (run twice) 2173 to get accurate timings 2174 2175 Synopsis: 2176 #include <petsclog.h> 2177 void PetscPreLoadEnd(void); 2178 2179 Not Collective 2180 2181 Usage: 2182 .vb 2183 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2184 lines of code 2185 PetscPreLoadStage("second stage"); 2186 lines of code 2187 PetscPreLoadEnd(); 2188 .ve 2189 2190 Notes: 2191 only works in C/C++ not fortran 2192 2193 Level: intermediate 2194 2195 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadStage() 2196 2197 M*/ 2198 2199 /*MC 2200 PetscPreLoadStage - Start a new segment of code to be timed separately. 2201 to get accurate timings 2202 2203 Synopsis: 2204 #include <petsclog.h> 2205 void PetscPreLoadStage(char *name); 2206 2207 Not Collective 2208 2209 Usage: 2210 .vb 2211 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2212 lines of code 2213 PetscPreLoadStage("second stage"); 2214 lines of code 2215 PetscPreLoadEnd(); 2216 .ve 2217 2218 Notes: 2219 only works in C/C++ not fortran 2220 2221 Level: intermediate 2222 2223 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd() 2224 2225 M*/ 2226 2227 #if PetscDefined(HAVE_DEVICE) 2228 #include <petsc/private/deviceimpl.h> 2229 2230 /*-------------------------------------------- GPU event Functions ----------------------------------------------*/ 2231 /*@C 2232 PetscLogGpuTimeBegin - Start timer for device 2233 2234 Notes: 2235 When CUDA or HIP is enabled, the timer is run on the GPU, it is a separate logging of time devoted to GPU computations (excluding kernel launch times). 2236 When CUDA or HIP is not available, the timer is run on the CPU, it is a separate logging of time devoted to GPU computations (including kernel launch times). 2237 There is no need to call WaitForCUDA() or WaitForHIP() between PetscLogGpuTimeBegin and PetscLogGpuTimeEnd 2238 This timer should NOT include times for data transfers between the GPU and CPU, nor setup actions such as allocating space. 2239 The regular logging captures the time for data transfers and any CPU activites during the event 2240 It is used to compute the flop rate on the GPU as it is actively engaged in running a kernel. 2241 2242 Developer Notes: 2243 The GPU event timer captures the execution time of all the kernels launched in the default stream by the CPU between PetscLogGpuTimeBegin() and PetsLogGpuTimeEnd(). 2244 PetscLogGpuTimeBegin() and PetsLogGpuTimeEnd() insert the begin and end events into the default stream (stream 0). The device will record a time stamp for the event when it reaches that event in the stream. The function xxxEventSynchronize() is called in PetsLogGpuTimeEnd() to block CPU execution, but not continued GPU excution, until the timer event is recorded. 2245 2246 Level: intermediate 2247 2248 .seealso: PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeEnd() 2249 @*/ 2250 PetscErrorCode PetscLogGpuTimeBegin(void) 2251 { 2252 PetscFunctionBegin; 2253 if (!PetscLogPLB) PetscFunctionReturn(0); 2254 if (PetscDefined(HAVE_CUDA) || PetscDefined(HAVE_HIP)) { 2255 PetscDeviceContext dctx; 2256 2257 PetscCall(PetscDeviceContextGetCurrentContext(&dctx)); 2258 PetscCall(PetscDeviceContextBeginTimer_Internal(dctx)); 2259 } else { 2260 PetscCall(PetscTimeSubtract(&petsc_gtime)); 2261 } 2262 PetscFunctionReturn(0); 2263 } 2264 2265 /*@C 2266 PetscLogGpuTimeEnd - Stop timer for device 2267 2268 Level: intermediate 2269 2270 .seealso: PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeBegin() 2271 @*/ 2272 PetscErrorCode PetscLogGpuTimeEnd(void) 2273 { 2274 PetscFunctionBegin; 2275 if (!PetscLogPLE) PetscFunctionReturn(0); 2276 if (PetscDefined(HAVE_CUDA) || PetscDefined(HAVE_HIP)) { 2277 PetscDeviceContext dctx; 2278 PetscLogDouble elapsed; 2279 2280 PetscCall(PetscDeviceContextGetCurrentContext(&dctx)); 2281 PetscCall(PetscDeviceContextEndTimer_Internal(dctx,&elapsed)); 2282 petsc_gtime += (elapsed/1000.0); 2283 } else { 2284 PetscCall(PetscTimeAdd(&petsc_gtime)); 2285 } 2286 PetscFunctionReturn(0); 2287 } 2288 #endif /* end of PETSC_HAVE_DEVICE */ 2289 2290 #else /* end of -DPETSC_USE_LOG section */ 2291 2292 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...) 2293 { 2294 PetscFunctionBegin; 2295 PetscFunctionReturn(0); 2296 } 2297 2298 #endif /* PETSC_USE_LOG*/ 2299 2300 PetscClassId PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID; 2301 PetscClassId PETSC_OBJECT_CLASSID = 0; 2302 2303 /*@C 2304 PetscClassIdRegister - Registers a new class name for objects and logging operations in an application code. 2305 2306 Not Collective 2307 2308 Input Parameter: 2309 . name - The class name 2310 2311 Output Parameter: 2312 . oclass - The class id or classid 2313 2314 Level: developer 2315 2316 @*/ 2317 PetscErrorCode PetscClassIdRegister(const char name[],PetscClassId *oclass) 2318 { 2319 #if defined(PETSC_USE_LOG) 2320 PetscStageLog stageLog; 2321 PetscInt stage; 2322 #endif 2323 2324 PetscFunctionBegin; 2325 *oclass = ++PETSC_LARGEST_CLASSID; 2326 #if defined(PETSC_USE_LOG) 2327 PetscCall(PetscLogGetStageLog(&stageLog)); 2328 PetscCall(PetscClassRegLogRegister(stageLog->classLog, name, *oclass)); 2329 for (stage = 0; stage < stageLog->numStages; stage++) { 2330 PetscCall(PetscClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses)); 2331 } 2332 #endif 2333 PetscFunctionReturn(0); 2334 } 2335 2336 #if defined(PETSC_USE_LOG) && defined(PETSC_HAVE_MPE) 2337 #include <mpe.h> 2338 2339 PetscBool PetscBeganMPE = PETSC_FALSE; 2340 2341 PETSC_INTERN PetscErrorCode PetscLogEventBeginMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject); 2342 PETSC_INTERN PetscErrorCode PetscLogEventEndMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject); 2343 2344 /*@C 2345 PetscLogMPEBegin - Turns on MPE logging of events. This creates large log files 2346 and slows the program down. 2347 2348 Collective over PETSC_COMM_WORLD 2349 2350 Options Database Keys: 2351 . -log_mpe - Prints extensive log information 2352 2353 Notes: 2354 A related routine is PetscLogDefaultBegin() (with the options key -log_view), which is 2355 intended for production runs since it logs only flop rates and object 2356 creation (and should not significantly slow the programs). 2357 2358 Level: advanced 2359 2360 .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogEventActivate(), 2361 PetscLogEventDeactivate() 2362 @*/ 2363 PetscErrorCode PetscLogMPEBegin(void) 2364 { 2365 PetscFunctionBegin; 2366 /* Do MPE initialization */ 2367 if (!MPE_Initialized_logging()) { /* This function exists in mpich 1.1.2 and higher */ 2368 PetscCall(PetscInfo(0,"Initializing MPE.\n")); 2369 PetscCall(MPE_Init_log()); 2370 2371 PetscBeganMPE = PETSC_TRUE; 2372 } else { 2373 PetscCall(PetscInfo(0,"MPE already initialized. Not attempting to reinitialize.\n")); 2374 } 2375 PetscCall(PetscLogSet(PetscLogEventBeginMPE, PetscLogEventEndMPE)); 2376 PetscFunctionReturn(0); 2377 } 2378 2379 /*@C 2380 PetscLogMPEDump - Dumps the MPE logging info to file for later use with Jumpshot. 2381 2382 Collective over PETSC_COMM_WORLD 2383 2384 Level: advanced 2385 2386 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogMPEBegin() 2387 @*/ 2388 PetscErrorCode PetscLogMPEDump(const char sname[]) 2389 { 2390 char name[PETSC_MAX_PATH_LEN]; 2391 2392 PetscFunctionBegin; 2393 if (PetscBeganMPE) { 2394 PetscCall(PetscInfo(0,"Finalizing MPE.\n")); 2395 if (sname) { 2396 PetscCall(PetscStrcpy(name,sname)); 2397 } else { 2398 PetscCall(PetscGetProgramName(name,sizeof(name))); 2399 } 2400 PetscCall(MPE_Finish_log(name)); 2401 } else { 2402 PetscCall(PetscInfo(0,"Not finalizing MPE (not started by PETSc).\n")); 2403 } 2404 PetscFunctionReturn(0); 2405 } 2406 2407 #define PETSC_RGB_COLORS_MAX 39 2408 static const char *PetscLogMPERGBColors[PETSC_RGB_COLORS_MAX] = { 2409 "OliveDrab: ", 2410 "BlueViolet: ", 2411 "CadetBlue: ", 2412 "CornflowerBlue: ", 2413 "DarkGoldenrod: ", 2414 "DarkGreen: ", 2415 "DarkKhaki: ", 2416 "DarkOliveGreen: ", 2417 "DarkOrange: ", 2418 "DarkOrchid: ", 2419 "DarkSeaGreen: ", 2420 "DarkSlateGray: ", 2421 "DarkTurquoise: ", 2422 "DeepPink: ", 2423 "DarkKhaki: ", 2424 "DimGray: ", 2425 "DodgerBlue: ", 2426 "GreenYellow: ", 2427 "HotPink: ", 2428 "IndianRed: ", 2429 "LavenderBlush: ", 2430 "LawnGreen: ", 2431 "LemonChiffon: ", 2432 "LightCoral: ", 2433 "LightCyan: ", 2434 "LightPink: ", 2435 "LightSalmon: ", 2436 "LightSlateGray: ", 2437 "LightYellow: ", 2438 "LimeGreen: ", 2439 "MediumPurple: ", 2440 "MediumSeaGreen: ", 2441 "MediumSlateBlue:", 2442 "MidnightBlue: ", 2443 "MintCream: ", 2444 "MistyRose: ", 2445 "NavajoWhite: ", 2446 "NavyBlue: ", 2447 "OliveDrab: " 2448 }; 2449 2450 /*@C 2451 PetscLogMPEGetRGBColor - This routine returns a rgb color useable with PetscLogEventRegister() 2452 2453 Not collective. Maybe it should be? 2454 2455 Output Parameter: 2456 . str - character string representing the color 2457 2458 Level: developer 2459 2460 .seealso: PetscLogEventRegister 2461 @*/ 2462 PetscErrorCode PetscLogMPEGetRGBColor(const char *str[]) 2463 { 2464 static int idx = 0; 2465 2466 PetscFunctionBegin; 2467 *str = PetscLogMPERGBColors[idx]; 2468 idx = (idx + 1)% PETSC_RGB_COLORS_MAX; 2469 PetscFunctionReturn(0); 2470 } 2471 2472 #endif /* PETSC_USE_LOG && PETSC_HAVE_MPE */ 2473