1 2 /* 3 PETSc code to log object creation and destruction and PETSc events. 4 5 This provides the public API used by the rest of PETSc and by users. 6 7 These routines use a private API that is not used elsewhere in PETSc and is not 8 accessible to users. The private API is defined in logimpl.h and the utils directory. 9 10 */ 11 #include <petsc/private/logimpl.h> /*I "petscsys.h" I*/ 12 #include <petsctime.h> 13 #include <petscviewer.h> 14 15 PetscErrorCode PetscLogObjectParent(PetscObject p,PetscObject c) 16 { 17 if (!c || !p) return 0; 18 c->parent = p; 19 c->parentid = p->id; 20 return 0; 21 } 22 23 /*@C 24 PetscLogObjectMemory - Adds to an object a count of additional amount of memory that is used by the object. 25 26 Not collective. 27 28 Input Parameters: 29 + obj - the PETSc object 30 - mem - the amount of memory that is being added to the object 31 32 Level: developer 33 34 Developer Notes: 35 Currently we do not always do a good job of associating all memory allocations with an object. 36 37 .seealso: PetscFinalize(), PetscInitializeFortran(), PetscGetArgs(), PetscInitializeNoArguments() 38 39 @*/ 40 PetscErrorCode PetscLogObjectMemory(PetscObject p,PetscLogDouble m) 41 { 42 if (!p) return 0; 43 p->mem += m; 44 return 0; 45 } 46 47 PetscLogEvent PETSC_LARGEST_EVENT = PETSC_EVENT; 48 49 #if defined(PETSC_USE_LOG) 50 #include <petscmachineinfo.h> 51 #include <petscconfiginfo.h> 52 53 /* used in the MPI_XXX() count macros in petsclog.h */ 54 55 /* Action and object logging variables */ 56 Action *petsc_actions = NULL; 57 Object *petsc_objects = NULL; 58 PetscBool petsc_logActions = PETSC_FALSE; 59 PetscBool petsc_logObjects = PETSC_FALSE; 60 int petsc_numActions = 0, petsc_maxActions = 100; 61 int petsc_numObjects = 0, petsc_maxObjects = 100; 62 int petsc_numObjectsDestroyed = 0; 63 64 /* Global counters */ 65 PetscLogDouble petsc_BaseTime = 0.0; 66 PetscLogDouble petsc_TotalFlops = 0.0; /* The number of flops */ 67 PetscLogDouble petsc_tmp_flops = 0.0; /* The incremental number of flops */ 68 PetscLogDouble petsc_send_ct = 0.0; /* The number of sends */ 69 PetscLogDouble petsc_recv_ct = 0.0; /* The number of receives */ 70 PetscLogDouble petsc_send_len = 0.0; /* The total length of all sent messages */ 71 PetscLogDouble petsc_recv_len = 0.0; /* The total length of all received messages */ 72 PetscLogDouble petsc_isend_ct = 0.0; /* The number of immediate sends */ 73 PetscLogDouble petsc_irecv_ct = 0.0; /* The number of immediate receives */ 74 PetscLogDouble petsc_isend_len = 0.0; /* The total length of all immediate send messages */ 75 PetscLogDouble petsc_irecv_len = 0.0; /* The total length of all immediate receive messages */ 76 PetscLogDouble petsc_wait_ct = 0.0; /* The number of waits */ 77 PetscLogDouble petsc_wait_any_ct = 0.0; /* The number of anywaits */ 78 PetscLogDouble petsc_wait_all_ct = 0.0; /* The number of waitalls */ 79 PetscLogDouble petsc_sum_of_waits_ct = 0.0; /* The total number of waits */ 80 PetscLogDouble petsc_allreduce_ct = 0.0; /* The number of reductions */ 81 PetscLogDouble petsc_gather_ct = 0.0; /* The number of gathers and gathervs */ 82 PetscLogDouble petsc_scatter_ct = 0.0; /* The number of scatters and scattervs */ 83 #if defined(PETSC_HAVE_DEVICE) 84 PetscLogDouble petsc_ctog_ct = 0.0; /* The total number of CPU to GPU copies */ 85 PetscLogDouble petsc_gtoc_ct = 0.0; /* The total number of GPU to CPU copies */ 86 PetscLogDouble petsc_ctog_sz = 0.0; /* The total size of CPU to GPU copies */ 87 PetscLogDouble petsc_gtoc_sz = 0.0; /* The total size of GPU to CPU copies */ 88 PetscLogDouble petsc_ctog_ct_scalar = 0.0; /* The total number of CPU to GPU copies */ 89 PetscLogDouble petsc_gtoc_ct_scalar = 0.0; /* The total number of GPU to CPU copies */ 90 PetscLogDouble petsc_ctog_sz_scalar = 0.0; /* The total size of CPU to GPU copies */ 91 PetscLogDouble petsc_gtoc_sz_scalar = 0.0; /* The total size of GPU to CPU copies */ 92 PetscLogDouble petsc_gflops = 0.0; /* The flops done on a GPU */ 93 PetscLogDouble petsc_gtime = 0.0; /* The time spent on a GPU */ 94 95 #if defined(PETSC_USE_DEBUG) 96 PetscBool petsc_gtime_inuse = PETSC_FALSE; 97 #endif 98 #endif 99 100 /* Logging functions */ 101 PetscErrorCode (*PetscLogPHC)(PetscObject) = NULL; 102 PetscErrorCode (*PetscLogPHD)(PetscObject) = NULL; 103 PetscErrorCode (*PetscLogPLB)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL; 104 PetscErrorCode (*PetscLogPLE)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL; 105 106 /* Tracing event logging variables */ 107 FILE *petsc_tracefile = NULL; 108 int petsc_tracelevel = 0; 109 const char *petsc_traceblanks = " "; 110 char petsc_tracespace[128] = " "; 111 PetscLogDouble petsc_tracetime = 0.0; 112 static PetscBool PetscLogInitializeCalled = PETSC_FALSE; 113 114 PETSC_INTERN PetscErrorCode PetscLogInitialize(void) 115 { 116 int stage; 117 PetscBool opt; 118 PetscErrorCode ierr; 119 120 PetscFunctionBegin; 121 if (PetscLogInitializeCalled) PetscFunctionReturn(0); 122 PetscLogInitializeCalled = PETSC_TRUE; 123 124 ierr = PetscOptionsHasName(NULL,NULL, "-log_exclude_actions", &opt);CHKERRQ(ierr); 125 if (opt) petsc_logActions = PETSC_FALSE; 126 ierr = PetscOptionsHasName(NULL,NULL, "-log_exclude_objects", &opt);CHKERRQ(ierr); 127 if (opt) petsc_logObjects = PETSC_FALSE; 128 if (petsc_logActions) { 129 ierr = PetscMalloc1(petsc_maxActions, &petsc_actions);CHKERRQ(ierr); 130 } 131 if (petsc_logObjects) { 132 ierr = PetscMalloc1(petsc_maxObjects, &petsc_objects);CHKERRQ(ierr); 133 } 134 PetscLogPHC = PetscLogObjCreateDefault; 135 PetscLogPHD = PetscLogObjDestroyDefault; 136 /* Setup default logging structures */ 137 ierr = PetscStageLogCreate(&petsc_stageLog);CHKERRQ(ierr); 138 ierr = PetscStageLogRegister(petsc_stageLog, "Main Stage", &stage);CHKERRQ(ierr); 139 140 /* All processors sync here for more consistent logging */ 141 ierr = MPI_Barrier(PETSC_COMM_WORLD);CHKERRMPI(ierr); 142 PetscTime(&petsc_BaseTime); 143 ierr = PetscLogStagePush(stage);CHKERRQ(ierr); 144 PetscFunctionReturn(0); 145 } 146 147 PETSC_INTERN PetscErrorCode PetscLogFinalize(void) 148 { 149 PetscStageLog stageLog; 150 PetscErrorCode ierr; 151 152 PetscFunctionBegin; 153 ierr = PetscFree(petsc_actions);CHKERRQ(ierr); 154 ierr = PetscFree(petsc_objects);CHKERRQ(ierr); 155 ierr = PetscLogNestedEnd();CHKERRQ(ierr); 156 ierr = PetscLogSet(NULL, NULL);CHKERRQ(ierr); 157 158 /* Resetting phase */ 159 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 160 ierr = PetscStageLogDestroy(stageLog);CHKERRQ(ierr); 161 162 petsc_TotalFlops = 0.0; 163 petsc_numActions = 0; 164 petsc_numObjects = 0; 165 petsc_numObjectsDestroyed = 0; 166 petsc_maxActions = 100; 167 petsc_maxObjects = 100; 168 petsc_actions = NULL; 169 petsc_objects = NULL; 170 petsc_logActions = PETSC_FALSE; 171 petsc_logObjects = PETSC_FALSE; 172 petsc_BaseTime = 0.0; 173 petsc_TotalFlops = 0.0; 174 petsc_tmp_flops = 0.0; 175 petsc_send_ct = 0.0; 176 petsc_recv_ct = 0.0; 177 petsc_send_len = 0.0; 178 petsc_recv_len = 0.0; 179 petsc_isend_ct = 0.0; 180 petsc_irecv_ct = 0.0; 181 petsc_isend_len = 0.0; 182 petsc_irecv_len = 0.0; 183 petsc_wait_ct = 0.0; 184 petsc_wait_any_ct = 0.0; 185 petsc_wait_all_ct = 0.0; 186 petsc_sum_of_waits_ct = 0.0; 187 petsc_allreduce_ct = 0.0; 188 petsc_gather_ct = 0.0; 189 petsc_scatter_ct = 0.0; 190 #if defined(PETSC_HAVE_DEVICE) 191 petsc_ctog_ct = 0.0; 192 petsc_gtoc_ct = 0.0; 193 petsc_ctog_sz = 0.0; 194 petsc_gtoc_sz = 0.0; 195 petsc_gflops = 0.0; 196 petsc_gtime = 0.0; 197 #endif 198 PETSC_LARGEST_EVENT = PETSC_EVENT; 199 PetscLogPHC = NULL; 200 PetscLogPHD = NULL; 201 petsc_tracefile = NULL; 202 petsc_tracelevel = 0; 203 petsc_traceblanks = " "; 204 petsc_tracespace[0] = ' '; petsc_tracespace[1] = 0; 205 petsc_tracetime = 0.0; 206 PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID; 207 PETSC_OBJECT_CLASSID = 0; 208 petsc_stageLog = NULL; 209 PetscLogInitializeCalled = PETSC_FALSE; 210 PetscFunctionReturn(0); 211 } 212 213 /*@C 214 PetscLogSet - Sets the logging functions called at the beginning and ending of every event. 215 216 Not Collective 217 218 Input Parameters: 219 + b - The function called at beginning of event 220 - e - The function called at end of event 221 222 Level: developer 223 224 .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogTraceBegin() 225 @*/ 226 PetscErrorCode PetscLogSet(PetscErrorCode (*b)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject), 227 PetscErrorCode (*e)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject)) 228 { 229 PetscFunctionBegin; 230 PetscLogPLB = b; 231 PetscLogPLE = e; 232 PetscFunctionReturn(0); 233 } 234 235 /*@C 236 PetscLogDefaultBegin - Turns on logging of objects and events. This logs flop 237 rates and object creation and should not slow programs down too much. 238 This routine may be called more than once. 239 240 Logically Collective over PETSC_COMM_WORLD 241 242 Options Database Keys: 243 . -log_view [viewertype:filename:viewerformat] - Prints summary of flop and timing information to the 244 screen (for code configured with --with-log=1 (which is the default)) 245 246 Usage: 247 .vb 248 PetscInitialize(...); 249 PetscLogDefaultBegin(); 250 ... code ... 251 PetscLogView(viewer); or PetscLogDump(); 252 PetscFinalize(); 253 .ve 254 255 Notes: 256 PetscLogView(viewer) or PetscLogDump() actually cause the printing of 257 the logging information. 258 259 Level: advanced 260 261 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogTraceBegin() 262 @*/ 263 PetscErrorCode PetscLogDefaultBegin(void) 264 { 265 PetscErrorCode ierr; 266 267 PetscFunctionBegin; 268 ierr = PetscLogSet(PetscLogEventBeginDefault, PetscLogEventEndDefault);CHKERRQ(ierr); 269 PetscFunctionReturn(0); 270 } 271 272 /*@C 273 PetscLogAllBegin - Turns on extensive logging of objects and events. Logs 274 all events. This creates large log files and slows the program down. 275 276 Logically Collective on PETSC_COMM_WORLD 277 278 Options Database Keys: 279 . -log_all - Prints extensive log information 280 281 Usage: 282 .vb 283 PetscInitialize(...); 284 PetscLogAllBegin(); 285 ... code ... 286 PetscLogDump(filename); 287 PetscFinalize(); 288 .ve 289 290 Notes: 291 A related routine is PetscLogDefaultBegin() (with the options key -log), which is 292 intended for production runs since it logs only flop rates and object 293 creation (and shouldn't significantly slow the programs). 294 295 Level: advanced 296 297 .seealso: PetscLogDump(), PetscLogDefaultBegin(), PetscLogTraceBegin() 298 @*/ 299 PetscErrorCode PetscLogAllBegin(void) 300 { 301 PetscErrorCode ierr; 302 303 PetscFunctionBegin; 304 ierr = PetscLogSet(PetscLogEventBeginComplete, PetscLogEventEndComplete);CHKERRQ(ierr); 305 PetscFunctionReturn(0); 306 } 307 308 /*@C 309 PetscLogTraceBegin - Activates trace logging. Every time a PETSc event 310 begins or ends, the event name is printed. 311 312 Logically Collective on PETSC_COMM_WORLD 313 314 Input Parameter: 315 . file - The file to print trace in (e.g. stdout) 316 317 Options Database Key: 318 . -log_trace [filename] - Activates PetscLogTraceBegin() 319 320 Notes: 321 PetscLogTraceBegin() prints the processor number, the execution time (sec), 322 then "Event begin:" or "Event end:" followed by the event name. 323 324 PetscLogTraceBegin() allows tracing of all PETSc calls, which is useful 325 to determine where a program is hanging without running in the 326 debugger. Can be used in conjunction with the -info option. 327 328 Level: intermediate 329 330 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogDefaultBegin() 331 @*/ 332 PetscErrorCode PetscLogTraceBegin(FILE *file) 333 { 334 PetscErrorCode ierr; 335 336 PetscFunctionBegin; 337 petsc_tracefile = file; 338 339 ierr = PetscLogSet(PetscLogEventBeginTrace, PetscLogEventEndTrace);CHKERRQ(ierr); 340 PetscFunctionReturn(0); 341 } 342 343 /*@ 344 PetscLogActions - Determines whether actions are logged for the graphical viewer. 345 346 Not Collective 347 348 Input Parameter: 349 . flag - PETSC_TRUE if actions are to be logged 350 351 Level: intermediate 352 353 Note: Logging of actions continues to consume more memory as the program 354 runs. Long running programs should consider turning this feature off. 355 356 Options Database Keys: 357 . -log_exclude_actions - Turns off actions logging 358 359 .seealso: PetscLogStagePush(), PetscLogStagePop() 360 @*/ 361 PetscErrorCode PetscLogActions(PetscBool flag) 362 { 363 PetscFunctionBegin; 364 petsc_logActions = flag; 365 PetscFunctionReturn(0); 366 } 367 368 /*@ 369 PetscLogObjects - Determines whether objects are logged for the graphical viewer. 370 371 Not Collective 372 373 Input Parameter: 374 . flag - PETSC_TRUE if objects are to be logged 375 376 Level: intermediate 377 378 Note: Logging of objects continues to consume more memory as the program 379 runs. Long running programs should consider turning this feature off. 380 381 Options Database Keys: 382 . -log_exclude_objects - Turns off objects logging 383 384 .seealso: PetscLogStagePush(), PetscLogStagePop() 385 @*/ 386 PetscErrorCode PetscLogObjects(PetscBool flag) 387 { 388 PetscFunctionBegin; 389 petsc_logObjects = flag; 390 PetscFunctionReturn(0); 391 } 392 393 /*------------------------------------------------ Stage Functions --------------------------------------------------*/ 394 /*@C 395 PetscLogStageRegister - Attaches a character string name to a logging stage. 396 397 Not Collective 398 399 Input Parameter: 400 . sname - The name to associate with that stage 401 402 Output Parameter: 403 . stage - The stage number 404 405 Level: intermediate 406 407 .seealso: PetscLogStagePush(), PetscLogStagePop() 408 @*/ 409 PetscErrorCode PetscLogStageRegister(const char sname[],PetscLogStage *stage) 410 { 411 PetscStageLog stageLog; 412 PetscLogEvent event; 413 PetscErrorCode ierr; 414 415 PetscFunctionBegin; 416 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 417 ierr = PetscStageLogRegister(stageLog, sname, stage);CHKERRQ(ierr); 418 /* Copy events already changed in the main stage, this sucks */ 419 ierr = PetscEventPerfLogEnsureSize(stageLog->stageInfo[*stage].eventLog, stageLog->eventLog->numEvents);CHKERRQ(ierr); 420 for (event = 0; event < stageLog->eventLog->numEvents; event++) { 421 ierr = PetscEventPerfInfoCopy(&stageLog->stageInfo[0].eventLog->eventInfo[event],&stageLog->stageInfo[*stage].eventLog->eventInfo[event]);CHKERRQ(ierr); 422 } 423 ierr = PetscClassPerfLogEnsureSize(stageLog->stageInfo[*stage].classLog, stageLog->classLog->numClasses);CHKERRQ(ierr); 424 PetscFunctionReturn(0); 425 } 426 427 /*@C 428 PetscLogStagePush - This function pushes a stage on the stack. 429 430 Not Collective 431 432 Input Parameter: 433 . stage - The stage on which to log 434 435 Usage: 436 If the option -log_sumary is used to run the program containing the 437 following code, then 2 sets of summary data will be printed during 438 PetscFinalize(). 439 .vb 440 PetscInitialize(int *argc,char ***args,0,0); 441 [stage 0 of code] 442 PetscLogStagePush(1); 443 [stage 1 of code] 444 PetscLogStagePop(); 445 PetscBarrier(...); 446 [more stage 0 of code] 447 PetscFinalize(); 448 .ve 449 450 Notes: 451 Use PetscLogStageRegister() to register a stage. 452 453 Level: intermediate 454 455 .seealso: PetscLogStagePop(), PetscLogStageRegister(), PetscBarrier() 456 @*/ 457 PetscErrorCode PetscLogStagePush(PetscLogStage stage) 458 { 459 PetscStageLog stageLog; 460 PetscErrorCode ierr; 461 462 PetscFunctionBegin; 463 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 464 ierr = PetscStageLogPush(stageLog, stage);CHKERRQ(ierr); 465 PetscFunctionReturn(0); 466 } 467 468 /*@C 469 PetscLogStagePop - This function pops a stage from the stack. 470 471 Not Collective 472 473 Usage: 474 If the option -log_sumary is used to run the program containing the 475 following code, then 2 sets of summary data will be printed during 476 PetscFinalize(). 477 .vb 478 PetscInitialize(int *argc,char ***args,0,0); 479 [stage 0 of code] 480 PetscLogStagePush(1); 481 [stage 1 of code] 482 PetscLogStagePop(); 483 PetscBarrier(...); 484 [more stage 0 of code] 485 PetscFinalize(); 486 .ve 487 488 Notes: 489 Use PetscLogStageRegister() to register a stage. 490 491 Level: intermediate 492 493 .seealso: PetscLogStagePush(), PetscLogStageRegister(), PetscBarrier() 494 @*/ 495 PetscErrorCode PetscLogStagePop(void) 496 { 497 PetscStageLog stageLog; 498 PetscErrorCode ierr; 499 500 PetscFunctionBegin; 501 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 502 ierr = PetscStageLogPop(stageLog);CHKERRQ(ierr); 503 PetscFunctionReturn(0); 504 } 505 506 /*@ 507 PetscLogStageSetActive - Determines stage activity for PetscLogEventBegin() and PetscLogEventEnd(). 508 509 Not Collective 510 511 Input Parameters: 512 + stage - The stage 513 - isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE) 514 515 Level: intermediate 516 517 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 518 @*/ 519 PetscErrorCode PetscLogStageSetActive(PetscLogStage stage, PetscBool isActive) 520 { 521 PetscStageLog stageLog; 522 PetscErrorCode ierr; 523 524 PetscFunctionBegin; 525 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 526 ierr = PetscStageLogSetActive(stageLog, stage, isActive);CHKERRQ(ierr); 527 PetscFunctionReturn(0); 528 } 529 530 /*@ 531 PetscLogStageGetActive - Returns stage activity for PetscLogEventBegin() and PetscLogEventEnd(). 532 533 Not Collective 534 535 Input Parameter: 536 . stage - The stage 537 538 Output Parameter: 539 . isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE) 540 541 Level: intermediate 542 543 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 544 @*/ 545 PetscErrorCode PetscLogStageGetActive(PetscLogStage stage, PetscBool *isActive) 546 { 547 PetscStageLog stageLog; 548 PetscErrorCode ierr; 549 550 PetscFunctionBegin; 551 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 552 ierr = PetscStageLogGetActive(stageLog, stage, isActive);CHKERRQ(ierr); 553 PetscFunctionReturn(0); 554 } 555 556 /*@ 557 PetscLogStageSetVisible - Determines stage visibility in PetscLogView() 558 559 Not Collective 560 561 Input Parameters: 562 + stage - The stage 563 - isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE) 564 565 Level: intermediate 566 567 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView() 568 @*/ 569 PetscErrorCode PetscLogStageSetVisible(PetscLogStage stage, PetscBool isVisible) 570 { 571 PetscStageLog stageLog; 572 PetscErrorCode ierr; 573 574 PetscFunctionBegin; 575 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 576 ierr = PetscStageLogSetVisible(stageLog, stage, isVisible);CHKERRQ(ierr); 577 PetscFunctionReturn(0); 578 } 579 580 /*@ 581 PetscLogStageGetVisible - Returns stage visibility in PetscLogView() 582 583 Not Collective 584 585 Input Parameter: 586 . stage - The stage 587 588 Output Parameter: 589 . isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE) 590 591 Level: intermediate 592 593 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView() 594 @*/ 595 PetscErrorCode PetscLogStageGetVisible(PetscLogStage stage, PetscBool *isVisible) 596 { 597 PetscStageLog stageLog; 598 PetscErrorCode ierr; 599 600 PetscFunctionBegin; 601 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 602 ierr = PetscStageLogGetVisible(stageLog, stage, isVisible);CHKERRQ(ierr); 603 PetscFunctionReturn(0); 604 } 605 606 /*@C 607 PetscLogStageGetId - Returns the stage id when given the stage name. 608 609 Not Collective 610 611 Input Parameter: 612 . name - The stage name 613 614 Output Parameter: 615 . stage - The stage, , or -1 if no stage with that name exists 616 617 Level: intermediate 618 619 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 620 @*/ 621 PetscErrorCode PetscLogStageGetId(const char name[], PetscLogStage *stage) 622 { 623 PetscStageLog stageLog; 624 PetscErrorCode ierr; 625 626 PetscFunctionBegin; 627 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 628 ierr = PetscStageLogGetStage(stageLog, name, stage);CHKERRQ(ierr); 629 PetscFunctionReturn(0); 630 } 631 632 /*------------------------------------------------ Event Functions --------------------------------------------------*/ 633 /*@C 634 PetscLogEventRegister - Registers an event name for logging operations in an application code. 635 636 Not Collective 637 638 Input Parameters: 639 + name - The name associated with the event 640 - classid - The classid associated to the class for this event, obtain either with 641 PetscClassIdRegister() or use a predefined one such as KSP_CLASSID, SNES_CLASSID, the predefined ones 642 are only available in C code 643 644 Output Parameter: 645 . event - The event id for use with PetscLogEventBegin() and PetscLogEventEnd(). 646 647 Example of Usage: 648 .vb 649 PetscLogEvent USER_EVENT; 650 PetscClassId classid; 651 PetscLogDouble user_event_flops; 652 PetscClassIdRegister("class name",&classid); 653 PetscLogEventRegister("User event name",classid,&USER_EVENT); 654 PetscLogEventBegin(USER_EVENT,0,0,0,0); 655 [code segment to monitor] 656 PetscLogFlops(user_event_flops); 657 PetscLogEventEnd(USER_EVENT,0,0,0,0); 658 .ve 659 660 Notes: 661 PETSc automatically logs library events if the code has been 662 configured with --with-log (which is the default) and 663 -log_view or -log_all is specified. PetscLogEventRegister() is 664 intended for logging user events to supplement this PETSc 665 information. 666 667 PETSc can gather data for use with the utilities Jumpshot 668 (part of the MPICH distribution). If PETSc has been compiled 669 with flag -DPETSC_HAVE_MPE (MPE is an additional utility within 670 MPICH), the user can employ another command line option, -log_mpe, 671 to create a logfile, "mpe.log", which can be visualized 672 Jumpshot. 673 674 The classid is associated with each event so that classes of events 675 can be disabled simultaneously, such as all matrix events. The user 676 can either use an existing classid, such as MAT_CLASSID, or create 677 their own as shown in the example. 678 679 If an existing event with the same name exists, its event handle is 680 returned instead of creating a new event. 681 682 Level: intermediate 683 684 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogFlops(), 685 PetscLogEventActivate(), PetscLogEventDeactivate(), PetscClassIdRegister() 686 @*/ 687 PetscErrorCode PetscLogEventRegister(const char name[],PetscClassId classid,PetscLogEvent *event) 688 { 689 PetscStageLog stageLog; 690 int stage; 691 PetscErrorCode ierr; 692 693 PetscFunctionBegin; 694 *event = PETSC_DECIDE; 695 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 696 ierr = PetscEventRegLogGetEvent(stageLog->eventLog, name, event);CHKERRQ(ierr); 697 if (*event > 0) PetscFunctionReturn(0); 698 ierr = PetscEventRegLogRegister(stageLog->eventLog, name, classid, event);CHKERRQ(ierr); 699 for (stage = 0; stage < stageLog->numStages; stage++) { 700 ierr = PetscEventPerfLogEnsureSize(stageLog->stageInfo[stage].eventLog, stageLog->eventLog->numEvents);CHKERRQ(ierr); 701 ierr = PetscClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses);CHKERRQ(ierr); 702 } 703 PetscFunctionReturn(0); 704 } 705 706 /*@ 707 PetscLogEventSetCollective - Indicates that a particular event is collective. 708 709 Not Collective 710 711 Input Parameters: 712 + event - The event id 713 - collective - Bolean flag indicating whether a particular event is collective 714 715 Note: 716 New events returned from PetscLogEventRegister() are collective by default. 717 718 Level: developer 719 720 .seealso: PetscLogEventRegister() 721 @*/ 722 PetscErrorCode PetscLogEventSetCollective(PetscLogEvent event,PetscBool collective) 723 { 724 PetscStageLog stageLog; 725 PetscEventRegLog eventRegLog; 726 PetscErrorCode ierr; 727 728 PetscFunctionBegin; 729 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 730 ierr = PetscStageLogGetEventRegLog(stageLog,&eventRegLog);CHKERRQ(ierr); 731 if (event < 0 || event > eventRegLog->numEvents) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Invalid event id"); 732 eventRegLog->eventInfo[event].collective = collective; 733 PetscFunctionReturn(0); 734 } 735 736 /*@ 737 PetscLogEventIncludeClass - Activates event logging for a PETSc object class in every stage. 738 739 Not Collective 740 741 Input Parameter: 742 . classid - The object class, for example MAT_CLASSID, SNES_CLASSID, etc. 743 744 Level: developer 745 746 .seealso: PetscLogEventActivateClass(),PetscLogEventDeactivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate() 747 @*/ 748 PetscErrorCode PetscLogEventIncludeClass(PetscClassId classid) 749 { 750 PetscStageLog stageLog; 751 int stage; 752 PetscErrorCode ierr; 753 754 PetscFunctionBegin; 755 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 756 for (stage = 0; stage < stageLog->numStages; stage++) { 757 ierr = PetscEventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);CHKERRQ(ierr); 758 } 759 PetscFunctionReturn(0); 760 } 761 762 /*@ 763 PetscLogEventExcludeClass - Deactivates event logging for a PETSc object class in every stage. 764 765 Not Collective 766 767 Input Parameter: 768 . classid - The object class, for example MAT_CLASSID, SNES_CLASSID, etc. 769 770 Level: developer 771 772 .seealso: PetscLogEventDeactivateClass(),PetscLogEventActivateClass(),PetscLogEventDeactivate(),PetscLogEventActivate() 773 @*/ 774 PetscErrorCode PetscLogEventExcludeClass(PetscClassId classid) 775 { 776 PetscStageLog stageLog; 777 int stage; 778 PetscErrorCode ierr; 779 780 PetscFunctionBegin; 781 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 782 for (stage = 0; stage < stageLog->numStages; stage++) { 783 ierr = PetscEventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);CHKERRQ(ierr); 784 } 785 PetscFunctionReturn(0); 786 } 787 788 /*@ 789 PetscLogEventActivate - Indicates that a particular event should be logged. 790 791 Not Collective 792 793 Input Parameter: 794 . event - The event id 795 796 Usage: 797 .vb 798 PetscLogEventDeactivate(VEC_SetValues); 799 [code where you do not want to log VecSetValues()] 800 PetscLogEventActivate(VEC_SetValues); 801 [code where you do want to log VecSetValues()] 802 .ve 803 804 Note: 805 The event may be either a pre-defined PETSc event (found in include/petsclog.h) 806 or an event number obtained with PetscLogEventRegister(). 807 808 Level: advanced 809 810 .seealso: PlogEventDeactivate(), PlogEventDeactivatePush(), PetscLogEventDeactivatePop() 811 @*/ 812 PetscErrorCode PetscLogEventActivate(PetscLogEvent event) 813 { 814 PetscStageLog stageLog; 815 int stage; 816 PetscErrorCode ierr; 817 818 PetscFunctionBegin; 819 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 820 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 821 ierr = PetscEventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 822 PetscFunctionReturn(0); 823 } 824 825 /*@ 826 PetscLogEventDeactivate - Indicates that a particular event should not be logged. 827 828 Not Collective 829 830 Input Parameter: 831 . event - The event id 832 833 Usage: 834 .vb 835 PetscLogEventDeactivate(VEC_SetValues); 836 [code where you do not want to log VecSetValues()] 837 PetscLogEventActivate(VEC_SetValues); 838 [code where you do want to log VecSetValues()] 839 .ve 840 841 Note: 842 The event may be either a pre-defined PETSc event (found in 843 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 844 845 Level: advanced 846 847 .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePush(), PetscLogEventDeactivatePop() 848 @*/ 849 PetscErrorCode PetscLogEventDeactivate(PetscLogEvent event) 850 { 851 PetscStageLog stageLog; 852 int stage; 853 PetscErrorCode ierr; 854 855 PetscFunctionBegin; 856 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 857 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 858 ierr = PetscEventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 859 PetscFunctionReturn(0); 860 } 861 862 /*@ 863 PetscLogEventDeactivatePush - Indicates that a particular event should not be logged. 864 865 Not Collective 866 867 Input Parameter: 868 . event - The event id 869 870 Usage: 871 .vb 872 PetscLogEventDeactivatePush(VEC_SetValues); 873 [code where you do not want to log VecSetValues()] 874 PetscLogEventDeactivatePop(VEC_SetValues); 875 [code where you do want to log VecSetValues()] 876 .ve 877 878 Note: 879 The event may be either a pre-defined PETSc event (found in 880 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 881 882 Level: advanced 883 884 .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePop() 885 @*/ 886 PetscErrorCode PetscLogEventDeactivatePush(PetscLogEvent event) 887 { 888 PetscStageLog stageLog; 889 int stage; 890 PetscErrorCode ierr; 891 892 PetscFunctionBegin; 893 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 894 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 895 ierr = PetscEventPerfLogDeactivatePush(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 896 PetscFunctionReturn(0); 897 } 898 899 /*@ 900 PetscLogEventDeactivatePop - Indicates that a particular event shouldbe logged. 901 902 Not Collective 903 904 Input Parameter: 905 . event - The event id 906 907 Usage: 908 .vb 909 PetscLogEventDeactivatePush(VEC_SetValues); 910 [code where you do not want to log VecSetValues()] 911 PetscLogEventDeactivatePop(VEC_SetValues); 912 [code where you do want to log VecSetValues()] 913 .ve 914 915 Note: 916 The event may be either a pre-defined PETSc event (found in 917 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 918 919 Level: advanced 920 921 .seealso: PetscLogEventActivate(), PetscLogEventDeactivatePush() 922 @*/ 923 PetscErrorCode PetscLogEventDeactivatePop(PetscLogEvent event) 924 { 925 PetscStageLog stageLog; 926 int stage; 927 PetscErrorCode ierr; 928 929 PetscFunctionBegin; 930 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 931 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 932 ierr = PetscEventPerfLogDeactivatePop(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 933 PetscFunctionReturn(0); 934 } 935 936 /*@ 937 PetscLogEventSetActiveAll - Sets the event activity in every stage. 938 939 Not Collective 940 941 Input Parameters: 942 + event - The event id 943 - isActive - The activity flag determining whether the event is logged 944 945 Level: advanced 946 947 .seealso: PlogEventActivate(),PlogEventDeactivate() 948 @*/ 949 PetscErrorCode PetscLogEventSetActiveAll(PetscLogEvent event, PetscBool isActive) 950 { 951 PetscStageLog stageLog; 952 int stage; 953 PetscErrorCode ierr; 954 955 PetscFunctionBegin; 956 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 957 for (stage = 0; stage < stageLog->numStages; stage++) { 958 if (isActive) { 959 ierr = PetscEventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 960 } else { 961 ierr = PetscEventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 962 } 963 } 964 PetscFunctionReturn(0); 965 } 966 967 /*@ 968 PetscLogEventActivateClass - Activates event logging for a PETSc object class. 969 970 Not Collective 971 972 Input Parameter: 973 . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc. 974 975 Level: developer 976 977 .seealso: PetscLogEventDeactivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate() 978 @*/ 979 PetscErrorCode PetscLogEventActivateClass(PetscClassId classid) 980 { 981 PetscStageLog stageLog; 982 int stage; 983 PetscErrorCode ierr; 984 985 PetscFunctionBegin; 986 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 987 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 988 ierr = PetscEventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);CHKERRQ(ierr); 989 PetscFunctionReturn(0); 990 } 991 992 /*@ 993 PetscLogEventDeactivateClass - Deactivates event logging for a PETSc object class. 994 995 Not Collective 996 997 Input Parameter: 998 . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc. 999 1000 Level: developer 1001 1002 .seealso: PetscLogEventActivateClass(),PetscLogEventActivate(),PetscLogEventDeactivate() 1003 @*/ 1004 PetscErrorCode PetscLogEventDeactivateClass(PetscClassId classid) 1005 { 1006 PetscStageLog stageLog; 1007 int stage; 1008 PetscErrorCode ierr; 1009 1010 PetscFunctionBegin; 1011 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1012 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 1013 ierr = PetscEventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);CHKERRQ(ierr); 1014 PetscFunctionReturn(0); 1015 } 1016 1017 /*MC 1018 PetscLogEventSync - Synchronizes the beginning of a user event. 1019 1020 Synopsis: 1021 #include <petsclog.h> 1022 PetscErrorCode PetscLogEventSync(int e,MPI_Comm comm) 1023 1024 Collective 1025 1026 Input Parameters: 1027 + e - integer associated with the event obtained from PetscLogEventRegister() 1028 - comm - an MPI communicator 1029 1030 Usage: 1031 .vb 1032 PetscLogEvent USER_EVENT; 1033 PetscLogEventRegister("User event",0,&USER_EVENT); 1034 PetscLogEventSync(USER_EVENT,PETSC_COMM_WORLD); 1035 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1036 [code segment to monitor] 1037 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1038 .ve 1039 1040 Notes: 1041 This routine should be called only if there is not a 1042 PetscObject available to pass to PetscLogEventBegin(). 1043 1044 Level: developer 1045 1046 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd() 1047 1048 M*/ 1049 1050 /*MC 1051 PetscLogEventBegin - Logs the beginning of a user event. 1052 1053 Synopsis: 1054 #include <petsclog.h> 1055 PetscErrorCode PetscLogEventBegin(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 1056 1057 Not Collective 1058 1059 Input Parameters: 1060 + e - integer associated with the event obtained from PetscLogEventRegister() 1061 - o1,o2,o3,o4 - objects associated with the event, or 0 1062 1063 Fortran Synopsis: 1064 void PetscLogEventBegin(int e,PetscErrorCode ierr) 1065 1066 Usage: 1067 .vb 1068 PetscLogEvent USER_EVENT; 1069 PetscLogDouble user_event_flops; 1070 PetscLogEventRegister("User event",0,&USER_EVENT); 1071 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1072 [code segment to monitor] 1073 PetscLogFlops(user_event_flops); 1074 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1075 .ve 1076 1077 Notes: 1078 You need to register each integer event with the command 1079 PetscLogEventRegister(). 1080 1081 Level: intermediate 1082 1083 .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops() 1084 1085 M*/ 1086 1087 /*MC 1088 PetscLogEventEnd - Log the end of a user event. 1089 1090 Synopsis: 1091 #include <petsclog.h> 1092 PetscErrorCode PetscLogEventEnd(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 1093 1094 Not Collective 1095 1096 Input Parameters: 1097 + e - integer associated with the event obtained with PetscLogEventRegister() 1098 - o1,o2,o3,o4 - objects associated with the event, or 0 1099 1100 Fortran Synopsis: 1101 void PetscLogEventEnd(int e,PetscErrorCode ierr) 1102 1103 Usage: 1104 .vb 1105 PetscLogEvent USER_EVENT; 1106 PetscLogDouble user_event_flops; 1107 PetscLogEventRegister("User event",0,&USER_EVENT,); 1108 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1109 [code segment to monitor] 1110 PetscLogFlops(user_event_flops); 1111 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1112 .ve 1113 1114 Notes: 1115 You should also register each additional integer event with the command 1116 PetscLogEventRegister(). 1117 1118 Level: intermediate 1119 1120 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogFlops() 1121 1122 M*/ 1123 1124 /*@C 1125 PetscLogEventGetId - Returns the event id when given the event name. 1126 1127 Not Collective 1128 1129 Input Parameter: 1130 . name - The event name 1131 1132 Output Parameter: 1133 . event - The event, or -1 if no event with that name exists 1134 1135 Level: intermediate 1136 1137 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogStageGetId() 1138 @*/ 1139 PetscErrorCode PetscLogEventGetId(const char name[], PetscLogEvent *event) 1140 { 1141 PetscStageLog stageLog; 1142 PetscErrorCode ierr; 1143 1144 PetscFunctionBegin; 1145 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1146 ierr = PetscEventRegLogGetEvent(stageLog->eventLog, name, event);CHKERRQ(ierr); 1147 PetscFunctionReturn(0); 1148 } 1149 1150 /*------------------------------------------------ Output Functions -------------------------------------------------*/ 1151 /*@C 1152 PetscLogDump - Dumps logs of objects to a file. This file is intended to 1153 be read by bin/petscview. This program no longer exists. 1154 1155 Collective on PETSC_COMM_WORLD 1156 1157 Input Parameter: 1158 . name - an optional file name 1159 1160 Usage: 1161 .vb 1162 PetscInitialize(...); 1163 PetscLogDefaultBegin(); or PetscLogAllBegin(); 1164 ... code ... 1165 PetscLogDump(filename); 1166 PetscFinalize(); 1167 .ve 1168 1169 Notes: 1170 The default file name is 1171 $ Log.<rank> 1172 where <rank> is the processor number. If no name is specified, 1173 this file will be used. 1174 1175 Level: advanced 1176 1177 .seealso: PetscLogDefaultBegin(), PetscLogAllBegin(), PetscLogView() 1178 @*/ 1179 PetscErrorCode PetscLogDump(const char sname[]) 1180 { 1181 PetscStageLog stageLog; 1182 PetscEventPerfInfo *eventInfo; 1183 FILE *fd; 1184 char file[PETSC_MAX_PATH_LEN], fname[PETSC_MAX_PATH_LEN]; 1185 PetscLogDouble flops, _TotalTime; 1186 PetscMPIInt rank; 1187 int action, object, curStage; 1188 PetscLogEvent event; 1189 PetscErrorCode ierr; 1190 1191 PetscFunctionBegin; 1192 /* Calculate the total elapsed time */ 1193 PetscTime(&_TotalTime); 1194 _TotalTime -= petsc_BaseTime; 1195 /* Open log file */ 1196 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank);CHKERRMPI(ierr); 1197 if (sname && sname[0]) sprintf(file, "%s.%d", sname, rank); 1198 else sprintf(file, "Log.%d", rank); 1199 ierr = PetscFixFilename(file, fname);CHKERRQ(ierr); 1200 ierr = PetscFOpen(PETSC_COMM_WORLD, fname, "w", &fd);CHKERRQ(ierr); 1201 if ((rank == 0) && (!fd)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN, "Cannot open file: %s", fname); 1202 /* Output totals */ 1203 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Total Flop %14e %16.8e\n", petsc_TotalFlops, _TotalTime);CHKERRQ(ierr); 1204 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Clock Resolution %g\n", 0.0);CHKERRQ(ierr); 1205 /* Output actions */ 1206 if (petsc_logActions) { 1207 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Actions accomplished %d\n", petsc_numActions);CHKERRQ(ierr); 1208 for (action = 0; action < petsc_numActions; action++) { 1209 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "%g %d %d %d %d %d %d %g %g %g\n", 1210 petsc_actions[action].time, petsc_actions[action].action, (int)petsc_actions[action].event, (int)petsc_actions[action].classid, petsc_actions[action].id1, 1211 petsc_actions[action].id2, petsc_actions[action].id3, petsc_actions[action].flops, petsc_actions[action].mem, petsc_actions[action].maxmem);CHKERRQ(ierr); 1212 } 1213 } 1214 /* Output objects */ 1215 if (petsc_logObjects) { 1216 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Objects created %d destroyed %d\n", petsc_numObjects, petsc_numObjectsDestroyed);CHKERRQ(ierr); 1217 for (object = 0; object < petsc_numObjects; object++) { 1218 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Parent ID: %d Memory: %d\n", petsc_objects[object].parent, (int) petsc_objects[object].mem);CHKERRQ(ierr); 1219 if (!petsc_objects[object].name[0]) { 1220 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd,"No Name\n");CHKERRQ(ierr); 1221 } else { 1222 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Name: %s\n", petsc_objects[object].name);CHKERRQ(ierr); 1223 } 1224 if (petsc_objects[object].info[0] != 0) { 1225 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "No Info\n");CHKERRQ(ierr); 1226 } else { 1227 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Info: %s\n", petsc_objects[object].info);CHKERRQ(ierr); 1228 } 1229 } 1230 } 1231 /* Output events */ 1232 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Event log:\n");CHKERRQ(ierr); 1233 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1234 ierr = PetscIntStackTop(stageLog->stack, &curStage);CHKERRQ(ierr); 1235 eventInfo = stageLog->stageInfo[curStage].eventLog->eventInfo; 1236 for (event = 0; event < stageLog->stageInfo[curStage].eventLog->numEvents; event++) { 1237 if (eventInfo[event].time != 0.0) flops = eventInfo[event].flops/eventInfo[event].time; 1238 else flops = 0.0; 1239 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "%d %16d %16g %16g %16g\n", event, eventInfo[event].count, 1240 eventInfo[event].flops, eventInfo[event].time, flops);CHKERRQ(ierr); 1241 } 1242 ierr = PetscFClose(PETSC_COMM_WORLD, fd);CHKERRQ(ierr); 1243 PetscFunctionReturn(0); 1244 } 1245 1246 /* 1247 PetscLogView_Detailed - Each process prints the times for its own events 1248 1249 */ 1250 PetscErrorCode PetscLogView_Detailed(PetscViewer viewer) 1251 { 1252 PetscStageLog stageLog; 1253 PetscEventPerfInfo *eventInfo = NULL, *stageInfo = NULL; 1254 PetscLogDouble locTotalTime, numRed, maxMem; 1255 int numStages,numEvents,stage,event; 1256 MPI_Comm comm = PetscObjectComm((PetscObject) viewer); 1257 PetscMPIInt rank,size; 1258 PetscErrorCode ierr; 1259 1260 PetscFunctionBegin; 1261 ierr = MPI_Comm_size(comm, &size);CHKERRMPI(ierr); 1262 ierr = MPI_Comm_rank(comm, &rank);CHKERRMPI(ierr); 1263 /* Must preserve reduction count before we go on */ 1264 numRed = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1265 /* Get the total elapsed time */ 1266 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1267 ierr = PetscViewerASCIIPrintf(viewer,"size = %d\n",size);CHKERRQ(ierr); 1268 ierr = PetscViewerASCIIPrintf(viewer,"LocalTimes = {}\n");CHKERRQ(ierr); 1269 ierr = PetscViewerASCIIPrintf(viewer,"LocalMessages = {}\n");CHKERRQ(ierr); 1270 ierr = PetscViewerASCIIPrintf(viewer,"LocalMessageLens = {}\n");CHKERRQ(ierr); 1271 ierr = PetscViewerASCIIPrintf(viewer,"LocalReductions = {}\n");CHKERRQ(ierr); 1272 ierr = PetscViewerASCIIPrintf(viewer,"LocalFlop = {}\n");CHKERRQ(ierr); 1273 ierr = PetscViewerASCIIPrintf(viewer,"LocalObjects = {}\n");CHKERRQ(ierr); 1274 ierr = PetscViewerASCIIPrintf(viewer,"LocalMemory = {}\n");CHKERRQ(ierr); 1275 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1276 ierr = MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);CHKERRMPI(ierr); 1277 ierr = PetscViewerASCIIPrintf(viewer,"Stages = {}\n");CHKERRQ(ierr); 1278 for (stage=0; stage<numStages; stage++) { 1279 ierr = PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"] = {}\n",stageLog->stageInfo[stage].name);CHKERRQ(ierr); 1280 ierr = PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"summary\"] = {}\n",stageLog->stageInfo[stage].name);CHKERRQ(ierr); 1281 ierr = MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);CHKERRMPI(ierr); 1282 for (event = 0; event < numEvents; event++) { 1283 ierr = PetscViewerASCIIPrintf(viewer,"Stages[\"%s\"][\"%s\"] = {}\n",stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name);CHKERRQ(ierr); 1284 } 1285 } 1286 ierr = PetscMallocGetMaximumUsage(&maxMem);CHKERRQ(ierr); 1287 ierr = PetscViewerASCIIPushSynchronized(viewer);CHKERRQ(ierr); 1288 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalTimes[%d] = %g\n",rank,locTotalTime);CHKERRQ(ierr); 1289 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessages[%d] = %g\n",rank,(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct));CHKERRQ(ierr); 1290 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalMessageLens[%d] = %g\n",rank,(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len));CHKERRQ(ierr); 1291 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalReductions[%d] = %g\n",rank,numRed);CHKERRQ(ierr); 1292 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalFlop[%d] = %g\n",rank,petsc_TotalFlops);CHKERRQ(ierr); 1293 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalObjects[%d] = %d\n",rank,petsc_numObjects);CHKERRQ(ierr); 1294 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"LocalMemory[%d] = %g\n",rank,maxMem);CHKERRQ(ierr); 1295 ierr = PetscViewerFlush(viewer);CHKERRQ(ierr); 1296 for (stage=0; stage<numStages; stage++) { 1297 stageInfo = &stageLog->stageInfo[stage].perfInfo; 1298 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"summary\"][%d] = {\"time\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g}\n", 1299 stageLog->stageInfo[stage].name,rank, 1300 stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops);CHKERRQ(ierr); 1301 ierr = MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);CHKERRMPI(ierr); 1302 for (event = 0; event < numEvents; event++) { 1303 eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event]; 1304 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"Stages[\"%s\"][\"%s\"][%d] = {\"count\" : %D, \"time\" : %g, \"syncTime\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g", 1305 stageLog->stageInfo[stage].name,stageLog->eventLog->eventInfo[event].name,rank, 1306 eventInfo->count,eventInfo->time,eventInfo->syncTime,eventInfo->numMessages,eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops);CHKERRQ(ierr); 1307 if (eventInfo->dof[0] >= 0.) { 1308 PetscInt d, e; 1309 1310 ierr = PetscViewerASCIISynchronizedPrintf(viewer, ", \"dof\" : [");CHKERRQ(ierr); 1311 for (d = 0; d < 8; ++d) { 1312 if (d > 0) {ierr = PetscViewerASCIISynchronizedPrintf(viewer, ", ");CHKERRQ(ierr);} 1313 ierr = PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->dof[d]);CHKERRQ(ierr); 1314 } 1315 ierr = PetscViewerASCIISynchronizedPrintf(viewer, "]");CHKERRQ(ierr); 1316 ierr = PetscViewerASCIISynchronizedPrintf(viewer, ", \"error\" : [");CHKERRQ(ierr); 1317 for (e = 0; e < 8; ++e) { 1318 if (e > 0) {ierr = PetscViewerASCIISynchronizedPrintf(viewer, ", ");CHKERRQ(ierr);} 1319 ierr = PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->errors[e]);CHKERRQ(ierr); 1320 } 1321 ierr = PetscViewerASCIISynchronizedPrintf(viewer, "]");CHKERRQ(ierr); 1322 } 1323 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"}\n");CHKERRQ(ierr); 1324 } 1325 } 1326 ierr = PetscViewerFlush(viewer);CHKERRQ(ierr); 1327 ierr = PetscViewerASCIIPopSynchronized(viewer);CHKERRQ(ierr); 1328 PetscFunctionReturn(0); 1329 } 1330 1331 /* 1332 PetscLogView_CSV - Each process prints the times for its own events in Comma-Separated Value Format 1333 */ 1334 PetscErrorCode PetscLogView_CSV(PetscViewer viewer) 1335 { 1336 PetscStageLog stageLog; 1337 PetscEventPerfInfo *eventInfo = NULL; 1338 PetscLogDouble locTotalTime, maxMem; 1339 int numStages,numEvents,stage,event; 1340 MPI_Comm comm = PetscObjectComm((PetscObject) viewer); 1341 PetscMPIInt rank,size; 1342 PetscErrorCode ierr; 1343 1344 PetscFunctionBegin; 1345 ierr = MPI_Comm_size(comm, &size);CHKERRMPI(ierr); 1346 ierr = MPI_Comm_rank(comm, &rank);CHKERRMPI(ierr); 1347 /* Must preserve reduction count before we go on */ 1348 /* Get the total elapsed time */ 1349 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1350 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1351 ierr = MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);CHKERRMPI(ierr); 1352 ierr = PetscMallocGetMaximumUsage(&maxMem);CHKERRQ(ierr); 1353 ierr = PetscViewerASCIIPushSynchronized(viewer);CHKERRQ(ierr); 1354 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); 1355 ierr = PetscViewerFlush(viewer);CHKERRQ(ierr); 1356 for (stage=0; stage<numStages; stage++) { 1357 PetscEventPerfInfo *stageInfo = &stageLog->stageInfo[stage].perfInfo; 1358 1359 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"%s,summary,%d,1,%g,%g,%g,%g,%g\n", 1360 stageLog->stageInfo[stage].name,rank,stageInfo->time,stageInfo->numMessages,stageInfo->messageLength,stageInfo->numReductions,stageInfo->flops);CHKERRQ(ierr); 1361 ierr = MPI_Allreduce(&stageLog->stageInfo[stage].eventLog->numEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);CHKERRMPI(ierr); 1362 for (event = 0; event < numEvents; event++) { 1363 eventInfo = &stageLog->stageInfo[stage].eventLog->eventInfo[event]; 1364 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"%s,%s,%d,%D,%g,%g,%g,%g,%g",stageLog->stageInfo[stage].name, 1365 stageLog->eventLog->eventInfo[event].name,rank,eventInfo->count,eventInfo->time,eventInfo->numMessages, 1366 eventInfo->messageLength,eventInfo->numReductions,eventInfo->flops);CHKERRQ(ierr); 1367 if (eventInfo->dof[0] >= 0.) { 1368 PetscInt d, e; 1369 1370 for (d = 0; d < 8; ++d) { 1371 ierr = PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->dof[d]);CHKERRQ(ierr); 1372 } 1373 for (e = 0; e < 8; ++e) { 1374 ierr = PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->errors[e]);CHKERRQ(ierr); 1375 } 1376 } 1377 ierr = PetscViewerASCIISynchronizedPrintf(viewer,"\n");CHKERRQ(ierr); 1378 } 1379 } 1380 ierr = PetscViewerFlush(viewer);CHKERRQ(ierr); 1381 ierr = PetscViewerASCIIPopSynchronized(viewer);CHKERRQ(ierr); 1382 PetscFunctionReturn(0); 1383 } 1384 1385 static PetscErrorCode PetscLogViewWarnSync(MPI_Comm comm,FILE *fd) 1386 { 1387 PetscErrorCode ierr; 1388 PetscFunctionBegin; 1389 if (!PetscLogSyncOn) PetscFunctionReturn(0); 1390 ierr = PetscFPrintf(comm, fd, "\n\n");CHKERRQ(ierr); 1391 ierr = PetscFPrintf(comm, fd, " ##########################################################\n");CHKERRQ(ierr); 1392 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1393 ierr = PetscFPrintf(comm, fd, " # WARNING!!! #\n");CHKERRQ(ierr); 1394 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1395 ierr = PetscFPrintf(comm, fd, " # This program was run with logging synchronization. #\n");CHKERRQ(ierr); 1396 ierr = PetscFPrintf(comm, fd, " # This option provides more meaningful imbalance #\n");CHKERRQ(ierr); 1397 ierr = PetscFPrintf(comm, fd, " # figures at the expense of slowing things down and #\n");CHKERRQ(ierr); 1398 ierr = PetscFPrintf(comm, fd, " # providing a distorted view of the overall runtime. #\n");CHKERRQ(ierr); 1399 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1400 ierr = PetscFPrintf(comm, fd, " ##########################################################\n\n\n");CHKERRQ(ierr); 1401 PetscFunctionReturn(0); 1402 } 1403 1404 static PetscErrorCode PetscLogViewWarnDebugging(MPI_Comm comm,FILE *fd) 1405 { 1406 PetscErrorCode ierr; 1407 1408 PetscFunctionBegin; 1409 if (PetscDefined(USE_DEBUG)) { 1410 ierr = PetscFPrintf(comm, fd, "\n\n");CHKERRQ(ierr); 1411 ierr = PetscFPrintf(comm, fd, " ##########################################################\n");CHKERRQ(ierr); 1412 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1413 ierr = PetscFPrintf(comm, fd, " # WARNING!!! #\n");CHKERRQ(ierr); 1414 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1415 ierr = PetscFPrintf(comm, fd, " # This code was compiled with a debugging option. #\n");CHKERRQ(ierr); 1416 ierr = PetscFPrintf(comm, fd, " # To get timing results run ./configure #\n");CHKERRQ(ierr); 1417 ierr = PetscFPrintf(comm, fd, " # using --with-debugging=no, the performance will #\n");CHKERRQ(ierr); 1418 ierr = PetscFPrintf(comm, fd, " # be generally two or three times faster. #\n");CHKERRQ(ierr); 1419 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1420 ierr = PetscFPrintf(comm, fd, " ##########################################################\n\n\n");CHKERRQ(ierr); 1421 } 1422 PetscFunctionReturn(0); 1423 } 1424 1425 static PetscErrorCode PetscLogViewWarnNoGpuAwareMpi(MPI_Comm comm,FILE *fd) 1426 { 1427 #if defined(PETSC_HAVE_CUDA) || defined(PETSC_HAVE_HIP) 1428 PetscErrorCode ierr; 1429 1430 PetscFunctionBegin; 1431 if (use_gpu_aware_mpi || !PetscCreatedGpuObjects) PetscFunctionReturn(0); 1432 ierr = PetscFPrintf(comm, fd, "\n\n");CHKERRQ(ierr); 1433 ierr = PetscFPrintf(comm, fd, " ##########################################################\n");CHKERRQ(ierr); 1434 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1435 ierr = PetscFPrintf(comm, fd, " # WARNING!!! #\n");CHKERRQ(ierr); 1436 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1437 ierr = PetscFPrintf(comm, fd, " # This code was compiled with GPU support and you've #\n");CHKERRQ(ierr); 1438 ierr = PetscFPrintf(comm, fd, " # created PETSc/GPU objects, but you intentionally used #\n");CHKERRQ(ierr); 1439 ierr = PetscFPrintf(comm, fd, " # -use_gpu_aware_mpi 0, such that PETSc had to copy data #\n");CHKERRQ(ierr); 1440 ierr = PetscFPrintf(comm, fd, " # from GPU to CPU for communication. To get meaningfull #\n");CHKERRQ(ierr); 1441 ierr = PetscFPrintf(comm, fd, " # timing results, please use GPU-aware MPI instead. #\n");CHKERRQ(ierr); 1442 ierr = PetscFPrintf(comm, fd, " ##########################################################\n\n\n");CHKERRQ(ierr); 1443 PetscFunctionReturn(0); 1444 #else 1445 return 0; 1446 #endif 1447 } 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 Parameters: 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 <petscdevice.h> 2245 PETSC_EXTERN cudaEvent_t petsc_gputimer_begin; 2246 PETSC_EXTERN cudaEvent_t petsc_gputimer_end; 2247 #endif 2248 2249 #if defined(PETSC_HAVE_HIP) 2250 #include <petscdevice.h> 2251 PETSC_EXTERN hipEvent_t petsc_gputimer_begin; 2252 PETSC_EXTERN hipEvent_t petsc_gputimer_end; 2253 #endif 2254 2255 /*-------------------------------------------- GPU event Functions ----------------------------------------------*/ 2256 /*@C 2257 PetscLogGpuTimeBegin - Start timer for device 2258 2259 Notes: 2260 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). 2261 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). 2262 There is no need to call WaitForCUDA() or WaitForHIP() between PetscLogGpuTimeBegin and PetscLogGpuTimeEnd 2263 This timer should NOT include times for data transfers between the GPU and CPU, nor setup actions such as allocating space. 2264 The regular logging captures the time for data transfers and any CPU activites during the event 2265 It is used to compute the flop rate on the GPU as it is actively engaged in running a kernel. 2266 2267 Developer Notes: 2268 The GPU event timer captures the execution time of all the kernels launched in the default stream by the CPU between PetscLogGpuTimeBegin() and PetsLogGpuTimeEnd(). 2269 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. 2270 2271 Level: intermediate 2272 2273 .seealso: PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeEnd() 2274 @*/ 2275 PetscErrorCode PetscLogGpuTimeBegin(void) 2276 { 2277 #if defined(PETSC_HAVE_CUDA) 2278 cudaError_t cerr; 2279 #elif defined(PETSC_HAVE_HIP) 2280 hipError_t cerr; 2281 #else 2282 PetscErrorCode ierr; 2283 #endif 2284 PetscFunctionBegin; 2285 #if defined(PETSC_USE_DEBUG) 2286 if (petsc_gtime_inuse) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Forgot to call PetscLogGpuTimeEnd()?"); 2287 petsc_gtime_inuse = PETSC_TRUE; 2288 #endif 2289 if (!PetscLogPLB) PetscFunctionReturn(0); 2290 #if defined(PETSC_HAVE_CUDA) 2291 cerr = cudaEventRecord(petsc_gputimer_begin,PetscDefaultCudaStream);CHKERRCUDA(cerr); 2292 #elif defined(PETSC_HAVE_HIP) 2293 cerr = hipEventRecord(petsc_gputimer_begin,PetscDefaultHipStream);CHKERRHIP(cerr); 2294 #else 2295 ierr = PetscTimeSubtract(&petsc_gtime);CHKERRQ(ierr); 2296 #endif 2297 PetscFunctionReturn(0); 2298 } 2299 2300 /*@C 2301 PetscLogGpuTimeEnd - Stop timer for device 2302 2303 Level: intermediate 2304 2305 .seealso: PetscLogView(), PetscLogGpuFlops(), PetscLogGpuTimeBegin() 2306 @*/ 2307 PetscErrorCode PetscLogGpuTimeEnd(void) 2308 { 2309 #if defined(PETSC_HAVE_CUDA) 2310 float gtime; 2311 cudaError_t cerr; 2312 #elif defined(PETSC_HAVE_HIP) 2313 float gtime; 2314 hipError_t cerr; 2315 #else 2316 PetscErrorCode ierr; 2317 #endif 2318 PetscFunctionBegin; 2319 #if defined(PETSC_USE_DEBUG) 2320 if (!petsc_gtime_inuse) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Forgot to call PetscLogGpuTimeBegin()?"); 2321 petsc_gtime_inuse = PETSC_FALSE; 2322 #endif 2323 if (!PetscLogPLE) PetscFunctionReturn(0); 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