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