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 Level: intermediate 1070 1071 .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops() 1072 1073 M*/ 1074 1075 /*MC 1076 PetscLogEventEnd - Log the end of a user event. 1077 1078 Synopsis: 1079 #include <petsclog.h> 1080 PetscErrorCode PetscLogEventEnd(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 1081 1082 Not Collective 1083 1084 Input Parameters: 1085 + e - integer associated with the event obtained with PetscLogEventRegister() 1086 - o1,o2,o3,o4 - objects associated with the event, or 0 1087 1088 Fortran Synopsis: 1089 void PetscLogEventEnd(int e,PetscErrorCode ierr) 1090 1091 Usage: 1092 .vb 1093 PetscLogEvent USER_EVENT; 1094 PetscLogDouble user_event_flops; 1095 PetscLogEventRegister("User event",0,&USER_EVENT,); 1096 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1097 [code segment to monitor] 1098 PetscLogFlops(user_event_flops); 1099 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1100 .ve 1101 1102 Notes: 1103 You should also register each additional integer event with the command 1104 PetscLogEventRegister(). 1105 1106 Level: intermediate 1107 1108 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogFlops() 1109 1110 M*/ 1111 1112 /*@C 1113 PetscLogEventGetId - Returns the event id when given the event name. 1114 1115 Not Collective 1116 1117 Input Parameter: 1118 . name - The event name 1119 1120 Output Parameter: 1121 . event - The event, or -1 if no event with that name exists 1122 1123 Level: intermediate 1124 1125 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogStageGetId() 1126 @*/ 1127 PetscErrorCode PetscLogEventGetId(const char name[], PetscLogEvent *event) 1128 { 1129 PetscStageLog stageLog; 1130 1131 PetscFunctionBegin; 1132 PetscCall(PetscLogGetStageLog(&stageLog)); 1133 PetscCall(PetscEventRegLogGetEvent(stageLog->eventLog, name, event)); 1134 PetscFunctionReturn(0); 1135 } 1136 1137 /*------------------------------------------------ Output Functions -------------------------------------------------*/ 1138 /*@C 1139 PetscLogDump - Dumps logs of objects to a file. This file is intended to 1140 be read by bin/petscview. This program no longer exists. 1141 1142 Collective on PETSC_COMM_WORLD 1143 1144 Input Parameter: 1145 . name - an optional file name 1146 1147 Usage: 1148 .vb 1149 PetscInitialize(...); 1150 PetscLogDefaultBegin(); or PetscLogAllBegin(); 1151 ... code ... 1152 PetscLogDump(filename); 1153 PetscFinalize(); 1154 .ve 1155 1156 Notes: 1157 The default file name is 1158 $ Log.<rank> 1159 where <rank> is the processor number. If no name is specified, 1160 this file will be used. 1161 1162 Level: advanced 1163 1164 .seealso: PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogView() 1165 @*/ 1166 PetscErrorCode PetscLogDump(const char sname[]) 1167 { 1168 PetscStageLog stageLog; 1169 PetscEventPerfInfo *eventInfo; 1170 FILE *fd; 1171 char file[PETSC_MAX_PATH_LEN], fname[PETSC_MAX_PATH_LEN]; 1172 PetscLogDouble flops, _TotalTime; 1173 PetscMPIInt rank; 1174 int action, object, curStage; 1175 PetscLogEvent event; 1176 PetscErrorCode ierr; 1177 1178 PetscFunctionBegin; 1179 /* Calculate the total elapsed time */ 1180 PetscTime(&_TotalTime); 1181 _TotalTime -= petsc_BaseTime; 1182 /* Open log file */ 1183 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 1184 if (sname && sname[0]) sprintf(file, "%s.%d", sname, rank); 1185 else sprintf(file, "Log.%d", rank); 1186 PetscCall(PetscFixFilename(file, fname)); 1187 PetscCall(PetscFOpen(PETSC_COMM_WORLD, fname, "w", &fd)); 1188 PetscCheck(!(rank == 0) || !(!fd),PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN, "Cannot open file: %s", fname); 1189 /* Output totals */ 1190 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Total Flop %14e %16.8e\n", petsc_TotalFlops, _TotalTime)); 1191 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Clock Resolution %g\n", 0.0)); 1192 /* Output actions */ 1193 if (petsc_logActions) { 1194 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Actions accomplished %d\n", petsc_numActions)); 1195 for (action = 0; action < petsc_numActions; action++) { 1196 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "%g %d %d %d %d %d %d %g %g %g\n", 1197 petsc_actions[action].time, petsc_actions[action].action, (int)petsc_actions[action].event, (int)petsc_actions[action].classid, petsc_actions[action].id1, 1198 petsc_actions[action].id2, petsc_actions[action].id3, petsc_actions[action].flops, petsc_actions[action].mem, petsc_actions[action].maxmem);PetscCall(ierr); 1199 } 1200 } 1201 /* Output objects */ 1202 if (petsc_logObjects) { 1203 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Objects created %d destroyed %d\n", petsc_numObjects, petsc_numObjectsDestroyed)); 1204 for (object = 0; object < petsc_numObjects; object++) { 1205 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Parent ID: %d Memory: %d\n", petsc_objects[object].parent, (int) petsc_objects[object].mem)); 1206 if (!petsc_objects[object].name[0]) { 1207 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd,"No Name\n")); 1208 } else { 1209 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Name: %s\n", petsc_objects[object].name)); 1210 } 1211 if (petsc_objects[object].info[0] != 0) { 1212 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "No Info\n")); 1213 } else { 1214 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Info: %s\n", petsc_objects[object].info)); 1215 } 1216 } 1217 } 1218 /* Output events */ 1219 PetscCall(PetscFPrintf(PETSC_COMM_WORLD, fd, "Event log:\n")); 1220 PetscCall(PetscLogGetStageLog(&stageLog)); 1221 PetscCall(PetscIntStackTop(stageLog->stack, &curStage)); 1222 eventInfo = stageLog->stageInfo[curStage].eventLog->eventInfo; 1223 for (event = 0; event < stageLog->stageInfo[curStage].eventLog->numEvents; event++) { 1224 if (eventInfo[event].time != 0.0) flops = eventInfo[event].flops/eventInfo[event].time; 1225 else flops = 0.0; 1226 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "%d %16d %16g %16g %16g\n", event, eventInfo[event].count, 1227 eventInfo[event].flops, eventInfo[event].time, flops);PetscCall(ierr); 1228 } 1229 PetscCall(PetscFClose(PETSC_COMM_WORLD, fd)); 1230 PetscFunctionReturn(0); 1231 } 1232 1233 /* 1234 PetscLogView_Detailed - Each process prints the times for its own events 1235 1236 */ 1237 PetscErrorCode PetscLogView_Detailed(PetscViewer viewer) 1238 { 1239 PetscStageLog stageLog; 1240 PetscEventPerfInfo *eventInfo = NULL, *stageInfo = NULL; 1241 PetscLogDouble locTotalTime, numRed, maxMem; 1242 int numStages,numEvents,stage,event; 1243 MPI_Comm comm = PetscObjectComm((PetscObject) viewer); 1244 PetscMPIInt rank,size; 1245 PetscErrorCode ierr; 1246 1247 PetscFunctionBegin; 1248 PetscCallMPI(MPI_Comm_size(comm, &size)); 1249 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1250 /* Must preserve reduction count before we go on */ 1251 numRed = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1252 /* Get the total elapsed time */ 1253 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1254 PetscCall(PetscViewerASCIIPrintf(viewer,"size = %d\n",size)); 1255 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalTimes = {}\n")); 1256 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalMessages = {}\n")); 1257 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalMessageLens = {}\n")); 1258 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalReductions = {}\n")); 1259 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalFlop = {}\n")); 1260 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalObjects = {}\n")); 1261 PetscCall(PetscViewerASCIIPrintf(viewer,"LocalMemory = {}\n")); 1262 PetscCall(PetscLogGetStageLog(&stageLog)); 1263 PetscCallMPI(MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm)); 1264 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages = {}\n")); 1265 for (stage=0; stage<numStages; stage++) { 1266 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"] = {}\n",stageLog->stageInfo[stage].name)); 1267 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"summary\"] = {}\n",stageLog->stageInfo[stage].name)); 1268 PetscCallMPI(MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1269 for (event = 0; event < numEvents; event++) { 1270 PetscCall(PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"%s\"] = {}\n",stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name)); 1271 } 1272 } 1273 PetscCall(PetscMallocGetMaximumUsage(&maxMem)); 1274 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 1275 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalTimes[%d] = %g\n",rank,locTotalTime)); 1276 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessages[%d] = %g\n",rank,(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct))); 1277 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessageLens[%d] = %g\n",rank,(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len))); 1278 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalReductions[%d] = %g\n",rank,numRed)); 1279 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalFlop[%d] = %g\n",rank,petsc_TotalFlops)); 1280 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalObjects[%d] = %d\n",rank,petsc_numObjects)); 1281 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"LocalMemory[%d] = %g\n",rank,maxMem)); 1282 PetscCall(PetscViewerFlush(viewer)); 1283 for (stage=0; stage<numStages; stage++) { 1284 stageInfo = &stageLog->stageInfo[stage].perfInfo; 1285 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"summary\"][%d] = {\"time\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g}\n", 1286 stageLog->stageInfo[stage].name,rank, 1287 stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops);PetscCall(ierr); 1288 PetscCallMPI(MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1289 for (event = 0; event < numEvents; event++) { 1290 eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event]; 1291 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"%s\"][%d] = {\"count\" : %d, \"time\" : %g, \"syncTime\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g", 1292 stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name,rank, 1293 eventInfo->count,eventInfo->time,eventInfo->syncTime,eventInfo->numMessages,eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops);PetscCall(ierr); 1294 if (eventInfo->dof[0] >= 0.) { 1295 PetscInt d, e; 1296 1297 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", \"dof\" : [")); 1298 for (d = 0; d < 8; ++d) { 1299 if (d > 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", ")); 1300 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->dof[d])); 1301 } 1302 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "]")); 1303 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", \"error\" : [")); 1304 for (e = 0; e < 8; ++e) { 1305 if (e > 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", ")); 1306 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->errors[e])); 1307 } 1308 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "]")); 1309 } 1310 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"}\n")); 1311 } 1312 } 1313 PetscCall(PetscViewerFlush(viewer)); 1314 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 1315 PetscFunctionReturn(0); 1316 } 1317 1318 /* 1319 PetscLogView_CSV - Each process prints the times for its own events in Comma-Separated Value Format 1320 */ 1321 PetscErrorCode PetscLogView_CSV(PetscViewer viewer) 1322 { 1323 PetscStageLog stageLog; 1324 PetscEventPerfInfo *eventInfo = NULL; 1325 PetscLogDouble locTotalTime, maxMem; 1326 int numStages,numEvents,stage,event; 1327 MPI_Comm comm = PetscObjectComm((PetscObject) viewer); 1328 PetscMPIInt rank,size; 1329 PetscErrorCode ierr; 1330 1331 PetscFunctionBegin; 1332 PetscCallMPI(MPI_Comm_size(comm, &size)); 1333 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1334 /* Must preserve reduction count before we go on */ 1335 /* Get the total elapsed time */ 1336 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1337 PetscCall(PetscLogGetStageLog(&stageLog)); 1338 PetscCallMPI(MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm)); 1339 PetscCall(PetscMallocGetMaximumUsage(&maxMem)); 1340 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 1341 ierr = 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); 1342 PetscCall(PetscViewerFlush(viewer)); 1343 for (stage=0; stage<numStages; stage++) { 1344 PetscEventPerfInfo *stageInfo = &stageLog->stageInfo[stage].perfInfo; 1345 1346 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"%s,summary,%d,1,%g,%g,%g,%g,%g\n", 1347 stageLog->stageInfo[stage].name,rank,stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops);PetscCall(ierr); 1348 PetscCallMPI(MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1349 for (event = 0; event < numEvents; event++) { 1350 eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event]; 1351 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"%s,%s,%d,%d,%g,%g,%g,%g,%g",stageLog->stageInfo[stage].name, 1352 stageLog->eventLog->eventInfo[event].name,rank,eventInfo->count,eventInfo->time,eventInfo->numMessages, 1353 eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops);PetscCall(ierr); 1354 if (eventInfo->dof[0] >= 0.) { 1355 PetscInt d, e; 1356 1357 for (d = 0; d < 8; ++d) { 1358 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->dof[d])); 1359 } 1360 for (e = 0; e < 8; ++e) { 1361 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->errors[e])); 1362 } 1363 } 1364 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer,"\n")); 1365 } 1366 } 1367 PetscCall(PetscViewerFlush(viewer)); 1368 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 1369 PetscFunctionReturn(0); 1370 } 1371 1372 static PetscErrorCode PetscLogViewWarnSync(MPI_Comm comm,FILE *fd) 1373 { 1374 PetscFunctionBegin; 1375 if (!PetscLogSyncOn) PetscFunctionReturn(0); 1376 PetscCall(PetscFPrintf(comm, fd, "\n\n")); 1377 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n")); 1378 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1379 PetscCall(PetscFPrintf(comm, fd, " # WARNING!!! #\n")); 1380 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1381 PetscCall(PetscFPrintf(comm, fd, " # This program was run with logging synchronization. #\n")); 1382 PetscCall(PetscFPrintf(comm, fd, " # This option provides more meaningful imbalance #\n")); 1383 PetscCall(PetscFPrintf(comm, fd, " # figures at the expense of slowing things down and #\n")); 1384 PetscCall(PetscFPrintf(comm, fd, " # providing a distorted view of the overall runtime. #\n")); 1385 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1386 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n\n\n")); 1387 PetscFunctionReturn(0); 1388 } 1389 1390 static PetscErrorCode PetscLogViewWarnDebugging(MPI_Comm comm,FILE *fd) 1391 { 1392 PetscFunctionBegin; 1393 if (PetscDefined(USE_DEBUG)) { 1394 PetscCall(PetscFPrintf(comm, fd, "\n\n")); 1395 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n")); 1396 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1397 PetscCall(PetscFPrintf(comm, fd, " # WARNING!!! #\n")); 1398 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1399 PetscCall(PetscFPrintf(comm, fd, " # This code was compiled with a debugging option. #\n")); 1400 PetscCall(PetscFPrintf(comm, fd, " # To get timing results run ./configure #\n")); 1401 PetscCall(PetscFPrintf(comm, fd, " # using --with-debugging=no, the performance will #\n")); 1402 PetscCall(PetscFPrintf(comm, fd, " # be generally two or three times faster. #\n")); 1403 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1404 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n\n\n")); 1405 } 1406 PetscFunctionReturn(0); 1407 } 1408 1409 static PetscErrorCode PetscLogViewWarnNoGpuAwareMpi(MPI_Comm comm,FILE *fd) 1410 { 1411 #if defined(PETSC_HAVE_DEVICE) 1412 PetscMPIInt size; 1413 1414 PetscFunctionBegin; 1415 PetscCallMPI(MPI_Comm_size(comm, &size)); 1416 if (use_gpu_aware_mpi || size == 1) PetscFunctionReturn(0); 1417 PetscCall(PetscFPrintf(comm, fd, "\n\n")); 1418 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n")); 1419 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1420 PetscCall(PetscFPrintf(comm, fd, " # WARNING!!! #\n")); 1421 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1422 PetscCall(PetscFPrintf(comm, fd, " # This code was compiled with GPU support and you've #\n")); 1423 PetscCall(PetscFPrintf(comm, fd, " # created PETSc/GPU objects, but you intentionally #\n")); 1424 PetscCall(PetscFPrintf(comm, fd, " # used -use_gpu_aware_mpi 0, requiring PETSc to copy #\n")); 1425 PetscCall(PetscFPrintf(comm, fd, " # additional data between the GPU and CPU. To obtain #\n")); 1426 PetscCall(PetscFPrintf(comm, fd, " # meaningful timing results on multi-rank runs, use #\n")); 1427 PetscCall(PetscFPrintf(comm, fd, " # GPU-aware MPI instead. #\n")); 1428 PetscCall(PetscFPrintf(comm, fd, " # #\n")); 1429 PetscCall(PetscFPrintf(comm, fd, " ##########################################################\n\n\n")); 1430 PetscFunctionReturn(0); 1431 #else 1432 return 0; 1433 #endif 1434 } 1435 1436 PetscErrorCode PetscLogView_Default(PetscViewer viewer) 1437 { 1438 FILE *fd; 1439 PetscLogDouble zero = 0.0; 1440 PetscStageLog stageLog; 1441 PetscStageInfo *stageInfo = NULL; 1442 PetscEventPerfInfo *eventInfo = NULL; 1443 PetscClassPerfInfo *classInfo; 1444 char arch[128],hostname[128],username[128],pname[PETSC_MAX_PATH_LEN],date[128]; 1445 const char *name; 1446 PetscLogDouble locTotalTime, TotalTime, TotalFlops; 1447 PetscLogDouble numMessages, messageLength, avgMessLen, numReductions; 1448 PetscLogDouble stageTime, flops, flopr, mem, mess, messLen, red; 1449 PetscLogDouble fracTime, fracFlops, fracMessages, fracLength, fracReductions, fracMess, fracMessLen, fracRed; 1450 PetscLogDouble fracStageTime, fracStageFlops, fracStageMess, fracStageMessLen, fracStageRed; 1451 PetscLogDouble min, max, tot, ratio, avg, x, y; 1452 PetscLogDouble minf, maxf, totf, ratf, mint, maxt, tott, ratt, ratC, totm, totml, totr, mal, malmax, emalmax; 1453 #if defined(PETSC_HAVE_DEVICE) 1454 PetscLogDouble cct, gct, csz, gsz, gmaxt, gflops, gflopr, fracgflops; 1455 #endif 1456 PetscMPIInt minC, maxC; 1457 PetscMPIInt size, rank; 1458 PetscBool *localStageUsed, *stageUsed; 1459 PetscBool *localStageVisible, *stageVisible; 1460 int numStages, localNumEvents, numEvents; 1461 int stage, oclass; 1462 PetscLogEvent event; 1463 PetscErrorCode ierr; 1464 char version[256]; 1465 MPI_Comm comm; 1466 1467 PetscFunctionBegin; 1468 PetscCall(PetscObjectGetComm((PetscObject)viewer,&comm)); 1469 PetscCall(PetscViewerASCIIGetPointer(viewer,&fd)); 1470 PetscCallMPI(MPI_Comm_size(comm, &size)); 1471 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1472 /* Get the total elapsed time */ 1473 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1474 1475 PetscCall(PetscFPrintf(comm, fd, "**************************************** ***********************************************************************************************************************\n")); 1476 PetscCall(PetscFPrintf(comm, fd, "*** WIDEN YOUR WINDOW TO 160 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***\n")); 1477 PetscCall(PetscFPrintf(comm, fd, "****************************************************************************************************************************************************************\n")); 1478 PetscCall(PetscFPrintf(comm, fd, "\n------------------------------------------------------------------ PETSc Performance Summary: -------------------------------------------------------------------\n\n")); 1479 PetscCall(PetscLogViewWarnSync(comm,fd)); 1480 PetscCall(PetscLogViewWarnDebugging(comm,fd)); 1481 PetscCall(PetscLogViewWarnNoGpuAwareMpi(comm,fd)); 1482 PetscCall(PetscGetArchType(arch,sizeof(arch))); 1483 PetscCall(PetscGetHostName(hostname,sizeof(hostname))); 1484 PetscCall(PetscGetUserName(username,sizeof(username))); 1485 PetscCall(PetscGetProgramName(pname,sizeof(pname))); 1486 PetscCall(PetscGetDate(date,sizeof(date))); 1487 PetscCall(PetscGetVersion(version,sizeof(version))); 1488 if (size == 1) { 1489 PetscCall(PetscFPrintf(comm,fd,"%s on a %s named %s with %d processor, by %s %s\n", pname, arch, hostname, size, username, date)); 1490 } else { 1491 PetscCall(PetscFPrintf(comm,fd,"%s on a %s named %s with %d processors, by %s %s\n", pname, arch, hostname, size, username, date)); 1492 } 1493 #if defined(PETSC_HAVE_OPENMP) 1494 PetscCall(PetscFPrintf(comm,fd,"Using %" PetscInt_FMT " OpenMP threads\n", PetscNumOMPThreads)); 1495 #endif 1496 PetscCall(PetscFPrintf(comm, fd, "Using %s\n", version)); 1497 1498 /* Must preserve reduction count before we go on */ 1499 red = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1500 1501 /* Calculate summary information */ 1502 PetscCall(PetscFPrintf(comm, fd, "\n Max Max/Min Avg Total\n")); 1503 /* Time */ 1504 PetscCallMPI(MPI_Allreduce(&locTotalTime, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1505 PetscCallMPI(MPI_Allreduce(&locTotalTime, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1506 PetscCallMPI(MPI_Allreduce(&locTotalTime, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1507 avg = tot/((PetscLogDouble) size); 1508 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1509 PetscCall(PetscFPrintf(comm, fd, "Time (sec): %5.3e %7.3f %5.3e\n", max, ratio, avg)); 1510 TotalTime = tot; 1511 /* Objects */ 1512 avg = (PetscLogDouble) petsc_numObjects; 1513 PetscCallMPI(MPI_Allreduce(&avg, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1514 PetscCallMPI(MPI_Allreduce(&avg, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1515 PetscCallMPI(MPI_Allreduce(&avg, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1516 avg = tot/((PetscLogDouble) size); 1517 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1518 PetscCall(PetscFPrintf(comm, fd, "Objects: %5.3e %7.3f %5.3e\n", max, ratio, avg)); 1519 /* Flops */ 1520 PetscCallMPI(MPI_Allreduce(&petsc_TotalFlops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1521 PetscCallMPI(MPI_Allreduce(&petsc_TotalFlops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1522 PetscCallMPI(MPI_Allreduce(&petsc_TotalFlops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1523 avg = tot/((PetscLogDouble) size); 1524 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1525 PetscCall(PetscFPrintf(comm, fd, "Flops: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1526 TotalFlops = tot; 1527 /* Flops/sec -- Must talk to Barry here */ 1528 if (locTotalTime != 0.0) flops = petsc_TotalFlops/locTotalTime; else flops = 0.0; 1529 PetscCallMPI(MPI_Allreduce(&flops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1530 PetscCallMPI(MPI_Allreduce(&flops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1531 PetscCallMPI(MPI_Allreduce(&flops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1532 avg = tot/((PetscLogDouble) size); 1533 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1534 PetscCall(PetscFPrintf(comm, fd, "Flops/sec: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1535 /* Memory */ 1536 PetscCall(PetscMallocGetMaximumUsage(&mem)); 1537 if (mem > 0.0) { 1538 PetscCallMPI(MPI_Allreduce(&mem, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1539 PetscCallMPI(MPI_Allreduce(&mem, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1540 PetscCallMPI(MPI_Allreduce(&mem, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1541 avg = tot/((PetscLogDouble) size); 1542 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1543 PetscCall(PetscFPrintf(comm, fd, "Memory (bytes): %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1544 } 1545 /* Messages */ 1546 mess = 0.5*(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct); 1547 PetscCallMPI(MPI_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1548 PetscCallMPI(MPI_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1549 PetscCallMPI(MPI_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1550 avg = tot/((PetscLogDouble) size); 1551 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1552 PetscCall(PetscFPrintf(comm, fd, "MPI Msg Count: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1553 numMessages = tot; 1554 /* Message Lengths */ 1555 mess = 0.5*(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len); 1556 PetscCallMPI(MPI_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1557 PetscCallMPI(MPI_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1558 PetscCallMPI(MPI_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1559 if (numMessages != 0) avg = tot/numMessages; else avg = 0.0; 1560 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1561 PetscCall(PetscFPrintf(comm, fd, "MPI Msg Len (bytes): %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1562 messageLength = tot; 1563 /* Reductions */ 1564 PetscCallMPI(MPI_Allreduce(&red, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1565 PetscCallMPI(MPI_Allreduce(&red, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1566 PetscCallMPI(MPI_Allreduce(&red, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1567 if (min != 0.0) ratio = max/min; else ratio = 0.0; 1568 PetscCall(PetscFPrintf(comm, fd, "MPI Reductions: %5.3e %7.3f\n", max, ratio)); 1569 numReductions = red; /* wrong because uses count from process zero */ 1570 PetscCall(PetscFPrintf(comm, fd, "\nFlop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)\n")); 1571 PetscCall(PetscFPrintf(comm, fd, " e.g., VecAXPY() for real vectors of length N --> 2N flops\n")); 1572 PetscCall(PetscFPrintf(comm, fd, " and VecAXPY() for complex vectors of length N --> 8N flops\n")); 1573 1574 /* Get total number of stages -- 1575 Currently, a single processor can register more stages than another, but stages must all be registered in order. 1576 We can removed this requirement if necessary by having a global stage numbering and indirection on the stage ID. 1577 This seems best accomplished by assoicating a communicator with each stage. 1578 */ 1579 PetscCall(PetscLogGetStageLog(&stageLog)); 1580 PetscCallMPI(MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm)); 1581 PetscCall(PetscMalloc1(numStages, &localStageUsed)); 1582 PetscCall(PetscMalloc1(numStages, &stageUsed)); 1583 PetscCall(PetscMalloc1(numStages, &localStageVisible)); 1584 PetscCall(PetscMalloc1(numStages, &stageVisible)); 1585 if (numStages > 0) { 1586 stageInfo = stageLog->stageInfo; 1587 for (stage = 0; stage < numStages; stage++) { 1588 if (stage < stageLog->numStages) { 1589 localStageUsed[stage] = stageInfo[stage].used; 1590 localStageVisible[stage] = stageInfo[stage].perfInfo.visible; 1591 } else { 1592 localStageUsed[stage] = PETSC_FALSE; 1593 localStageVisible[stage] = PETSC_TRUE; 1594 } 1595 } 1596 PetscCallMPI(MPI_Allreduce(localStageUsed, stageUsed, numStages, MPIU_BOOL, MPI_LOR, comm)); 1597 PetscCallMPI(MPI_Allreduce(localStageVisible, stageVisible, numStages, MPIU_BOOL, MPI_LAND, comm)); 1598 for (stage = 0; stage < numStages; stage++) { 1599 if (stageUsed[stage]) { 1600 PetscCall(PetscFPrintf(comm, fd, "\nSummary of Stages: ----- Time ------ ----- Flop ------ --- Messages --- -- Message Lengths -- -- Reductions --\n")); 1601 PetscCall(PetscFPrintf(comm, fd, " Avg %%Total Avg %%Total Count %%Total Avg %%Total Count %%Total\n")); 1602 break; 1603 } 1604 } 1605 for (stage = 0; stage < numStages; stage++) { 1606 if (!stageUsed[stage]) continue; 1607 /* CANNOT use MPI_Allreduce() since it might fail the line number check */ 1608 if (localStageUsed[stage]) { 1609 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1610 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1611 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1612 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1613 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1614 name = stageInfo[stage].name; 1615 } else { 1616 PetscCallMPI(MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1617 PetscCallMPI(MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1618 PetscCallMPI(MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1619 PetscCallMPI(MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1620 PetscCallMPI(MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1621 name = ""; 1622 } 1623 mess *= 0.5; messLen *= 0.5; red /= size; 1624 if (TotalTime != 0.0) fracTime = stageTime/TotalTime; else fracTime = 0.0; 1625 if (TotalFlops != 0.0) fracFlops = flops/TotalFlops; else fracFlops = 0.0; 1626 /* Talk to Barry if (stageTime != 0.0) flops = (size*flops)/stageTime; else flops = 0.0; */ 1627 if (numMessages != 0.0) fracMessages = mess/numMessages; else fracMessages = 0.0; 1628 if (mess != 0.0) avgMessLen = messLen/mess; else avgMessLen = 0.0; 1629 if (messageLength != 0.0) fracLength = messLen/messageLength; else fracLength = 0.0; 1630 if (numReductions != 0.0) fracReductions = red/numReductions; else fracReductions = 0.0; 1631 ierr = 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", 1632 stage, name, stageTime/size, 100.0*fracTime, flops, 100.0*fracFlops, 1633 mess, 100.0*fracMessages, avgMessLen, 100.0*fracLength, red, 100.0*fracReductions);PetscCall(ierr); 1634 } 1635 } 1636 1637 PetscCall(PetscFPrintf(comm, fd,"\n------------------------------------------------------------------------------------------------------------------------\n")); 1638 PetscCall(PetscFPrintf(comm, fd, "See the 'Profiling' chapter of the users' manual for details on interpreting output.\n")); 1639 PetscCall(PetscFPrintf(comm, fd, "Phase summary info:\n")); 1640 PetscCall(PetscFPrintf(comm, fd, " Count: number of times phase was executed\n")); 1641 PetscCall(PetscFPrintf(comm, fd, " Time and Flop: Max - maximum over all processors\n")); 1642 PetscCall(PetscFPrintf(comm, fd, " Ratio - ratio of maximum to minimum over all processors\n")); 1643 PetscCall(PetscFPrintf(comm, fd, " Mess: number of messages sent\n")); 1644 PetscCall(PetscFPrintf(comm, fd, " AvgLen: average message length (bytes)\n")); 1645 PetscCall(PetscFPrintf(comm, fd, " Reduct: number of global reductions\n")); 1646 PetscCall(PetscFPrintf(comm, fd, " Global: entire computation\n")); 1647 PetscCall(PetscFPrintf(comm, fd, " Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().\n")); 1648 PetscCall(PetscFPrintf(comm, fd, " %%T - percent time in this phase %%F - percent flop in this phase\n")); 1649 PetscCall(PetscFPrintf(comm, fd, " %%M - percent messages in this phase %%L - percent message lengths in this phase\n")); 1650 PetscCall(PetscFPrintf(comm, fd, " %%R - percent reductions in this phase\n")); 1651 PetscCall(PetscFPrintf(comm, fd, " Total Mflop/s: 10e-6 * (sum of flop over all processors)/(max time over all processors)\n")); 1652 if (PetscLogMemory) { 1653 PetscCall(PetscFPrintf(comm, fd, " Malloc Mbytes: Memory allocated and kept during event (sum over all calls to event)\n")); 1654 PetscCall(PetscFPrintf(comm, fd, " EMalloc Mbytes: extra memory allocated during event and then freed (maximum over all calls to events)\n")); 1655 PetscCall(PetscFPrintf(comm, fd, " MMalloc Mbytes: Increase in high water mark of allocated memory (sum over all calls to event)\n")); 1656 PetscCall(PetscFPrintf(comm, fd, " RMI Mbytes: Increase in resident memory (sum over all calls to event)\n")); 1657 } 1658 #if defined(PETSC_HAVE_DEVICE) 1659 PetscCall(PetscFPrintf(comm, fd, " GPU Mflop/s: 10e-6 * (sum of flop on GPU over all processors)/(max GPU time over all processors)\n")); 1660 PetscCall(PetscFPrintf(comm, fd, " CpuToGpu Count: total number of CPU to GPU copies per processor\n")); 1661 PetscCall(PetscFPrintf(comm, fd, " CpuToGpu Size (Mbytes): 10e-6 * (total size of CPU to GPU copies per processor)\n")); 1662 PetscCall(PetscFPrintf(comm, fd, " GpuToCpu Count: total number of GPU to CPU copies per processor\n")); 1663 PetscCall(PetscFPrintf(comm, fd, " GpuToCpu Size (Mbytes): 10e-6 * (total size of GPU to CPU copies per processor)\n")); 1664 PetscCall(PetscFPrintf(comm, fd, " GPU %%F: percent flops on GPU in this event\n")); 1665 #endif 1666 PetscCall(PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------\n")); 1667 1668 PetscCall(PetscLogViewWarnDebugging(comm,fd)); 1669 1670 /* Report events */ 1671 PetscCall(PetscFPrintf(comm, fd,"Event Count Time (sec) Flop --- Global --- --- Stage ---- Total")); 1672 if (PetscLogMemory) { 1673 PetscCall(PetscFPrintf(comm, fd," Malloc EMalloc MMalloc RMI")); 1674 } 1675 #if defined(PETSC_HAVE_DEVICE) 1676 PetscCall(PetscFPrintf(comm, fd," GPU - CpuToGpu - - GpuToCpu - GPU")); 1677 #endif 1678 PetscCall(PetscFPrintf(comm, fd,"\n")); 1679 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")); 1680 if (PetscLogMemory) { 1681 PetscCall(PetscFPrintf(comm, fd," Mbytes Mbytes Mbytes Mbytes")); 1682 } 1683 #if defined(PETSC_HAVE_DEVICE) 1684 PetscCall(PetscFPrintf(comm, fd," Mflop/s Count Size Count Size %%F")); 1685 #endif 1686 PetscCall(PetscFPrintf(comm, fd,"\n")); 1687 PetscCall(PetscFPrintf(comm, fd,"------------------------------------------------------------------------------------------------------------------------")); 1688 if (PetscLogMemory) { 1689 PetscCall(PetscFPrintf(comm, fd,"-----------------------------")); 1690 } 1691 #if defined(PETSC_HAVE_DEVICE) 1692 PetscCall(PetscFPrintf(comm, fd,"---------------------------------------")); 1693 #endif 1694 PetscCall(PetscFPrintf(comm, fd,"\n")); 1695 1696 /* Problem: The stage name will not show up unless the stage executed on proc 1 */ 1697 for (stage = 0; stage < numStages; stage++) { 1698 if (!stageVisible[stage]) continue; 1699 /* CANNOT use MPI_Allreduce() since it might fail the line number check */ 1700 if (localStageUsed[stage]) { 1701 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name)); 1702 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1703 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1704 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1705 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1706 PetscCallMPI(MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1707 } else { 1708 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage)); 1709 PetscCallMPI(MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1710 PetscCallMPI(MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1711 PetscCallMPI(MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1712 PetscCallMPI(MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1713 PetscCallMPI(MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1714 } 1715 mess *= 0.5; messLen *= 0.5; red /= size; 1716 1717 /* Get total number of events in this stage -- 1718 Currently, a single processor can register more events than another, but events must all be registered in order, 1719 just like stages. We can removed this requirement if necessary by having a global event numbering and indirection 1720 on the event ID. This seems best accomplished by associating a communicator with each stage. 1721 1722 Problem: If the event did not happen on proc 1, its name will not be available. 1723 Problem: Event visibility is not implemented 1724 */ 1725 if (localStageUsed[stage]) { 1726 eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo; 1727 localNumEvents = stageLog->stageInfo[stage].eventLog->numEvents; 1728 } else localNumEvents = 0; 1729 PetscCallMPI(MPI_Allreduce(&localNumEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm)); 1730 for (event = 0; event < numEvents; event++) { 1731 /* CANNOT use MPI_Allreduce() since it might fail the line number check */ 1732 if (localStageUsed[stage] && (event < stageLog->stageInfo[stage].eventLog->numEvents) && (eventInfo[event].depth == 0)) { 1733 if ((eventInfo[event].count > 0) && (eventInfo[event].time > 0.0)) flopr = eventInfo[event].flops; else flopr = 0.0; 1734 PetscCallMPI(MPI_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1735 PetscCallMPI(MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1736 PetscCallMPI(MPI_Allreduce(&eventInfo[event].flops, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1737 PetscCallMPI(MPI_Allreduce(&eventInfo[event].time, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1738 PetscCallMPI(MPI_Allreduce(&eventInfo[event].time, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1739 PetscCallMPI(MPI_Allreduce(&eventInfo[event].time, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1740 PetscCallMPI(MPI_Allreduce(&eventInfo[event].numMessages, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1741 PetscCallMPI(MPI_Allreduce(&eventInfo[event].messageLength, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1742 PetscCallMPI(MPI_Allreduce(&eventInfo[event].numReductions, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1743 PetscCallMPI(MPI_Allreduce(&eventInfo[event].count, &minC, 1, MPI_INT, MPI_MIN, comm)); 1744 PetscCallMPI(MPI_Allreduce(&eventInfo[event].count, &maxC, 1, MPI_INT, MPI_MAX, comm)); 1745 if (PetscLogMemory) { 1746 PetscCallMPI(MPI_Allreduce(&eventInfo[event].memIncrease, &mem, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1747 PetscCallMPI(MPI_Allreduce(&eventInfo[event].mallocSpace, &mal, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1748 PetscCallMPI(MPI_Allreduce(&eventInfo[event].mallocIncrease, &malmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1749 PetscCallMPI(MPI_Allreduce(&eventInfo[event].mallocIncreaseEvent, &emalmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1750 } 1751 #if defined(PETSC_HAVE_DEVICE) 1752 PetscCallMPI(MPI_Allreduce(&eventInfo[event].CpuToGpuCount, &cct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1753 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuToCpuCount, &gct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1754 PetscCallMPI(MPI_Allreduce(&eventInfo[event].CpuToGpuSize, &csz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1755 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuToCpuSize, &gsz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1756 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuFlops, &gflops,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1757 PetscCallMPI(MPI_Allreduce(&eventInfo[event].GpuTime, &gmaxt ,1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1758 #endif 1759 name = stageLog->eventLog->eventInfo[event].name; 1760 } else { 1761 flopr = 0.0; 1762 PetscCallMPI(MPI_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1763 PetscCallMPI(MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1764 PetscCallMPI(MPI_Allreduce(&zero, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1765 PetscCallMPI(MPI_Allreduce(&zero, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1766 PetscCallMPI(MPI_Allreduce(&zero, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1767 PetscCallMPI(MPI_Allreduce(&zero, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1768 PetscCallMPI(MPI_Allreduce(&zero, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1769 PetscCallMPI(MPI_Allreduce(&zero, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1770 PetscCallMPI(MPI_Allreduce(&zero, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1771 PetscCallMPI(MPI_Allreduce(&ierr, &minC, 1, MPI_INT, MPI_MIN, comm)); 1772 PetscCallMPI(MPI_Allreduce(&ierr, &maxC, 1, MPI_INT, MPI_MAX, comm)); 1773 if (PetscLogMemory) { 1774 PetscCallMPI(MPI_Allreduce(&zero, &mem, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1775 PetscCallMPI(MPI_Allreduce(&zero, &mal, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1776 PetscCallMPI(MPI_Allreduce(&zero, &malmax, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1777 PetscCallMPI(MPI_Allreduce(&zero, &emalmax,1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1778 } 1779 #if defined(PETSC_HAVE_DEVICE) 1780 PetscCallMPI(MPI_Allreduce(&zero, &cct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1781 PetscCallMPI(MPI_Allreduce(&zero, &gct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1782 PetscCallMPI(MPI_Allreduce(&zero, &csz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1783 PetscCallMPI(MPI_Allreduce(&zero, &gsz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1784 PetscCallMPI(MPI_Allreduce(&zero, &gflops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1785 PetscCallMPI(MPI_Allreduce(&zero, &gmaxt , 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1786 #endif 1787 name = ""; 1788 } 1789 if (mint < 0.0) { 1790 ierr = 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); 1791 mint = 0; 1792 } 1793 PetscCheck(minf >= 0.0,PETSC_COMM_SELF,PETSC_ERR_PLIB,"Minimum flop %g over all processors for %s is negative! Not possible!",minf,name); 1794 totm *= 0.5; totml *= 0.5; totr /= size; 1795 1796 if (maxC != 0) { 1797 if (minC != 0) ratC = ((PetscLogDouble)maxC)/minC;else ratC = 0.0; 1798 if (mint != 0.0) ratt = maxt/mint; else ratt = 0.0; 1799 if (minf != 0.0) ratf = maxf/minf; else ratf = 0.0; 1800 if (TotalTime != 0.0) fracTime = tott/TotalTime; else fracTime = 0.0; 1801 if (TotalFlops != 0.0) fracFlops = totf/TotalFlops; else fracFlops = 0.0; 1802 if (stageTime != 0.0) fracStageTime = tott/stageTime; else fracStageTime = 0.0; 1803 if (flops != 0.0) fracStageFlops = totf/flops; else fracStageFlops = 0.0; 1804 if (numMessages != 0.0) fracMess = totm/numMessages; else fracMess = 0.0; 1805 if (messageLength != 0.0) fracMessLen = totml/messageLength; else fracMessLen = 0.0; 1806 if (numReductions != 0.0) fracRed = totr/numReductions; else fracRed = 0.0; 1807 if (mess != 0.0) fracStageMess = totm/mess; else fracStageMess = 0.0; 1808 if (messLen != 0.0) fracStageMessLen = totml/messLen; else fracStageMessLen = 0.0; 1809 if (red != 0.0) fracStageRed = totr/red; else fracStageRed = 0.0; 1810 if (totm != 0.0) totml /= totm; else totml = 0.0; 1811 if (maxt != 0.0) flopr = totf/maxt; else flopr = 0.0; 1812 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")); 1813 ierr = PetscFPrintf(comm, fd, 1814 "%-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", 1815 name, maxC, ratC, maxt, ratt, maxf, ratf, totm, totml, totr, 1816 100.0*fracTime, 100.0*fracFlops, 100.0*fracMess, 100.0*fracMessLen, 100.0*fracRed, 1817 100.0*fracStageTime, 100.0*fracStageFlops, 100.0*fracStageMess, 100.0*fracStageMessLen, 100.0*fracStageRed, 1818 PetscAbs(flopr)/1.0e6);PetscCall(ierr); 1819 if (PetscLogMemory) { 1820 PetscCall(PetscFPrintf(comm, fd," %5.0f %5.0f %5.0f %5.0f",mal/1.0e6,emalmax/1.0e6,malmax/1.0e6,mem/1.0e6)); 1821 } 1822 #if defined(PETSC_HAVE_DEVICE) 1823 if (totf != 0.0) fracgflops = gflops/totf; else fracgflops = 0.0; 1824 if (gmaxt != 0.0) gflopr = gflops/gmaxt; else gflopr = 0.0; 1825 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)); 1826 #endif 1827 PetscCall(PetscFPrintf(comm, fd,"\n")); 1828 } 1829 } 1830 } 1831 1832 /* Memory usage and object creation */ 1833 PetscCall(PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------")); 1834 if (PetscLogMemory) { 1835 PetscCall(PetscFPrintf(comm, fd, "-----------------------------")); 1836 } 1837 #if defined(PETSC_HAVE_DEVICE) 1838 PetscCall(PetscFPrintf(comm, fd, "---------------------------------------")); 1839 #endif 1840 PetscCall(PetscFPrintf(comm, fd, "\n")); 1841 PetscCall(PetscFPrintf(comm, fd, "\n")); 1842 PetscCall(PetscFPrintf(comm, fd, "Memory usage is given in bytes:\n\n")); 1843 1844 /* Right now, only stages on the first processor are reported here, meaning only objects associated with 1845 the global communicator, or MPI_COMM_SELF for proc 1. We really should report global stats and then 1846 stats for stages local to processor sets. 1847 */ 1848 /* We should figure out the longest object name here (now 20 characters) */ 1849 PetscCall(PetscFPrintf(comm, fd, "Object Type Creations Destructions Memory Descendants' Mem.\n")); 1850 PetscCall(PetscFPrintf(comm, fd, "Reports information only for process 0.\n")); 1851 for (stage = 0; stage < numStages; stage++) { 1852 if (localStageUsed[stage]) { 1853 classInfo = stageLog->stageInfo[stage].classLog->classInfo; 1854 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name)); 1855 for (oclass = 0; oclass < stageLog->stageInfo[stage].classLog->numClasses; oclass++) { 1856 if ((classInfo[oclass].creations > 0) || (classInfo[oclass].destructions > 0)) { 1857 ierr = PetscFPrintf(comm, fd, "%20s %5d %5d %11.0f %g\n", stageLog->classLog->classInfo[oclass].name, 1858 classInfo[oclass].creations, classInfo[oclass].destructions, classInfo[oclass].mem, 1859 classInfo[oclass].descMem);PetscCall(ierr); 1860 } 1861 } 1862 } else { 1863 if (!localStageVisible[stage]) continue; 1864 PetscCall(PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage)); 1865 } 1866 } 1867 1868 PetscCall(PetscFree(localStageUsed)); 1869 PetscCall(PetscFree(stageUsed)); 1870 PetscCall(PetscFree(localStageVisible)); 1871 PetscCall(PetscFree(stageVisible)); 1872 1873 /* Information unrelated to this particular run */ 1874 PetscCall(PetscFPrintf(comm, fd, "========================================================================================================================\n")); 1875 PetscTime(&y); 1876 PetscTime(&x); 1877 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 1878 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 1879 PetscCall(PetscFPrintf(comm,fd,"Average time to get PetscTime(): %g\n", (y-x)/10.0)); 1880 /* MPI information */ 1881 if (size > 1) { 1882 MPI_Status status; 1883 PetscMPIInt tag; 1884 MPI_Comm newcomm; 1885 1886 PetscCallMPI(MPI_Barrier(comm)); 1887 PetscTime(&x); 1888 PetscCallMPI(MPI_Barrier(comm)); 1889 PetscCallMPI(MPI_Barrier(comm)); 1890 PetscCallMPI(MPI_Barrier(comm)); 1891 PetscCallMPI(MPI_Barrier(comm)); 1892 PetscCallMPI(MPI_Barrier(comm)); 1893 PetscTime(&y); 1894 PetscCall(PetscFPrintf(comm, fd, "Average time for MPI_Barrier(): %g\n", (y-x)/5.0)); 1895 PetscCall(PetscCommDuplicate(comm,&newcomm, &tag)); 1896 PetscCallMPI(MPI_Barrier(comm)); 1897 if (rank) { 1898 PetscCallMPI(MPI_Recv(NULL, 0, MPI_INT, rank-1, tag, newcomm, &status)); 1899 PetscCallMPI(MPI_Send(NULL, 0, MPI_INT, (rank+1)%size, tag, newcomm)); 1900 } else { 1901 PetscTime(&x); 1902 PetscCallMPI(MPI_Send(NULL, 0, MPI_INT, 1, tag, newcomm)); 1903 PetscCallMPI(MPI_Recv(NULL, 0, MPI_INT, size-1, tag, newcomm, &status)); 1904 PetscTime(&y); 1905 PetscCall(PetscFPrintf(comm,fd,"Average time for zero size MPI_Send(): %g\n", (y-x)/size)); 1906 } 1907 PetscCall(PetscCommDestroy(&newcomm)); 1908 } 1909 PetscCall(PetscOptionsView(NULL,viewer)); 1910 1911 /* Machine and compile information */ 1912 #if defined(PETSC_USE_FORTRAN_KERNELS) 1913 PetscCall(PetscFPrintf(comm, fd, "Compiled with FORTRAN kernels\n")); 1914 #else 1915 PetscCall(PetscFPrintf(comm, fd, "Compiled without FORTRAN kernels\n")); 1916 #endif 1917 #if defined(PETSC_USE_64BIT_INDICES) 1918 PetscCall(PetscFPrintf(comm, fd, "Compiled with 64 bit PetscInt\n")); 1919 #elif defined(PETSC_USE___FLOAT128) 1920 PetscCall(PetscFPrintf(comm, fd, "Compiled with 32 bit PetscInt\n")); 1921 #endif 1922 #if defined(PETSC_USE_REAL_SINGLE) 1923 PetscCall(PetscFPrintf(comm, fd, "Compiled with single precision PetscScalar and PetscReal\n")); 1924 #elif defined(PETSC_USE___FLOAT128) 1925 PetscCall(PetscFPrintf(comm, fd, "Compiled with 128 bit precision PetscScalar and PetscReal\n")); 1926 #endif 1927 #if defined(PETSC_USE_REAL_MAT_SINGLE) 1928 PetscCall(PetscFPrintf(comm, fd, "Compiled with single precision matrices\n")); 1929 #else 1930 PetscCall(PetscFPrintf(comm, fd, "Compiled with full precision matrices (default)\n")); 1931 #endif 1932 ierr = PetscFPrintf(comm, fd, "sizeof(short) %d sizeof(int) %d sizeof(long) %d sizeof(void*) %d sizeof(PetscScalar) %d sizeof(PetscInt) %d\n", 1933 (int) sizeof(short), (int) sizeof(int), (int) sizeof(long), (int) sizeof(void*),(int) sizeof(PetscScalar),(int) sizeof(PetscInt));PetscCall(ierr); 1934 1935 PetscCall(PetscFPrintf(comm, fd, "Configure options: %s",petscconfigureoptions)); 1936 PetscCall(PetscFPrintf(comm, fd, "%s", petscmachineinfo)); 1937 PetscCall(PetscFPrintf(comm, fd, "%s", petsccompilerinfo)); 1938 PetscCall(PetscFPrintf(comm, fd, "%s", petsccompilerflagsinfo)); 1939 PetscCall(PetscFPrintf(comm, fd, "%s", petsclinkerinfo)); 1940 1941 /* Cleanup */ 1942 PetscCall(PetscFPrintf(comm, fd, "\n")); 1943 PetscCall(PetscLogViewWarnNoGpuAwareMpi(comm,fd)); 1944 PetscCall(PetscLogViewWarnDebugging(comm,fd)); 1945 PetscFunctionReturn(0); 1946 } 1947 1948 /*@C 1949 PetscLogView - Prints a summary of the logging. 1950 1951 Collective over MPI_Comm 1952 1953 Input Parameter: 1954 . viewer - an ASCII viewer 1955 1956 Options Database Keys: 1957 + -log_view [:filename] - Prints summary of log information 1958 . -log_view :filename.py:ascii_info_detail - Saves logging information from each process as a Python file 1959 . -log_view :filename.xml:ascii_xml - Saves a summary of the logging information in a nested format (see below for how to view it) 1960 . -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) 1961 . -log_all - Saves a file Log.rank for each MPI process with details of each step of the computation 1962 - -log_trace [filename] - Displays a trace of what each process is doing 1963 1964 Notes: 1965 It is possible to control the logging programatically but we recommend using the options database approach whenever possible 1966 By default the summary is printed to stdout. 1967 1968 Before calling this routine you must have called either PetscLogDefaultBegin() or PetscLogNestedBegin() 1969 1970 If PETSc is configured with --with-logging=0 then this functionality is not available 1971 1972 To view the nested XML format filename.xml first copy ${PETSC_DIR}/share/petsc/xml/performance_xml2html.xsl to the current 1973 directory then open filename.xml with your browser. Specific notes for certain browsers 1974 $ Firefox and Internet explorer - simply open the file 1975 $ Google Chrome - you must start up Chrome with the option --allow-file-access-from-files 1976 $ Safari - see https://ccm.net/faq/36342-safari-how-to-enable-local-file-access 1977 or one can use the package http://xmlsoft.org/XSLT/xsltproc2.html to translate the xml file to html and then open it with 1978 your browser. 1979 Alternatively, use the script ${PETSC_DIR}/lib/petsc/bin/petsc-performance-view to automatically open a new browser 1980 window and render the XML log file contents. 1981 1982 The nested XML format was kindly donated by Koos Huijssen and Christiaan M. Klaij MARITIME RESEARCH INSTITUTE NETHERLANDS 1983 1984 The Flame Graph output can be visualised using either the original Flame Graph script (https://github.com/brendangregg/FlameGraph) 1985 or using speedscope (https://www.speedscope.app). 1986 Old XML profiles may be converted into this format using the script ${PETSC_DIR}/lib/petsc/bin/xml2flamegraph.py. 1987 1988 Level: beginner 1989 1990 .seealso: PetscLogDefaultBegin(), PetscLogDump() 1991 @*/ 1992 PetscErrorCode PetscLogView(PetscViewer viewer) 1993 { 1994 PetscBool isascii; 1995 PetscViewerFormat format; 1996 int stage, lastStage; 1997 PetscStageLog stageLog; 1998 1999 PetscFunctionBegin; 2000 PetscCheck(PetscLogPLB,PETSC_COMM_SELF,PETSC_ERR_SUP,"Must use -log_view or PetscLogDefaultBegin() before calling this routine"); 2001 /* Pop off any stages the user forgot to remove */ 2002 lastStage = 0; 2003 PetscCall(PetscLogGetStageLog(&stageLog)); 2004 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 2005 while (stage >= 0) { 2006 lastStage = stage; 2007 PetscCall(PetscStageLogPop(stageLog)); 2008 PetscCall(PetscStageLogGetCurrent(stageLog, &stage)); 2009 } 2010 PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii)); 2011 PetscCheck(isascii,PetscObjectComm((PetscObject)viewer),PETSC_ERR_SUP,"Currently can only view logging to ASCII"); 2012 PetscCall(PetscViewerGetFormat(viewer,&format)); 2013 if (format == PETSC_VIEWER_DEFAULT || format == PETSC_VIEWER_ASCII_INFO) { 2014 PetscCall(PetscLogView_Default(viewer)); 2015 } else if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) { 2016 PetscCall(PetscLogView_Detailed(viewer)); 2017 } else if (format == PETSC_VIEWER_ASCII_CSV) { 2018 PetscCall(PetscLogView_CSV(viewer)); 2019 } else if (format == PETSC_VIEWER_ASCII_XML) { 2020 PetscCall(PetscLogView_Nested(viewer)); 2021 } else if (format == PETSC_VIEWER_ASCII_FLAMEGRAPH) { 2022 PetscCall(PetscLogView_Flamegraph(viewer)); 2023 } 2024 PetscCall(PetscStageLogPush(stageLog, lastStage)); 2025 PetscFunctionReturn(0); 2026 } 2027 2028 /*@C 2029 PetscLogViewFromOptions - Processes command line options to determine if/how a PetscLog is to be viewed. 2030 2031 Collective on PETSC_COMM_WORLD 2032 2033 Not normally called by user 2034 2035 Level: intermediate 2036 2037 @*/ 2038 PetscErrorCode PetscLogViewFromOptions(void) 2039 { 2040 PetscViewer viewer; 2041 PetscBool flg; 2042 PetscViewerFormat format; 2043 2044 PetscFunctionBegin; 2045 PetscCall(PetscOptionsGetViewer(PETSC_COMM_WORLD,NULL,NULL,"-log_view",&viewer,&format,&flg)); 2046 if (flg) { 2047 PetscCall(PetscViewerPushFormat(viewer,format)); 2048 PetscCall(PetscLogView(viewer)); 2049 PetscCall(PetscViewerPopFormat(viewer)); 2050 PetscCall(PetscViewerDestroy(&viewer)); 2051 } 2052 PetscFunctionReturn(0); 2053 } 2054 2055 /*----------------------------------------------- Counter Functions -------------------------------------------------*/ 2056 /*@C 2057 PetscGetFlops - Returns the number of flops used on this processor 2058 since the program began. 2059 2060 Not Collective 2061 2062 Output Parameter: 2063 flops - number of floating point operations 2064 2065 Notes: 2066 A global counter logs all PETSc flop counts. The user can use 2067 PetscLogFlops() to increment this counter to include flops for the 2068 application code. 2069 2070 Level: intermediate 2071 2072 .seealso: PetscTime(), PetscLogFlops() 2073 @*/ 2074 PetscErrorCode PetscGetFlops(PetscLogDouble *flops) 2075 { 2076 PetscFunctionBegin; 2077 *flops = petsc_TotalFlops; 2078 PetscFunctionReturn(0); 2079 } 2080 2081 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...) 2082 { 2083 size_t fullLength; 2084 va_list Argp; 2085 2086 PetscFunctionBegin; 2087 if (!petsc_logObjects) PetscFunctionReturn(0); 2088 va_start(Argp, format); 2089 PetscCall(PetscVSNPrintf(petsc_objects[obj->id].info, 64,format,&fullLength, Argp)); 2090 va_end(Argp); 2091 PetscFunctionReturn(0); 2092 } 2093 2094 /*MC 2095 PetscLogFlops - Adds floating point operations to the global counter. 2096 2097 Synopsis: 2098 #include <petsclog.h> 2099 PetscErrorCode PetscLogFlops(PetscLogDouble f) 2100 2101 Not Collective 2102 2103 Input Parameter: 2104 . f - flop counter 2105 2106 Usage: 2107 .vb 2108 PetscLogEvent USER_EVENT; 2109 PetscLogEventRegister("User event",0,&USER_EVENT); 2110 PetscLogEventBegin(USER_EVENT,0,0,0,0); 2111 [code segment to monitor] 2112 PetscLogFlops(user_flops) 2113 PetscLogEventEnd(USER_EVENT,0,0,0,0); 2114 .ve 2115 2116 Notes: 2117 A global counter logs all PETSc flop counts. The user can use 2118 PetscLogFlops() to increment this counter to include flops for the 2119 application code. 2120 2121 Level: intermediate 2122 2123 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscGetFlops() 2124 2125 M*/ 2126 2127 /*MC 2128 PetscPreLoadBegin - Begin a segment of code that may be preloaded (run twice) 2129 to get accurate timings 2130 2131 Synopsis: 2132 #include <petsclog.h> 2133 void PetscPreLoadBegin(PetscBool flag,char *name); 2134 2135 Not Collective 2136 2137 Input Parameters: 2138 + flag - PETSC_TRUE to run twice, PETSC_FALSE to run once, may be overridden 2139 with command line option -preload true or -preload false 2140 - name - name of first stage (lines of code timed separately with -log_view) to 2141 be preloaded 2142 2143 Usage: 2144 .vb 2145 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2146 lines of code 2147 PetscPreLoadStage("second stage"); 2148 lines of code 2149 PetscPreLoadEnd(); 2150 .ve 2151 2152 Notes: 2153 Only works in C/C++, not Fortran 2154 2155 Flags available within the macro. 2156 + PetscPreLoadingUsed - true if we are or have done preloading 2157 . PetscPreLoadingOn - true if it is CURRENTLY doing preload 2158 . PetscPreLoadIt - 0 for the first computation (with preloading turned off it is only 0) 1 for the second 2159 - PetscPreLoadMax - number of times it will do the computation, only one when preloading is turned on 2160 The first two variables are available throughout the program, the second two only between the PetscPreLoadBegin() 2161 and PetscPreLoadEnd() 2162 2163 Level: intermediate 2164 2165 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadEnd(), PetscPreLoadStage() 2166 2167 M*/ 2168 2169 /*MC 2170 PetscPreLoadEnd - End a segment of code that may be preloaded (run twice) 2171 to get accurate timings 2172 2173 Synopsis: 2174 #include <petsclog.h> 2175 void PetscPreLoadEnd(void); 2176 2177 Not Collective 2178 2179 Usage: 2180 .vb 2181 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2182 lines of code 2183 PetscPreLoadStage("second stage"); 2184 lines of code 2185 PetscPreLoadEnd(); 2186 .ve 2187 2188 Notes: 2189 only works in C/C++ not fortran 2190 2191 Level: intermediate 2192 2193 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadStage() 2194 2195 M*/ 2196 2197 /*MC 2198 PetscPreLoadStage - Start a new segment of code to be timed separately. 2199 to get accurate timings 2200 2201 Synopsis: 2202 #include <petsclog.h> 2203 void PetscPreLoadStage(char *name); 2204 2205 Not Collective 2206 2207 Usage: 2208 .vb 2209 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2210 lines of code 2211 PetscPreLoadStage("second stage"); 2212 lines of code 2213 PetscPreLoadEnd(); 2214 .ve 2215 2216 Notes: 2217 only works in C/C++ not fortran 2218 2219 Level: intermediate 2220 2221 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd() 2222 2223 M*/ 2224 2225 #if PetscDefined(HAVE_DEVICE) 2226 #include <petsc/private/deviceimpl.h> 2227 2228 /*-------------------------------------------- GPU event Functions ----------------------------------------------*/ 2229 /*@C 2230 PetscLogGpuTimeBegin - Start timer for device 2231 2232 Notes: 2233 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). 2234 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). 2235 There is no need to call WaitForCUDA() or WaitForHIP() between PetscLogGpuTimeBegin and PetscLogGpuTimeEnd 2236 This timer should NOT include times for data transfers between the GPU and CPU, nor setup actions such as allocating space. 2237 The regular logging captures the time for data transfers and any CPU activites during the event 2238 It is used to compute the flop rate on the GPU as it is actively engaged in running a kernel. 2239 2240 Developer Notes: 2241 The GPU event timer captures the execution time of all the kernels launched in the default stream by the CPU between PetscLogGpuTimeBegin() and PetsLogGpuTimeEnd(). 2242 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. 2243 2244 Level: intermediate 2245 2246 .seealso: PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeEnd() 2247 @*/ 2248 PetscErrorCode PetscLogGpuTimeBegin(void) 2249 { 2250 PetscFunctionBegin; 2251 if (!PetscLogPLB) PetscFunctionReturn(0); 2252 if (PetscDefined(HAVE_CUDA) || PetscDefined(HAVE_HIP)) { 2253 PetscDeviceContext dctx; 2254 2255 PetscCall(PetscDeviceContextGetCurrentContext(&dctx)); 2256 PetscCall(PetscDeviceContextBeginTimer_Internal(dctx)); 2257 } else { 2258 PetscCall(PetscTimeSubtract(&petsc_gtime)); 2259 } 2260 PetscFunctionReturn(0); 2261 } 2262 2263 /*@C 2264 PetscLogGpuTimeEnd - Stop timer for device 2265 2266 Level: intermediate 2267 2268 .seealso: PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeBegin() 2269 @*/ 2270 PetscErrorCode PetscLogGpuTimeEnd(void) 2271 { 2272 PetscFunctionBegin; 2273 if (!PetscLogPLE) PetscFunctionReturn(0); 2274 if (PetscDefined(HAVE_CUDA) || PetscDefined(HAVE_HIP)) { 2275 PetscDeviceContext dctx; 2276 PetscLogDouble elapsed; 2277 2278 PetscCall(PetscDeviceContextGetCurrentContext(&dctx)); 2279 PetscCall(PetscDeviceContextEndTimer_Internal(dctx,&elapsed)); 2280 petsc_gtime += (elapsed/1000.0); 2281 } else { 2282 PetscCall(PetscTimeAdd(&petsc_gtime)); 2283 } 2284 PetscFunctionReturn(0); 2285 } 2286 #endif /* end of PETSC_HAVE_DEVICE */ 2287 2288 #else /* end of -DPETSC_USE_LOG section */ 2289 2290 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...) 2291 { 2292 PetscFunctionBegin; 2293 PetscFunctionReturn(0); 2294 } 2295 2296 #endif /* PETSC_USE_LOG*/ 2297 2298 PetscClassId PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID; 2299 PetscClassId PETSC_OBJECT_CLASSID = 0; 2300 2301 /*@C 2302 PetscClassIdRegister - Registers a new class name for objects and logging operations in an application code. 2303 2304 Not Collective 2305 2306 Input Parameter: 2307 . name - The class name 2308 2309 Output Parameter: 2310 . oclass - The class id or classid 2311 2312 Level: developer 2313 2314 @*/ 2315 PetscErrorCode PetscClassIdRegister(const char name[],PetscClassId *oclass) 2316 { 2317 #if defined(PETSC_USE_LOG) 2318 PetscStageLog stageLog; 2319 PetscInt stage; 2320 #endif 2321 2322 PetscFunctionBegin; 2323 *oclass = ++PETSC_LARGEST_CLASSID; 2324 #if defined(PETSC_USE_LOG) 2325 PetscCall(PetscLogGetStageLog(&stageLog)); 2326 PetscCall(PetscClassRegLogRegister(stageLog->classLog, name, *oclass)); 2327 for (stage = 0; stage < stageLog->numStages; stage++) { 2328 PetscCall(PetscClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses)); 2329 } 2330 #endif 2331 PetscFunctionReturn(0); 2332 } 2333 2334 #if defined(PETSC_USE_LOG) && defined(PETSC_HAVE_MPE) 2335 #include <mpe.h> 2336 2337 PetscBool PetscBeganMPE = PETSC_FALSE; 2338 2339 PETSC_INTERN PetscErrorCode PetscLogEventBeginMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject); 2340 PETSC_INTERN PetscErrorCode PetscLogEventEndMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject); 2341 2342 /*@C 2343 PetscLogMPEBegin - Turns on MPE logging of events. This creates large log files 2344 and slows the program down. 2345 2346 Collective over PETSC_COMM_WORLD 2347 2348 Options Database Keys: 2349 . -log_mpe - Prints extensive log information 2350 2351 Notes: 2352 A related routine is PetscLogDefaultBegin() (with the options key -log_view), which is 2353 intended for production runs since it logs only flop rates and object 2354 creation (and should not significantly slow the programs). 2355 2356 Level: advanced 2357 2358 .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogEventActivate(), 2359 PetscLogEventDeactivate() 2360 @*/ 2361 PetscErrorCode PetscLogMPEBegin(void) 2362 { 2363 PetscFunctionBegin; 2364 /* Do MPE initialization */ 2365 if (!MPE_Initialized_logging()) { /* This function exists in mpich 1.1.2 and higher */ 2366 PetscCall(PetscInfo(0,"Initializing MPE.\n")); 2367 PetscCall(MPE_Init_log()); 2368 2369 PetscBeganMPE = PETSC_TRUE; 2370 } else { 2371 PetscCall(PetscInfo(0,"MPE already initialized. Not attempting to reinitialize.\n")); 2372 } 2373 PetscCall(PetscLogSet(PetscLogEventBeginMPE, PetscLogEventEndMPE)); 2374 PetscFunctionReturn(0); 2375 } 2376 2377 /*@C 2378 PetscLogMPEDump - Dumps the MPE logging info to file for later use with Jumpshot. 2379 2380 Collective over PETSC_COMM_WORLD 2381 2382 Level: advanced 2383 2384 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogMPEBegin() 2385 @*/ 2386 PetscErrorCode PetscLogMPEDump(const char sname[]) 2387 { 2388 char name[PETSC_MAX_PATH_LEN]; 2389 2390 PetscFunctionBegin; 2391 if (PetscBeganMPE) { 2392 PetscCall(PetscInfo(0,"Finalizing MPE.\n")); 2393 if (sname) { 2394 PetscCall(PetscStrcpy(name,sname)); 2395 } else { 2396 PetscCall(PetscGetProgramName(name,sizeof(name))); 2397 } 2398 PetscCall(MPE_Finish_log(name)); 2399 } else { 2400 PetscCall(PetscInfo(0,"Not finalizing MPE (not started by PETSc).\n")); 2401 } 2402 PetscFunctionReturn(0); 2403 } 2404 2405 #define PETSC_RGB_COLORS_MAX 39 2406 static const char *PetscLogMPERGBColors[PETSC_RGB_COLORS_MAX] = { 2407 "OliveDrab: ", 2408 "BlueViolet: ", 2409 "CadetBlue: ", 2410 "CornflowerBlue: ", 2411 "DarkGoldenrod: ", 2412 "DarkGreen: ", 2413 "DarkKhaki: ", 2414 "DarkOliveGreen: ", 2415 "DarkOrange: ", 2416 "DarkOrchid: ", 2417 "DarkSeaGreen: ", 2418 "DarkSlateGray: ", 2419 "DarkTurquoise: ", 2420 "DeepPink: ", 2421 "DarkKhaki: ", 2422 "DimGray: ", 2423 "DodgerBlue: ", 2424 "GreenYellow: ", 2425 "HotPink: ", 2426 "IndianRed: ", 2427 "LavenderBlush: ", 2428 "LawnGreen: ", 2429 "LemonChiffon: ", 2430 "LightCoral: ", 2431 "LightCyan: ", 2432 "LightPink: ", 2433 "LightSalmon: ", 2434 "LightSlateGray: ", 2435 "LightYellow: ", 2436 "LimeGreen: ", 2437 "MediumPurple: ", 2438 "MediumSeaGreen: ", 2439 "MediumSlateBlue:", 2440 "MidnightBlue: ", 2441 "MintCream: ", 2442 "MistyRose: ", 2443 "NavajoWhite: ", 2444 "NavyBlue: ", 2445 "OliveDrab: " 2446 }; 2447 2448 /*@C 2449 PetscLogMPEGetRGBColor - This routine returns a rgb color useable with PetscLogEventRegister() 2450 2451 Not collective. Maybe it should be? 2452 2453 Output Parameter: 2454 . str - character string representing the color 2455 2456 Level: developer 2457 2458 .seealso: PetscLogEventRegister 2459 @*/ 2460 PetscErrorCode PetscLogMPEGetRGBColor(const char *str[]) 2461 { 2462 static int idx = 0; 2463 2464 PetscFunctionBegin; 2465 *str = PetscLogMPERGBColors[idx]; 2466 idx = (idx + 1)% PETSC_RGB_COLORS_MAX; 2467 PetscFunctionReturn(0); 2468 } 2469 2470 #endif /* PETSC_USE_LOG && PETSC_HAVE_MPE */ 2471