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