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