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