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 #include <petscthreadcomm.h> 15 16 PetscErrorCode PetscLogObjectParent(PetscObject p,PetscObject c) 17 { 18 if (c && p) { 19 PetscObject pp = p; 20 while (!c->parent && pp) { 21 /* if not credited elsewhere credit all childs memory to all new ancestors */ 22 pp->mem += c->mem; 23 pp = pp->parent; 24 } 25 c->parent = p; 26 c->parentid = p->id; 27 } 28 return 0; 29 } 30 31 PetscErrorCode PetscLogObjectMemory(PetscObject p,PetscLogDouble m) 32 { 33 while (p) { 34 /* Create all ancestors with the memory */ 35 p->mem += m; 36 p = p->parent; 37 } 38 return 0; 39 } 40 41 PetscLogEvent PETSC_LARGEST_EVENT = PETSC_EVENT; 42 43 #if defined(PETSC_CLANGUAGE_CXX) 44 std::map<std::string,PETSc::LogEvent> PETSc::Log::event_registry; 45 std::map<std::string,PETSc::LogStage> PETSc::Log::stage_registry; 46 #endif 47 48 #if defined(PETSC_USE_LOG) 49 #include <petscmachineinfo.h> 50 #include <petscconfiginfo.h> 51 52 /* used in the MPI_XXX() count macros in petsclog.h */ 53 54 /* Action and object logging variables */ 55 Action *petsc_actions = NULL; 56 Object *petsc_objects = NULL; 57 PetscBool petsc_logActions = PETSC_FALSE; 58 PetscBool petsc_logObjects = PETSC_FALSE; 59 int petsc_numActions = 0, petsc_maxActions = 100; 60 int petsc_numObjects = 0, petsc_maxObjects = 100; 61 int petsc_numObjectsDestroyed = 0; 62 63 /* Global counters */ 64 PetscLogDouble petsc_BaseTime = 0.0; 65 PetscLogDouble petsc_TotalFlops = 0.0; /* The number of flops */ 66 PetscLogDouble petsc_tmp_flops = 0.0; /* The incremental number of flops */ 67 PetscLogDouble petsc_send_ct = 0.0; /* The number of sends */ 68 PetscLogDouble petsc_recv_ct = 0.0; /* The number of receives */ 69 PetscLogDouble petsc_send_len = 0.0; /* The total length of all sent messages */ 70 PetscLogDouble petsc_recv_len = 0.0; /* The total length of all received messages */ 71 PetscLogDouble petsc_isend_ct = 0.0; /* The number of immediate sends */ 72 PetscLogDouble petsc_irecv_ct = 0.0; /* The number of immediate receives */ 73 PetscLogDouble petsc_isend_len = 0.0; /* The total length of all immediate send messages */ 74 PetscLogDouble petsc_irecv_len = 0.0; /* The total length of all immediate receive messages */ 75 PetscLogDouble petsc_wait_ct = 0.0; /* The number of waits */ 76 PetscLogDouble petsc_wait_any_ct = 0.0; /* The number of anywaits */ 77 PetscLogDouble petsc_wait_all_ct = 0.0; /* The number of waitalls */ 78 PetscLogDouble petsc_sum_of_waits_ct = 0.0; /* The total number of waits */ 79 PetscLogDouble petsc_allreduce_ct = 0.0; /* The number of reductions */ 80 PetscLogDouble petsc_gather_ct = 0.0; /* The number of gathers and gathervs */ 81 PetscLogDouble petsc_scatter_ct = 0.0; /* The number of scatters and scattervs */ 82 83 /* Logging functions */ 84 PetscErrorCode (*PetscLogPHC)(PetscObject) = NULL; 85 PetscErrorCode (*PetscLogPHD)(PetscObject) = NULL; 86 PetscErrorCode (*PetscLogPLB)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL; 87 PetscErrorCode (*PetscLogPLE)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject) = NULL; 88 89 /* Tracing event logging variables */ 90 FILE *petsc_tracefile = NULL; 91 int petsc_tracelevel = 0; 92 const char *petsc_traceblanks = " "; 93 char petsc_tracespace[128] = " "; 94 PetscLogDouble petsc_tracetime = 0.0; 95 static PetscBool PetscLogBegin_PrivateCalled = PETSC_FALSE; 96 97 /*---------------------------------------------- General Functions --------------------------------------------------*/ 98 #undef __FUNCT__ 99 #define __FUNCT__ "PetscLogDestroy" 100 /*@C 101 PetscLogDestroy - Destroys the object and event logging data and resets the global counters. 102 103 Not Collective 104 105 Notes: 106 This routine should not usually be used by programmers. Instead employ 107 PetscLogStagePush() and PetscLogStagePop(). 108 109 Level: developer 110 111 .keywords: log, destroy 112 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogStagePush(), PlogStagePop() 113 @*/ 114 PetscErrorCode PetscLogDestroy(void) 115 { 116 PetscStageLog stageLog; 117 PetscErrorCode ierr; 118 119 PetscFunctionBegin; 120 ierr = PetscFree(petsc_actions);CHKERRQ(ierr); 121 ierr = PetscFree(petsc_objects);CHKERRQ(ierr); 122 ierr = PetscLogSet(NULL, NULL);CHKERRQ(ierr); 123 124 /* Resetting phase */ 125 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 126 ierr = PetscStageLogDestroy(stageLog);CHKERRQ(ierr); 127 128 petsc_TotalFlops = 0.0; 129 petsc_numActions = 0; 130 petsc_numObjects = 0; 131 petsc_numObjectsDestroyed = 0; 132 petsc_maxActions = 100; 133 petsc_maxObjects = 100; 134 petsc_actions = NULL; 135 petsc_objects = NULL; 136 petsc_logActions = PETSC_FALSE; 137 petsc_logObjects = PETSC_FALSE; 138 petsc_BaseTime = 0.0; 139 petsc_TotalFlops = 0.0; 140 petsc_tmp_flops = 0.0; 141 petsc_send_ct = 0.0; 142 petsc_recv_ct = 0.0; 143 petsc_send_len = 0.0; 144 petsc_recv_len = 0.0; 145 petsc_isend_ct = 0.0; 146 petsc_irecv_ct = 0.0; 147 petsc_isend_len = 0.0; 148 petsc_irecv_len = 0.0; 149 petsc_wait_ct = 0.0; 150 petsc_wait_any_ct = 0.0; 151 petsc_wait_all_ct = 0.0; 152 petsc_sum_of_waits_ct = 0.0; 153 petsc_allreduce_ct = 0.0; 154 petsc_gather_ct = 0.0; 155 petsc_scatter_ct = 0.0; 156 PETSC_LARGEST_EVENT = PETSC_EVENT; 157 PetscLogPHC = NULL; 158 PetscLogPHD = NULL; 159 petsc_tracefile = NULL; 160 petsc_tracelevel = 0; 161 petsc_traceblanks = " "; 162 petsc_tracespace[0] = ' '; petsc_tracespace[1] = 0; 163 petsc_tracetime = 0.0; 164 PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID; 165 PETSC_OBJECT_CLASSID = 0; 166 petsc_stageLog = 0; 167 PetscLogBegin_PrivateCalled = PETSC_FALSE; 168 PetscFunctionReturn(0); 169 } 170 171 #undef __FUNCT__ 172 #define __FUNCT__ "PetscLogSet" 173 /*@C 174 PetscLogSet - Sets the logging functions called at the beginning and ending of every event. 175 176 Not Collective 177 178 Input Parameters: 179 + b - The function called at beginning of event 180 - e - The function called at end of event 181 182 Level: developer 183 184 .seealso: PetscLogDump(), PetscLogBegin(), PetscLogAllBegin(), PetscLogTraceBegin() 185 @*/ 186 PetscErrorCode PetscLogSet(PetscErrorCode (*b)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject), 187 PetscErrorCode (*e)(PetscLogEvent, int, PetscObject, PetscObject, PetscObject, PetscObject)) 188 { 189 PetscFunctionBegin; 190 PetscLogPLB = b; 191 PetscLogPLE = e; 192 PetscFunctionReturn(0); 193 } 194 195 #if defined(PETSC_HAVE_CHUD) 196 #include <CHUD/CHUD.h> 197 #endif 198 #if defined(PETSC_HAVE_PAPI) 199 #include <papi.h> 200 int PAPIEventSet = PAPI_NULL; 201 #endif 202 203 /*------------------------------------------- Initialization Functions ----------------------------------------------*/ 204 #undef __FUNCT__ 205 #define __FUNCT__ "PetscLogBegin_Private" 206 PetscErrorCode PetscLogBegin_Private(void) 207 { 208 int stage; 209 PetscBool opt; 210 PetscErrorCode ierr; 211 212 PetscFunctionBegin; 213 if (PetscLogBegin_PrivateCalled) PetscFunctionReturn(0); 214 PetscLogBegin_PrivateCalled = PETSC_TRUE; 215 216 ierr = PetscOptionsHasName(NULL, "-log_exclude_actions", &opt);CHKERRQ(ierr); 217 if (opt) petsc_logActions = PETSC_FALSE; 218 ierr = PetscOptionsHasName(NULL, "-log_exclude_objects", &opt);CHKERRQ(ierr); 219 if (opt) petsc_logObjects = PETSC_FALSE; 220 if (petsc_logActions) { 221 ierr = PetscMalloc(petsc_maxActions * sizeof(Action), &petsc_actions);CHKERRQ(ierr); 222 } 223 if (petsc_logObjects) { 224 ierr = PetscMalloc(petsc_maxObjects * sizeof(Object), &petsc_objects);CHKERRQ(ierr); 225 } 226 PetscLogPHC = PetscLogObjCreateDefault; 227 PetscLogPHD = PetscLogObjDestroyDefault; 228 /* Setup default logging structures */ 229 ierr = PetscStageLogCreate(&petsc_stageLog);CHKERRQ(ierr); 230 ierr = PetscStageLogRegister(petsc_stageLog, "Main Stage", &stage);CHKERRQ(ierr); 231 #if defined(PETSC_HAVE_CHUD) 232 ierr = chudInitialize();CHKERRQ(ierr); 233 ierr = chudAcquireSamplingFacility(CHUD_BLOCKING);CHKERRQ(ierr); 234 ierr = chudSetSamplingDevice(chudCPU1Dev);CHKERRQ(ierr); 235 ierr = chudSetStartDelay(0,chudNanoSeconds);CHKERRQ(ierr); 236 ierr = chudClearPMCMode(chudCPU1Dev,chudUnused);CHKERRQ(ierr); 237 ierr = chudClearPMCs();CHKERRQ(ierr); 238 /* ierr = chudSetPMCMuxPosition(chudCPU1Dev,0,0);CHKERRQ(ierr); */ 239 printf("%s\n",chudGetEventName(chudCPU1Dev,PMC_1,193)); 240 printf("%s\n",chudGetEventDescription(chudCPU1Dev,PMC_1,193)); 241 printf("%s\n",chudGetEventNotes(chudCPU1Dev,PMC_1,193)); 242 ierr = chudSetPMCEvent(chudCPU1Dev,PMC_1,193);CHKERRQ(ierr); 243 ierr = chudSetPMCMode(chudCPU1Dev,PMC_1,chudCounter);CHKERRQ(ierr); 244 ierr = chudSetPrivilegeFilter(chudCPU1Dev,PMC_1,chudCountUserEvents);CHKERRQ(ierr); 245 ierr = chudSetPMCEventMask(chudCPU1Dev,PMC_1,0xFE);CHKERRQ(ierr); 246 if (!chudIsEventValid(chudCPU1Dev,PMC_1,193)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_SUP,"Event is not valid %d",193); 247 ierr = chudStartPMCs();CHKERRQ(ierr); 248 #endif 249 #if defined(PETSC_HAVE_PAPI) 250 ierr = PAPI_library_init(PAPI_VER_CURRENT); 251 if (ierr != PAPI_VER_CURRENT) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Cannot initialize PAPI"); 252 ierr = PAPI_query_event(PAPI_FP_INS);CHKERRQ(ierr); 253 ierr = PAPI_create_eventset(&PAPIEventSet);CHKERRQ(ierr); 254 ierr = PAPI_add_event(PAPIEventSet,PAPI_FP_INS);CHKERRQ(ierr); 255 ierr = PAPI_start(PAPIEventSet);CHKERRQ(ierr); 256 #endif 257 258 /* All processors sync here for more consistent logging */ 259 ierr = MPI_Barrier(PETSC_COMM_WORLD);CHKERRQ(ierr); 260 PetscTime(&petsc_BaseTime); 261 ierr = PetscLogStagePush(stage);CHKERRQ(ierr); 262 PetscFunctionReturn(0); 263 } 264 265 #undef __FUNCT__ 266 #define __FUNCT__ "PetscLogBegin" 267 /*@C 268 PetscLogBegin - Turns on logging of objects and events. This logs flop 269 rates and object creation and should not slow programs down too much. 270 This routine may be called more than once. 271 272 Logically Collective over PETSC_COMM_WORLD 273 274 Options Database Keys: 275 + -log_summary - Prints summary of flop and timing information to the 276 screen (for code compiled with PETSC_USE_LOG) 277 - -log - Prints detailed log information (for code compiled with PETSC_USE_LOG) 278 279 Usage: 280 .vb 281 PetscInitialize(...); 282 PetscLogBegin(); 283 ... code ... 284 PetscLogView(viewer); or PetscLogDump(); 285 PetscFinalize(); 286 .ve 287 288 Notes: 289 PetscLogView(viewer) or PetscLogDump() actually cause the printing of 290 the logging information. 291 292 Level: advanced 293 294 .keywords: log, begin 295 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogTraceBegin() 296 @*/ 297 PetscErrorCode PetscLogBegin(void) 298 { 299 PetscErrorCode ierr; 300 301 PetscFunctionBegin; 302 ierr = PetscLogSet(PetscLogEventBeginDefault, PetscLogEventEndDefault);CHKERRQ(ierr); 303 ierr = PetscLogBegin_Private();CHKERRQ(ierr); 304 PetscFunctionReturn(0); 305 } 306 307 #undef __FUNCT__ 308 #define __FUNCT__ "PetscLogAllBegin" 309 /*@C 310 PetscLogAllBegin - Turns on extensive logging of objects and events. Logs 311 all events. This creates large log files and slows the program down. 312 313 Logically Collective on PETSC_COMM_WORLD 314 315 Options Database Keys: 316 . -log_all - Prints extensive log information (for code compiled with PETSC_USE_LOG) 317 318 Usage: 319 .vb 320 PetscInitialize(...); 321 PetscLogAllBegin(); 322 ... code ... 323 PetscLogDump(filename); 324 PetscFinalize(); 325 .ve 326 327 Notes: 328 A related routine is PetscLogBegin() (with the options key -log), which is 329 intended for production runs since it logs only flop rates and object 330 creation (and shouldn't significantly slow the programs). 331 332 Level: advanced 333 334 .keywords: log, all, begin 335 .seealso: PetscLogDump(), PetscLogBegin(), PetscLogTraceBegin() 336 @*/ 337 PetscErrorCode PetscLogAllBegin(void) 338 { 339 PetscErrorCode ierr; 340 341 PetscFunctionBegin; 342 ierr = PetscLogSet(PetscLogEventBeginComplete, PetscLogEventEndComplete);CHKERRQ(ierr); 343 ierr = PetscLogBegin_Private();CHKERRQ(ierr); 344 PetscFunctionReturn(0); 345 } 346 347 #undef __FUNCT__ 348 #define __FUNCT__ "PetscLogTraceBegin" 349 /*@ 350 PetscLogTraceBegin - Activates trace logging. Every time a PETSc event 351 begins or ends, the event name is printed. 352 353 Logically Collective on PETSC_COMM_WORLD 354 355 Input Parameter: 356 . file - The file to print trace in (e.g. stdout) 357 358 Options Database Key: 359 . -log_trace [filename] - Activates PetscLogTraceBegin() 360 361 Notes: 362 PetscLogTraceBegin() prints the processor number, the execution time (sec), 363 then "Event begin:" or "Event end:" followed by the event name. 364 365 PetscLogTraceBegin() allows tracing of all PETSc calls, which is useful 366 to determine where a program is hanging without running in the 367 debugger. Can be used in conjunction with the -info option. 368 369 Level: intermediate 370 371 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogView(), PetscLogBegin() 372 @*/ 373 PetscErrorCode PetscLogTraceBegin(FILE *file) 374 { 375 PetscErrorCode ierr; 376 377 PetscFunctionBegin; 378 petsc_tracefile = file; 379 380 ierr = PetscLogSet(PetscLogEventBeginTrace, PetscLogEventEndTrace);CHKERRQ(ierr); 381 ierr = PetscLogBegin_Private();CHKERRQ(ierr); 382 PetscFunctionReturn(0); 383 } 384 385 #undef __FUNCT__ 386 #define __FUNCT__ "PetscLogActions" 387 /*@ 388 PetscLogActions - Determines whether actions are logged for the graphical viewer. 389 390 Not Collective 391 392 Input Parameter: 393 . flag - PETSC_TRUE if actions are to be logged 394 395 Level: intermediate 396 397 Note: Logging of actions continues to consume more memory as the program 398 runs. Long running programs should consider turning this feature off. 399 400 Options Database Keys: 401 . -log_exclude_actions - Turns off actions logging 402 403 .keywords: log, stage, register 404 .seealso: PetscLogStagePush(), PetscLogStagePop() 405 @*/ 406 PetscErrorCode PetscLogActions(PetscBool flag) 407 { 408 PetscFunctionBegin; 409 petsc_logActions = flag; 410 PetscFunctionReturn(0); 411 } 412 413 #undef __FUNCT__ 414 #define __FUNCT__ "PetscLogObjects" 415 /*@ 416 PetscLogObjects - Determines whether objects are logged for the graphical viewer. 417 418 Not Collective 419 420 Input Parameter: 421 . flag - PETSC_TRUE if objects are to be logged 422 423 Level: intermediate 424 425 Note: Logging of objects continues to consume more memory as the program 426 runs. Long running programs should consider turning this feature off. 427 428 Options Database Keys: 429 . -log_exclude_objects - Turns off objects logging 430 431 .keywords: log, stage, register 432 .seealso: PetscLogStagePush(), PetscLogStagePop() 433 @*/ 434 PetscErrorCode PetscLogObjects(PetscBool flag) 435 { 436 PetscFunctionBegin; 437 petsc_logObjects = flag; 438 PetscFunctionReturn(0); 439 } 440 441 /*------------------------------------------------ Stage Functions --------------------------------------------------*/ 442 #undef __FUNCT__ 443 #define __FUNCT__ "PetscLogStageRegister" 444 /*@C 445 PetscLogStageRegister - Attaches a charactor string name to a logging stage. 446 447 Not Collective 448 449 Input Parameter: 450 . sname - The name to associate with that stage 451 452 Output Parameter: 453 . stage - The stage number 454 455 Level: intermediate 456 457 .keywords: log, stage, register 458 .seealso: PetscLogStagePush(), PetscLogStagePop() 459 @*/ 460 PetscErrorCode PetscLogStageRegister(const char sname[],PetscLogStage *stage) 461 { 462 PetscStageLog stageLog; 463 PetscLogEvent event; 464 PetscErrorCode ierr; 465 466 PetscFunctionBegin; 467 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 468 ierr = PetscStageLogRegister(stageLog, sname, stage);CHKERRQ(ierr); 469 /* Copy events already changed in the main stage, this sucks */ 470 ierr = EventPerfLogEnsureSize(stageLog->stageInfo[*stage].eventLog, stageLog->eventLog->numEvents);CHKERRQ(ierr); 471 for (event = 0; event < stageLog->eventLog->numEvents; event++) { 472 ierr = EventPerfInfoCopy(&stageLog->stageInfo[0].eventLog->eventInfo[event],&stageLog->stageInfo[*stage].eventLog->eventInfo[event]);CHKERRQ(ierr); 473 } 474 ierr = ClassPerfLogEnsureSize(stageLog->stageInfo[*stage].classLog, stageLog->classLog->numClasses);CHKERRQ(ierr); 475 PetscFunctionReturn(0); 476 } 477 478 #undef __FUNCT__ 479 #define __FUNCT__ "PetscLogStagePush" 480 /*@C 481 PetscLogStagePush - This function pushes a stage on the stack. 482 483 Not Collective 484 485 Input Parameter: 486 . stage - The stage on which to log 487 488 Usage: 489 If the option -log_sumary is used to run the program containing the 490 following code, then 2 sets of summary data will be printed during 491 PetscFinalize(). 492 .vb 493 PetscInitialize(int *argc,char ***args,0,0); 494 [stage 0 of code] 495 PetscLogStagePush(1); 496 [stage 1 of code] 497 PetscLogStagePop(); 498 PetscBarrier(...); 499 [more stage 0 of code] 500 PetscFinalize(); 501 .ve 502 503 Notes: 504 Use PetscLogStageRegister() to register a stage. 505 506 Level: intermediate 507 508 .keywords: log, push, stage 509 .seealso: PetscLogStagePop(), PetscLogStageRegister(), PetscBarrier() 510 @*/ 511 PetscErrorCode PetscLogStagePush(PetscLogStage stage) 512 { 513 PetscStageLog stageLog; 514 PetscErrorCode ierr; 515 516 PetscFunctionBegin; 517 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 518 ierr = PetscStageLogPush(stageLog, stage);CHKERRQ(ierr); 519 PetscFunctionReturn(0); 520 } 521 522 #undef __FUNCT__ 523 #define __FUNCT__ "PetscLogStagePop" 524 /*@C 525 PetscLogStagePop - This function pops a stage from the stack. 526 527 Not Collective 528 529 Usage: 530 If the option -log_sumary is used to run the program containing the 531 following code, then 2 sets of summary data will be printed during 532 PetscFinalize(). 533 .vb 534 PetscInitialize(int *argc,char ***args,0,0); 535 [stage 0 of code] 536 PetscLogStagePush(1); 537 [stage 1 of code] 538 PetscLogStagePop(); 539 PetscBarrier(...); 540 [more stage 0 of code] 541 PetscFinalize(); 542 .ve 543 544 Notes: 545 Use PetscLogStageRegister() to register a stage. 546 547 Level: intermediate 548 549 .keywords: log, pop, stage 550 .seealso: PetscLogStagePush(), PetscLogStageRegister(), PetscBarrier() 551 @*/ 552 PetscErrorCode PetscLogStagePop(void) 553 { 554 PetscStageLog stageLog; 555 PetscErrorCode ierr; 556 557 PetscFunctionBegin; 558 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 559 ierr = PetscStageLogPop(stageLog);CHKERRQ(ierr); 560 PetscFunctionReturn(0); 561 } 562 563 #undef __FUNCT__ 564 #define __FUNCT__ "PetscLogStageSetActive" 565 /*@ 566 PetscLogStageSetActive - Determines stage activity for PetscLogEventBegin() and PetscLogEventEnd(). 567 568 Not Collective 569 570 Input Parameters: 571 + stage - The stage 572 - isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE) 573 574 Level: intermediate 575 576 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 577 @*/ 578 PetscErrorCode PetscLogStageSetActive(PetscLogStage stage, PetscBool isActive) 579 { 580 PetscStageLog stageLog; 581 PetscErrorCode ierr; 582 583 PetscFunctionBegin; 584 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 585 ierr = PetscStageLogSetActive(stageLog, stage, isActive);CHKERRQ(ierr); 586 PetscFunctionReturn(0); 587 } 588 589 #undef __FUNCT__ 590 #define __FUNCT__ "PetscLogStageGetActive" 591 /*@ 592 PetscLogStageGetActive - Returns stage activity for PetscLogEventBegin() and PetscLogEventEnd(). 593 594 Not Collective 595 596 Input Parameter: 597 . stage - The stage 598 599 Output Parameter: 600 . isActive - The activity flag, PETSC_TRUE for logging, else PETSC_FALSE (defaults to PETSC_TRUE) 601 602 Level: intermediate 603 604 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 605 @*/ 606 PetscErrorCode PetscLogStageGetActive(PetscLogStage stage, PetscBool *isActive) 607 { 608 PetscStageLog stageLog; 609 PetscErrorCode ierr; 610 611 PetscFunctionBegin; 612 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 613 ierr = PetscStageLogGetActive(stageLog, stage, isActive);CHKERRQ(ierr); 614 PetscFunctionReturn(0); 615 } 616 617 #undef __FUNCT__ 618 #define __FUNCT__ "PetscLogStageSetVisible" 619 /*@ 620 PetscLogStageSetVisible - Determines stage visibility in PetscLogView() 621 622 Not Collective 623 624 Input Parameters: 625 + stage - The stage 626 - isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE) 627 628 Level: intermediate 629 630 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView() 631 @*/ 632 PetscErrorCode PetscLogStageSetVisible(PetscLogStage stage, PetscBool isVisible) 633 { 634 PetscStageLog stageLog; 635 PetscErrorCode ierr; 636 637 PetscFunctionBegin; 638 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 639 ierr = PetscStageLogSetVisible(stageLog, stage, isVisible);CHKERRQ(ierr); 640 PetscFunctionReturn(0); 641 } 642 643 #undef __FUNCT__ 644 #define __FUNCT__ "PetscLogStageGetVisible" 645 /*@ 646 PetscLogStageGetVisible - Returns stage visibility in PetscLogView() 647 648 Not Collective 649 650 Input Parameter: 651 . stage - The stage 652 653 Output Parameter: 654 . isVisible - The visibility flag, PETSC_TRUE to print, else PETSC_FALSE (defaults to PETSC_TRUE) 655 656 Level: intermediate 657 658 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscLogView() 659 @*/ 660 PetscErrorCode PetscLogStageGetVisible(PetscLogStage stage, PetscBool *isVisible) 661 { 662 PetscStageLog stageLog; 663 PetscErrorCode ierr; 664 665 PetscFunctionBegin; 666 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 667 ierr = PetscStageLogGetVisible(stageLog, stage, isVisible);CHKERRQ(ierr); 668 PetscFunctionReturn(0); 669 } 670 671 #undef __FUNCT__ 672 #define __FUNCT__ "PetscLogStageGetId" 673 /*@C 674 PetscLogStageGetId - Returns the stage id when given the stage name. 675 676 Not Collective 677 678 Input Parameter: 679 . name - The stage name 680 681 Output Parameter: 682 . stage - The stage 683 684 Level: intermediate 685 686 .seealso: PetscLogStagePush(), PetscLogStagePop(), PetscPreLoadBegin(), PetscPreLoadEnd(), PetscPreLoadStage() 687 @*/ 688 PetscErrorCode PetscLogStageGetId(const char name[], PetscLogStage *stage) 689 { 690 PetscStageLog stageLog; 691 PetscErrorCode ierr; 692 693 PetscFunctionBegin; 694 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 695 ierr = PetscStageLogGetStage(stageLog, name, stage);CHKERRQ(ierr); 696 PetscFunctionReturn(0); 697 } 698 699 /*------------------------------------------------ Event Functions --------------------------------------------------*/ 700 #undef __FUNCT__ 701 #define __FUNCT__ "PetscLogEventRegister" 702 /*@C 703 PetscLogEventRegister - Registers an event name for logging operations in an application code. 704 705 Not Collective 706 707 Input Parameter: 708 + name - The name associated with the event 709 - classid - The classid associated to the class for this event, obtain either with 710 PetscClassIdRegister() or use a predefined one such as KSP_CLASSID, SNES_CLASSID, the predefined ones 711 are only available in C code 712 713 Output Parameter: 714 . event - The event id for use with PetscLogEventBegin() and PetscLogEventEnd(). 715 716 Example of Usage: 717 .vb 718 PetscLogEvent USER_EVENT; 719 PetscClassId classid; 720 PetscLogDouble user_event_flops; 721 PetscClassIdRegister("class name",&classid); 722 PetscLogEventRegister("User event name",classid,&USER_EVENT); 723 PetscLogEventBegin(USER_EVENT,0,0,0,0); 724 [code segment to monitor] 725 PetscLogFlops(user_event_flops); 726 PetscLogEventEnd(USER_EVENT,0,0,0,0); 727 .ve 728 729 Notes: 730 PETSc automatically logs library events if the code has been 731 compiled with -DPETSC_USE_LOG (which is the default) and -log, 732 -log_summary, or -log_all are specified. PetscLogEventRegister() is 733 intended for logging user events to supplement this PETSc 734 information. 735 736 PETSc can gather data for use with the utilities Jumpshot 737 (part of the MPICH distribution). If PETSc has been compiled 738 with flag -DPETSC_HAVE_MPE (MPE is an additional utility within 739 MPICH), the user can employ another command line option, -log_mpe, 740 to create a logfile, "mpe.log", which can be visualized 741 Jumpshot. 742 743 The classid is associated with each event so that classes of events 744 can be disabled simultaneously, such as all matrix events. The user 745 can either use an existing classid, such as MAT_CLASSID, or create 746 their own as shown in the example. 747 748 Level: intermediate 749 750 .keywords: log, event, register 751 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogFlops(), 752 PetscLogEventMPEActivate(), PetscLogEventMPEDeactivate(), 753 PetscLogEventActivate(), PetscLogEventDeactivate(), PetscClassIdRegister() 754 @*/ 755 PetscErrorCode PetscLogEventRegister(const char name[],PetscClassId classid,PetscLogEvent *event) 756 { 757 PetscStageLog stageLog; 758 int stage; 759 PetscErrorCode ierr; 760 761 PetscFunctionBegin; 762 *event = PETSC_DECIDE; 763 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 764 ierr = EventRegLogRegister(stageLog->eventLog, name, classid, event);CHKERRQ(ierr); 765 for (stage = 0; stage < stageLog->numStages; stage++) { 766 ierr = EventPerfLogEnsureSize(stageLog->stageInfo[stage].eventLog, stageLog->eventLog->numEvents);CHKERRQ(ierr); 767 ierr = ClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses);CHKERRQ(ierr); 768 } 769 PetscFunctionReturn(0); 770 } 771 772 #undef __FUNCT__ 773 #define __FUNCT__ "PetscLogEventActivate" 774 /*@ 775 PetscLogEventActivate - Indicates that a particular event should be logged. 776 777 Not Collective 778 779 Input Parameter: 780 . event - The event id 781 782 Usage: 783 .vb 784 PetscLogEventDeactivate(VEC_SetValues); 785 [code where you do not want to log VecSetValues()] 786 PetscLogEventActivate(VEC_SetValues); 787 [code where you do want to log VecSetValues()] 788 .ve 789 790 Note: 791 The event may be either a pre-defined PETSc event (found in include/petsclog.h) 792 or an event number obtained with PetscLogEventRegister(). 793 794 Level: advanced 795 796 .keywords: log, event, activate 797 .seealso: PetscLogEventMPEDeactivate(),PetscLogEventMPEActivate(),PlogEventDeactivate() 798 @*/ 799 PetscErrorCode PetscLogEventActivate(PetscLogEvent event) 800 { 801 PetscStageLog stageLog; 802 int stage; 803 PetscErrorCode ierr; 804 805 PetscFunctionBegin; 806 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 807 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 808 ierr = EventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 809 PetscFunctionReturn(0); 810 } 811 812 #undef __FUNCT__ 813 #define __FUNCT__ "PetscLogEventDeactivate" 814 /*@ 815 PetscLogEventDeactivate - Indicates that a particular event should not be logged. 816 817 Not Collective 818 819 Input Parameter: 820 . event - The event id 821 822 Usage: 823 .vb 824 PetscLogEventDeactivate(VEC_SetValues); 825 [code where you do not want to log VecSetValues()] 826 PetscLogEventActivate(VEC_SetValues); 827 [code where you do want to log VecSetValues()] 828 .ve 829 830 Note: 831 The event may be either a pre-defined PETSc event (found in 832 include/petsclog.h) or an event number obtained with PetscLogEventRegister()). 833 834 Level: advanced 835 836 .keywords: log, event, deactivate 837 .seealso: PetscLogEventMPEDeactivate(),PetscLogEventMPEActivate(),PlogEventActivate() 838 @*/ 839 PetscErrorCode PetscLogEventDeactivate(PetscLogEvent event) 840 { 841 PetscStageLog stageLog; 842 int stage; 843 PetscErrorCode ierr; 844 845 PetscFunctionBegin; 846 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 847 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 848 ierr = EventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 849 PetscFunctionReturn(0); 850 } 851 852 #undef __FUNCT__ 853 #define __FUNCT__ "PetscLogEventSetActiveAll" 854 /*@ 855 PetscLogEventSetActiveAll - Sets the event activity in every stage. 856 857 Not Collective 858 859 Input Parameters: 860 + event - The event id 861 - isActive - The activity flag determining whether the event is logged 862 863 Level: advanced 864 865 .keywords: log, event, activate 866 .seealso: PetscLogEventMPEDeactivate(),PetscLogEventMPEActivate(),PlogEventActivate(),PlogEventDeactivate() 867 @*/ 868 PetscErrorCode PetscLogEventSetActiveAll(PetscLogEvent event, PetscBool isActive) 869 { 870 PetscStageLog stageLog; 871 int stage; 872 PetscErrorCode ierr; 873 874 PetscFunctionBegin; 875 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 876 for (stage = 0; stage < stageLog->numStages; stage++) { 877 if (isActive) { 878 ierr = EventPerfLogActivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 879 } else { 880 ierr = EventPerfLogDeactivate(stageLog->stageInfo[stage].eventLog, event);CHKERRQ(ierr); 881 } 882 } 883 PetscFunctionReturn(0); 884 } 885 886 #undef __FUNCT__ 887 #define __FUNCT__ "PetscLogEventActivateClass" 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 .keywords: log, event, activate, class 899 .seealso: PetscInfoActivate(),PetscInfo(),PetscInfoAllow(),PetscLogEventDeactivateClass(), PetscLogEventActivate(),PetscLogEventDeactivate() 900 @*/ 901 PetscErrorCode PetscLogEventActivateClass(PetscClassId classid) 902 { 903 PetscStageLog stageLog; 904 int stage; 905 PetscErrorCode ierr; 906 907 PetscFunctionBegin; 908 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 909 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 910 ierr = EventPerfLogActivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);CHKERRQ(ierr); 911 PetscFunctionReturn(0); 912 } 913 914 #undef __FUNCT__ 915 #define __FUNCT__ "PetscLogEventDeactivateClass" 916 /*@ 917 PetscLogEventDeactivateClass - Deactivates event logging for a PETSc object class. 918 919 Not Collective 920 921 Input Parameter: 922 . classid - The event class, for example MAT_CLASSID, SNES_CLASSID, etc. 923 924 Level: developer 925 926 .keywords: log, event, deactivate, class 927 .seealso: PetscInfoActivate(),PetscInfo(),PetscInfoAllow(),PetscLogEventActivateClass(), PetscLogEventActivate(),PetscLogEventDeactivate() 928 @*/ 929 PetscErrorCode PetscLogEventDeactivateClass(PetscClassId classid) 930 { 931 PetscStageLog stageLog; 932 int stage; 933 PetscErrorCode ierr; 934 935 PetscFunctionBegin; 936 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 937 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 938 ierr = EventPerfLogDeactivateClass(stageLog->stageInfo[stage].eventLog, stageLog->eventLog, classid);CHKERRQ(ierr); 939 PetscFunctionReturn(0); 940 } 941 942 /*MC 943 PetscLogEventBegin - Logs the beginning of a user event. 944 945 Synopsis: 946 #include "petsclog.h" 947 PetscErrorCode PetscLogEventBegin(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 948 949 Not Collective 950 951 Input Parameters: 952 + e - integer associated with the event obtained from PetscLogEventRegister() 953 - o1,o2,o3,o4 - objects associated with the event, or 0 954 955 956 Fortran Synopsis: 957 void PetscLogEventBegin(int e,PetscErrorCode ierr) 958 959 Usage: 960 .vb 961 PetscLogEvent USER_EVENT; 962 PetscLogDouble user_event_flops; 963 PetscLogEventRegister("User event",0,&USER_EVENT); 964 PetscLogEventBegin(USER_EVENT,0,0,0,0); 965 [code segment to monitor] 966 PetscLogFlops(user_event_flops); 967 PetscLogEventEnd(USER_EVENT,0,0,0,0); 968 .ve 969 970 Notes: 971 You need to register each integer event with the command 972 PetscLogEventRegister(). The source code must be compiled with 973 -DPETSC_USE_LOG, which is the default. 974 975 PETSc automatically logs library events if the code has been 976 compiled with -DPETSC_USE_LOG, and -log, -log_summary, or -log_all are 977 specified. PetscLogEventBegin() is intended for logging user events 978 to supplement this PETSc information. 979 980 Level: intermediate 981 982 .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops() 983 984 .keywords: log, event, begin 985 M*/ 986 987 /*MC 988 PetscLogEventEnd - Log the end of a user event. 989 990 Synopsis: 991 #include "petsclog.h" 992 PetscErrorCode PetscLogEventEnd(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4) 993 994 Not Collective 995 996 Input Parameters: 997 + e - integer associated with the event obtained with PetscLogEventRegister() 998 - o1,o2,o3,o4 - objects associated with the event, or 0 999 1000 1001 Fortran Synopsis: 1002 void PetscLogEventEnd(int e,PetscErrorCode ierr) 1003 1004 Usage: 1005 .vb 1006 PetscLogEvent USER_EVENT; 1007 PetscLogDouble user_event_flops; 1008 PetscLogEventRegister("User event",0,&USER_EVENT,); 1009 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1010 [code segment to monitor] 1011 PetscLogFlops(user_event_flops); 1012 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1013 .ve 1014 1015 Notes: 1016 You should also register each additional integer event with the command 1017 PetscLogEventRegister(). Source code must be compiled with 1018 -DPETSC_USE_LOG, which is the default. 1019 1020 PETSc automatically logs library events if the code has been 1021 compiled with -DPETSC_USE_LOG, and -log, -log_summary, or -log_all are 1022 specified. PetscLogEventEnd() is intended for logging user events 1023 to supplement this PETSc information. 1024 1025 Level: intermediate 1026 1027 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogFlops() 1028 1029 .keywords: log, event, end 1030 M*/ 1031 1032 /*MC 1033 PetscLogEventBarrierBegin - Logs the time in a barrier before an event. 1034 1035 Synopsis: 1036 #include "petsclog.h" 1037 PetscErrorCode PetscLogEventBarrierBegin(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4,MPI_Comm comm) 1038 1039 Not Collective 1040 1041 Input Parameters: 1042 . e - integer associated with the event obtained from PetscLogEventRegister() 1043 . o1,o2,o3,o4 - objects associated with the event, or 0 1044 . comm - communicator the barrier takes place over 1045 1046 1047 Usage: 1048 .vb 1049 PetscLogEventBarrierBegin(VEC_NormBarrier,0,0,0,0,comm); 1050 MPI_Allreduce() 1051 PetscLogEventBarrierEnd(VEC_NormBarrier,0,0,0,0,comm); 1052 .ve 1053 1054 Notes: 1055 This is for logging the amount of time spent in a barrier for an event 1056 that requires synchronization. 1057 1058 Additional Notes: 1059 Synchronization events always come in pairs; for example, VEC_NormBarrier and 1060 VEC_NormComm = VEC_NormBarrier + 1 1061 1062 Level: advanced 1063 1064 .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops(), PetscLogEventBegin(), 1065 PetscLogEventBarrierEnd() 1066 1067 .keywords: log, event, begin, barrier 1068 M*/ 1069 1070 /*MC 1071 PetscLogEventBarrierEnd - Logs the time in a barrier before an event. 1072 1073 Synopsis: 1074 #include "petsclog.h" 1075 PetscErrorCode PetscLogEventBarrierEnd(int e,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4,MPI_Comm comm) 1076 1077 Logically Collective on MPI_Comm 1078 1079 Input Parameters: 1080 . e - integer associated with the event obtained from PetscLogEventRegister() 1081 . o1,o2,o3,o4 - objects associated with the event, or 0 1082 . comm - communicator the barrier takes place over 1083 1084 1085 Usage: 1086 .vb 1087 PetscLogEventBarrierBegin(VEC_NormBarrier,0,0,0,0,comm); 1088 MPI_Allreduce() 1089 PetscLogEventBarrierEnd(VEC_NormBarrier,0,0,0,0,comm); 1090 .ve 1091 1092 Notes: 1093 This is for logging the amount of time spent in a barrier for an event 1094 that requires synchronization. 1095 1096 Additional Notes: 1097 Synchronization events always come in pairs; for example, VEC_NormBarrier and 1098 VEC_NormComm = VEC_NormBarrier + 1 1099 1100 Level: advanced 1101 1102 .seealso: PetscLogEventRegister(), PetscLogEventEnd(), PetscLogFlops(), PetscLogEventBegin(), 1103 PetscLogEventBarrierBegin() 1104 1105 .keywords: log, event, begin, barrier 1106 M*/ 1107 1108 #undef __FUNCT__ 1109 #define __FUNCT__ "PetscLogEventGetId" 1110 /*@C 1111 PetscLogEventGetId - Returns the event id when given the event name. 1112 1113 Not Collective 1114 1115 Input Parameter: 1116 . name - The event name 1117 1118 Output Parameter: 1119 . event - The event 1120 1121 Level: intermediate 1122 1123 .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogStageGetId() 1124 @*/ 1125 PetscErrorCode PetscLogEventGetId(const char name[], PetscLogEvent *event) 1126 { 1127 PetscStageLog stageLog; 1128 PetscErrorCode ierr; 1129 1130 PetscFunctionBegin; 1131 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1132 ierr = EventRegLogGetEvent(stageLog->eventLog, name, event);CHKERRQ(ierr); 1133 PetscFunctionReturn(0); 1134 } 1135 1136 1137 /*------------------------------------------------ Output Functions -------------------------------------------------*/ 1138 #undef __FUNCT__ 1139 #define __FUNCT__ "PetscLogDump" 1140 /*@C 1141 PetscLogDump - Dumps logs of objects to a file. This file is intended to 1142 be read by bin/petscview. This program no longer exists. 1143 1144 Collective on PETSC_COMM_WORLD 1145 1146 Input Parameter: 1147 . name - an optional file name 1148 1149 Options Database Keys: 1150 + -log - Prints basic log information (for code compiled with PETSC_USE_LOG) 1151 - -log_all - Prints extensive log information (for code compiled with PETSC_USE_LOG) 1152 1153 Usage: 1154 .vb 1155 PetscInitialize(...); 1156 PetscLogBegin(); or PetscLogAllBegin(); 1157 ... code ... 1158 PetscLogDump(filename); 1159 PetscFinalize(); 1160 .ve 1161 1162 Notes: 1163 The default file name is 1164 $ Log.<rank> 1165 where <rank> is the processor number. If no name is specified, 1166 this file will be used. 1167 1168 Level: advanced 1169 1170 .keywords: log, dump 1171 .seealso: PetscLogBegin(), PetscLogAllBegin(), PetscLogView() 1172 @*/ 1173 PetscErrorCode PetscLogDump(const char sname[]) 1174 { 1175 PetscStageLog stageLog; 1176 PetscEventPerfInfo *eventInfo; 1177 FILE *fd; 1178 char file[PETSC_MAX_PATH_LEN], fname[PETSC_MAX_PATH_LEN]; 1179 PetscLogDouble flops, _TotalTime; 1180 PetscMPIInt rank; 1181 int action, object, curStage; 1182 PetscLogEvent event; 1183 PetscErrorCode ierr; 1184 1185 PetscFunctionBegin; 1186 /* Calculate the total elapsed time */ 1187 PetscTime(&_TotalTime); 1188 _TotalTime -= petsc_BaseTime; 1189 /* Open log file */ 1190 ierr = MPI_Comm_rank(PETSC_COMM_WORLD, &rank);CHKERRQ(ierr); 1191 if (sname) sprintf(file, "%s.%d", sname, rank); 1192 else sprintf(file, "Log.%d", rank); 1193 ierr = PetscFixFilename(file, fname);CHKERRQ(ierr); 1194 ierr = PetscFOpen(PETSC_COMM_WORLD, fname, "w", &fd);CHKERRQ(ierr); 1195 if ((!rank) && (!fd)) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_FILE_OPEN, "Cannot open file: %s", fname); 1196 /* Output totals */ 1197 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Total Flops %14e %16.8e\n", petsc_TotalFlops, _TotalTime);CHKERRQ(ierr); 1198 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Clock Resolution %g\n", 0.0);CHKERRQ(ierr); 1199 /* Output actions */ 1200 if (petsc_logActions) { 1201 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Actions accomplished %d\n", petsc_numActions);CHKERRQ(ierr); 1202 for (action = 0; action < petsc_numActions; action++) { 1203 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "%g %d %d %d %d %d %d %g %g %g\n", 1204 petsc_actions[action].time, petsc_actions[action].action, (int)petsc_actions[action].event, (int)petsc_actions[action].classid, petsc_actions[action].id1, 1205 petsc_actions[action].id2, petsc_actions[action].id3, petsc_actions[action].flops, petsc_actions[action].mem, petsc_actions[action].maxmem);CHKERRQ(ierr); 1206 } 1207 } 1208 /* Output objects */ 1209 if (petsc_logObjects) { 1210 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Objects created %d destroyed %d\n", petsc_numObjects, petsc_numObjectsDestroyed);CHKERRQ(ierr); 1211 for (object = 0; object < petsc_numObjects; object++) { 1212 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Parent ID: %d Memory: %d\n", petsc_objects[object].parent, (int) petsc_objects[object].mem);CHKERRQ(ierr); 1213 if (!petsc_objects[object].name[0]) { 1214 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd,"No Name\n");CHKERRQ(ierr); 1215 } else { 1216 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Name: %s\n", petsc_objects[object].name);CHKERRQ(ierr); 1217 } 1218 if (petsc_objects[object].info[0] != 0) { 1219 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "No Info\n");CHKERRQ(ierr); 1220 } else { 1221 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Info: %s\n", petsc_objects[object].info);CHKERRQ(ierr); 1222 } 1223 } 1224 } 1225 /* Output events */ 1226 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "Event log:\n");CHKERRQ(ierr); 1227 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1228 ierr = PetscIntStackTop(stageLog->stack, &curStage);CHKERRQ(ierr); 1229 eventInfo = stageLog->stageInfo[curStage].eventLog->eventInfo; 1230 for (event = 0; event < stageLog->stageInfo[curStage].eventLog->numEvents; event++) { 1231 if (eventInfo[event].time != 0.0) flops = eventInfo[event].flops/eventInfo[event].time; 1232 else flops = 0.0; 1233 ierr = PetscFPrintf(PETSC_COMM_WORLD, fd, "%d %16d %16g %16g %16g\n", event, eventInfo[event].count, 1234 eventInfo[event].flops, eventInfo[event].time, flops);CHKERRQ(ierr); 1235 } 1236 ierr = PetscFClose(PETSC_COMM_WORLD, fd);CHKERRQ(ierr); 1237 PetscFunctionReturn(0); 1238 } 1239 1240 #undef __FUNCT__ 1241 #define __FUNCT__ "PetscLogView" 1242 /*@C 1243 PetscLogViewer - Prints a summary of the logging. 1244 1245 Collective over MPI_Comm 1246 1247 Input Parameter: 1248 . viewer - an ASCII viewer 1249 1250 Options Database Keys: 1251 . -log_summary - Prints summary of log information (for code compiled with PETSC_USE_LOG) 1252 1253 Usage: 1254 .vb 1255 PetscInitialize(...); 1256 PetscLogBegin(); 1257 ... code ... 1258 PetscLogView(PetscViewer); 1259 PetscFinalize(...); 1260 .ve 1261 1262 Notes: 1263 By default the summary is printed to stdout. 1264 1265 Level: beginner 1266 1267 .keywords: log, dump, print 1268 .seealso: PetscLogBegin(), PetscLogDump() 1269 @*/ 1270 PetscErrorCode PetscLogView(PetscViewer viewer) 1271 { 1272 FILE *fd; 1273 PetscLogDouble zero = 0.0; 1274 PetscStageLog stageLog; 1275 PetscStageInfo *stageInfo = NULL; 1276 PetscEventPerfInfo *eventInfo = NULL; 1277 PetscClassPerfInfo *classInfo; 1278 char arch[128],hostname[128],username[128],pname[PETSC_MAX_PATH_LEN],date[128]; 1279 const char *name; 1280 PetscLogDouble locTotalTime, TotalTime, TotalFlops; 1281 PetscLogDouble numMessages, messageLength, avgMessLen, numReductions; 1282 PetscLogDouble stageTime, flops, flopr, mem, mess, messLen, red; 1283 PetscLogDouble fracTime, fracFlops, fracMessages, fracLength, fracReductions, fracMess, fracMessLen, fracRed; 1284 PetscLogDouble fracStageTime, fracStageFlops, fracStageMess, fracStageMessLen, fracStageRed; 1285 PetscLogDouble min, max, tot, ratio, avg, x, y; 1286 PetscLogDouble minf, maxf, totf, ratf, mint, maxt, tott, ratt, ratCt, totm, totml, totr; 1287 PetscMPIInt minCt, maxCt; 1288 PetscMPIInt size, rank; 1289 PetscBool *localStageUsed, *stageUsed; 1290 PetscBool *localStageVisible, *stageVisible; 1291 int numStages, localNumEvents, numEvents; 1292 int stage, lastStage, oclass; 1293 PetscLogEvent event; 1294 PetscErrorCode ierr; 1295 char version[256]; 1296 MPI_Comm comm; 1297 PetscInt nthreads; 1298 1299 PetscFunctionBegin; 1300 ierr = PetscObjectGetComm((PetscObject)viewer,&comm);CHKERRQ(ierr); 1301 if (!PetscLogBegin_PrivateCalled) SETERRQ(comm, PETSC_ERR_ORDER, "No call to PetscLogBegin() before PetscLogView()"); 1302 ierr = PetscViewerASCIIGetPointer(viewer,&fd);CHKERRQ(ierr); 1303 ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 1304 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 1305 /* Pop off any stages the user forgot to remove */ 1306 lastStage = 0; 1307 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1308 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 1309 while (stage >= 0) { 1310 lastStage = stage; 1311 ierr = PetscStageLogPop(stageLog);CHKERRQ(ierr); 1312 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 1313 } 1314 /* Get the total elapsed time */ 1315 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 1316 1317 ierr = PetscFPrintf(comm, fd, "************************************************************************************************************************\n");CHKERRQ(ierr); 1318 ierr = PetscFPrintf(comm, fd, "*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***\n");CHKERRQ(ierr); 1319 ierr = PetscFPrintf(comm, fd, "************************************************************************************************************************\n");CHKERRQ(ierr); 1320 ierr = PetscFPrintf(comm, fd, "\n---------------------------------------------- PETSc Performance Summary: ----------------------------------------------\n\n");CHKERRQ(ierr); 1321 ierr = PetscGetArchType(arch,sizeof(arch));CHKERRQ(ierr); 1322 ierr = PetscGetHostName(hostname,sizeof(hostname));CHKERRQ(ierr); 1323 ierr = PetscGetUserName(username,sizeof(username));CHKERRQ(ierr); 1324 ierr = PetscGetProgramName(pname,sizeof(pname));CHKERRQ(ierr); 1325 ierr = PetscGetDate(date,sizeof(date));CHKERRQ(ierr); 1326 ierr = PetscGetVersion(version,sizeof(version));CHKERRQ(ierr); 1327 if (size == 1) { 1328 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); 1329 } else { 1330 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); 1331 } 1332 ierr = PetscThreadCommGetNThreads(PETSC_COMM_WORLD,&nthreads);CHKERRQ(ierr); 1333 if (nthreads > 1) { 1334 ierr = PetscFPrintf(comm,fd,"With %d threads per MPI_Comm\n", (int)nthreads);CHKERRQ(ierr); 1335 } 1336 1337 ierr = PetscFPrintf(comm, fd, "Using %s\n", version);CHKERRQ(ierr); 1338 1339 /* Must preserve reduction count before we go on */ 1340 red = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1341 1342 /* Calculate summary information */ 1343 ierr = PetscFPrintf(comm, fd, "\n Max Max/Min Avg Total \n");CHKERRQ(ierr); 1344 /* Time */ 1345 ierr = MPI_Allreduce(&locTotalTime, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1346 ierr = MPI_Allreduce(&locTotalTime, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1347 ierr = MPI_Allreduce(&locTotalTime, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1348 avg = (tot)/((PetscLogDouble) size); 1349 if (min != 0.0) ratio = max/min; 1350 else ratio = 0.0; 1351 ierr = PetscFPrintf(comm, fd, "Time (sec): %5.3e %10.5f %5.3e\n", max, ratio, avg);CHKERRQ(ierr); 1352 TotalTime = tot; 1353 /* Objects */ 1354 avg = (PetscLogDouble) petsc_numObjects; 1355 ierr = MPI_Allreduce(&avg, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1356 ierr = MPI_Allreduce(&avg, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1357 ierr = MPI_Allreduce(&avg, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1358 avg = (tot)/((PetscLogDouble) size); 1359 if (min != 0.0) ratio = max/min; 1360 else ratio = 0.0; 1361 ierr = PetscFPrintf(comm, fd, "Objects: %5.3e %10.5f %5.3e\n", max, ratio, avg);CHKERRQ(ierr); 1362 /* Flops */ 1363 ierr = MPI_Allreduce(&petsc_TotalFlops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1364 ierr = MPI_Allreduce(&petsc_TotalFlops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1365 ierr = MPI_Allreduce(&petsc_TotalFlops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1366 avg = (tot)/((PetscLogDouble) size); 1367 if (min != 0.0) ratio = max/min; 1368 else ratio = 0.0; 1369 ierr = PetscFPrintf(comm, fd, "Flops: %5.3e %10.5f %5.3e %5.3e\n", max, ratio, avg, tot);CHKERRQ(ierr); 1370 TotalFlops = tot; 1371 /* Flops/sec -- Must talk to Barry here */ 1372 if (locTotalTime != 0.0) flops = petsc_TotalFlops/locTotalTime; 1373 else flops = 0.0; 1374 ierr = MPI_Allreduce(&flops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1375 ierr = MPI_Allreduce(&flops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1376 ierr = MPI_Allreduce(&flops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1377 avg = (tot)/((PetscLogDouble) size); 1378 if (min != 0.0) ratio = max/min; 1379 else ratio = 0.0; 1380 ierr = PetscFPrintf(comm, fd, "Flops/sec: %5.3e %10.5f %5.3e %5.3e\n", max, ratio, avg, tot);CHKERRQ(ierr); 1381 /* Memory */ 1382 ierr = PetscMallocGetMaximumUsage(&mem);CHKERRQ(ierr); 1383 if (mem > 0.0) { 1384 ierr = MPI_Allreduce(&mem, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1385 ierr = MPI_Allreduce(&mem, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1386 ierr = MPI_Allreduce(&mem, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1387 avg = (tot)/((PetscLogDouble) size); 1388 if (min != 0.0) ratio = max/min; 1389 else ratio = 0.0; 1390 ierr = PetscFPrintf(comm, fd, "Memory: %5.3e %10.5f %5.3e\n", max, ratio, tot);CHKERRQ(ierr); 1391 } 1392 /* Messages */ 1393 mess = 0.5*(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct); 1394 ierr = MPI_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1395 ierr = MPI_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1396 ierr = MPI_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1397 avg = (tot)/((PetscLogDouble) size); 1398 if (min != 0.0) ratio = max/min; 1399 else ratio = 0.0; 1400 ierr = PetscFPrintf(comm, fd, "MPI Messages: %5.3e %10.5f %5.3e %5.3e\n", max, ratio, avg, tot);CHKERRQ(ierr); 1401 numMessages = tot; 1402 /* Message Lengths */ 1403 mess = 0.5*(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len); 1404 ierr = MPI_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1405 ierr = MPI_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1406 ierr = MPI_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1407 if (numMessages != 0) avg = (tot)/(numMessages); 1408 else avg = 0.0; 1409 if (min != 0.0) ratio = max/min; 1410 else ratio = 0.0; 1411 ierr = PetscFPrintf(comm, fd, "MPI Message Lengths: %5.3e %10.5f %5.3e %5.3e\n", max, ratio, avg, tot);CHKERRQ(ierr); 1412 messageLength = tot; 1413 /* Reductions */ 1414 ierr = MPI_Allreduce(&red, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1415 ierr = MPI_Allreduce(&red, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1416 ierr = MPI_Allreduce(&red, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1417 if (min != 0.0) ratio = max/min; 1418 else ratio = 0.0; 1419 ierr = PetscFPrintf(comm, fd, "MPI Reductions: %5.3e %10.5f\n", max, ratio);CHKERRQ(ierr); 1420 numReductions = red; /* wrong because uses count from process zero */ 1421 ierr = PetscFPrintf(comm, fd, "\nFlop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)\n");CHKERRQ(ierr); 1422 ierr = PetscFPrintf(comm, fd, " e.g., VecAXPY() for real vectors of length N --> 2N flops\n");CHKERRQ(ierr); 1423 ierr = PetscFPrintf(comm, fd, " and VecAXPY() for complex vectors of length N --> 8N flops\n");CHKERRQ(ierr); 1424 1425 /* Get total number of stages -- 1426 Currently, a single processor can register more stages than another, but stages must all be registered in order. 1427 We can removed this requirement if necessary by having a global stage numbering and indirection on the stage ID. 1428 This seems best accomplished by assoicating a communicator with each stage. 1429 */ 1430 ierr = MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 1431 ierr = PetscMalloc(numStages * sizeof(PetscBool), &localStageUsed);CHKERRQ(ierr); 1432 ierr = PetscMalloc(numStages * sizeof(PetscBool), &stageUsed);CHKERRQ(ierr); 1433 ierr = PetscMalloc(numStages * sizeof(PetscBool), &localStageVisible);CHKERRQ(ierr); 1434 ierr = PetscMalloc(numStages * sizeof(PetscBool), &stageVisible);CHKERRQ(ierr); 1435 if (numStages > 0) { 1436 stageInfo = stageLog->stageInfo; 1437 for (stage = 0; stage < numStages; stage++) { 1438 if (stage < stageLog->numStages) { 1439 localStageUsed[stage] = stageInfo[stage].used; 1440 localStageVisible[stage] = stageInfo[stage].perfInfo.visible; 1441 } else { 1442 localStageUsed[stage] = PETSC_FALSE; 1443 localStageVisible[stage] = PETSC_TRUE; 1444 } 1445 } 1446 ierr = MPI_Allreduce(localStageUsed, stageUsed, numStages, MPIU_BOOL, MPI_LOR, comm);CHKERRQ(ierr); 1447 ierr = MPI_Allreduce(localStageVisible, stageVisible, numStages, MPIU_BOOL, MPI_LAND, comm);CHKERRQ(ierr); 1448 for (stage = 0; stage < numStages; stage++) { 1449 if (stageUsed[stage]) { 1450 ierr = PetscFPrintf(comm, fd, "\nSummary of Stages: ----- Time ------ ----- Flops ----- --- Messages --- -- Message Lengths -- -- Reductions --\n");CHKERRQ(ierr); 1451 ierr = PetscFPrintf(comm, fd, " Avg %%Total Avg %%Total counts %%Total Avg %%Total counts %%Total \n");CHKERRQ(ierr); 1452 break; 1453 } 1454 } 1455 for (stage = 0; stage < numStages; stage++) { 1456 if (!stageUsed[stage]) continue; 1457 if (localStageUsed[stage]) { 1458 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1459 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1460 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1461 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1462 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1463 name = stageInfo[stage].name; 1464 } else { 1465 ierr = MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1466 ierr = MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1467 ierr = MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1468 ierr = MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1469 ierr = MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1470 name = ""; 1471 } 1472 mess *= 0.5; messLen *= 0.5; red /= size; 1473 if (TotalTime != 0.0) fracTime = stageTime/TotalTime; else fracTime = 0.0; 1474 if (TotalFlops != 0.0) fracFlops = flops/TotalFlops; else fracFlops = 0.0; 1475 /* Talk to Barry if (stageTime != 0.0) flops = (size*flops)/stageTime; else flops = 0.0; */ 1476 if (numMessages != 0.0) fracMessages = mess/numMessages; else fracMessages = 0.0; 1477 if (numMessages != 0.0) avgMessLen = messLen/numMessages; else avgMessLen = 0.0; 1478 if (messageLength != 0.0) fracLength = messLen/messageLength; else fracLength = 0.0; 1479 if (numReductions != 0.0) fracReductions = red/numReductions; else fracReductions = 0.0; 1480 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", 1481 stage, name, stageTime/size, 100.0*fracTime, flops, 100.0*fracFlops, 1482 mess, 100.0*fracMessages, avgMessLen, 100.0*fracLength, red, 100.0*fracReductions);CHKERRQ(ierr); 1483 } 1484 } 1485 1486 ierr = PetscFPrintf(comm, fd,"\n------------------------------------------------------------------------------------------------------------------------\n");CHKERRQ(ierr); 1487 ierr = PetscFPrintf(comm, fd, "See the 'Profiling' chapter of the users' manual for details on interpreting output.\n");CHKERRQ(ierr); 1488 ierr = PetscFPrintf(comm, fd, "Phase summary info:\n");CHKERRQ(ierr); 1489 ierr = PetscFPrintf(comm, fd, " Count: number of times phase was executed\n");CHKERRQ(ierr); 1490 ierr = PetscFPrintf(comm, fd, " Time and Flops: Max - maximum over all processors\n");CHKERRQ(ierr); 1491 ierr = PetscFPrintf(comm, fd, " Ratio - ratio of maximum to minimum over all processors\n");CHKERRQ(ierr); 1492 ierr = PetscFPrintf(comm, fd, " Mess: number of messages sent\n");CHKERRQ(ierr); 1493 ierr = PetscFPrintf(comm, fd, " Avg. len: average message length (bytes)\n");CHKERRQ(ierr); 1494 ierr = PetscFPrintf(comm, fd, " Reduct: number of global reductions\n");CHKERRQ(ierr); 1495 ierr = PetscFPrintf(comm, fd, " Global: entire computation\n");CHKERRQ(ierr); 1496 ierr = PetscFPrintf(comm, fd, " Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().\n");CHKERRQ(ierr); 1497 ierr = PetscFPrintf(comm, fd, " %%T - percent time in this phase %%F - percent flops in this phase\n");CHKERRQ(ierr); 1498 ierr = PetscFPrintf(comm, fd, " %%M - percent messages in this phase %%L - percent message lengths in this phase\n");CHKERRQ(ierr); 1499 ierr = PetscFPrintf(comm, fd, " %%R - percent reductions in this phase\n");CHKERRQ(ierr); 1500 ierr = PetscFPrintf(comm, fd, " Total Mflop/s: 10e-6 * (sum of flops over all processors)/(max time over all processors)\n");CHKERRQ(ierr); 1501 ierr = PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------\n");CHKERRQ(ierr); 1502 1503 #if defined(PETSC_USE_DEBUG) 1504 ierr = PetscFPrintf(comm, fd, "\n\n");CHKERRQ(ierr); 1505 ierr = PetscFPrintf(comm, fd, " ##########################################################\n");CHKERRQ(ierr); 1506 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1507 ierr = PetscFPrintf(comm, fd, " # WARNING!!! #\n");CHKERRQ(ierr); 1508 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1509 ierr = PetscFPrintf(comm, fd, " # This code was compiled with a debugging option, #\n");CHKERRQ(ierr); 1510 ierr = PetscFPrintf(comm, fd, " # To get timing results run ./configure #\n");CHKERRQ(ierr); 1511 ierr = PetscFPrintf(comm, fd, " # using --with-debugging=no, the performance will #\n");CHKERRQ(ierr); 1512 ierr = PetscFPrintf(comm, fd, " # be generally two or three times faster. #\n");CHKERRQ(ierr); 1513 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1514 ierr = PetscFPrintf(comm, fd, " ##########################################################\n\n\n");CHKERRQ(ierr); 1515 #endif 1516 #if defined(PETSC_USE_COMPLEX) && !defined(PETSC_USE_FORTRAN_KERNELS) 1517 ierr = PetscFPrintf(comm, fd, "\n\n");CHKERRQ(ierr); 1518 ierr = PetscFPrintf(comm, fd, " ##########################################################\n");CHKERRQ(ierr); 1519 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1520 ierr = PetscFPrintf(comm, fd, " # WARNING!!! #\n");CHKERRQ(ierr); 1521 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1522 ierr = PetscFPrintf(comm, fd, " # The code for various complex numbers numerical #\n");CHKERRQ(ierr); 1523 ierr = PetscFPrintf(comm, fd, " # kernels uses C++, which generally is not well #\n");CHKERRQ(ierr); 1524 ierr = PetscFPrintf(comm, fd, " # optimized. For performance that is about 4-5 times #\n");CHKERRQ(ierr); 1525 ierr = PetscFPrintf(comm, fd, " # faster, specify --with-fortran-kernels=1 #\n");CHKERRQ(ierr); 1526 ierr = PetscFPrintf(comm, fd, " # when running ./configure.py. #\n");CHKERRQ(ierr); 1527 ierr = PetscFPrintf(comm, fd, " # #\n");CHKERRQ(ierr); 1528 ierr = PetscFPrintf(comm, fd, " ##########################################################\n\n\n");CHKERRQ(ierr); 1529 #endif 1530 1531 /* Report events */ 1532 ierr = PetscFPrintf(comm, fd,"Event Count Time (sec) Flops --- Global --- --- Stage --- Total\n");CHKERRQ(ierr); 1533 ierr = PetscFPrintf(comm, fd," Max Ratio Max Ratio Max Ratio Mess Avg len Reduct %%T %%F %%M %%L %%R %%T %%F %%M %%L %%R Mflop/s\n");CHKERRQ(ierr); 1534 ierr = PetscFPrintf(comm,fd,"------------------------------------------------------------------------------------------------------------------------\n");CHKERRQ(ierr); 1535 1536 /* Problem: The stage name will not show up unless the stage executed on proc 1 */ 1537 for (stage = 0; stage < numStages; stage++) { 1538 if (!stageVisible[stage]) continue; 1539 if (localStageUsed[stage]) { 1540 ierr = PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name);CHKERRQ(ierr); 1541 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1542 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1543 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1544 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1545 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1546 } else { 1547 ierr = PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage);CHKERRQ(ierr); 1548 ierr = MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1549 ierr = MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1550 ierr = MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1551 ierr = MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1552 ierr = MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1553 } 1554 mess *= 0.5; messLen *= 0.5; red /= size; 1555 1556 /* Get total number of events in this stage -- 1557 Currently, a single processor can register more events than another, but events must all be registered in order, 1558 just like stages. We can removed this requirement if necessary by having a global event numbering and indirection 1559 on the event ID. This seems best accomplished by assoicating a communicator with each stage. 1560 1561 Problem: If the event did not happen on proc 1, its name will not be available. 1562 Problem: Event visibility is not implemented 1563 */ 1564 if (localStageUsed[stage]) { 1565 eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo; 1566 localNumEvents = stageLog->stageInfo[stage].eventLog->numEvents; 1567 } else localNumEvents = 0; 1568 ierr = MPI_Allreduce(&localNumEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 1569 for (event = 0; event < numEvents; event++) { 1570 if (localStageUsed[stage] && (event < stageLog->stageInfo[stage].eventLog->numEvents) && (eventInfo[event].depth == 0)) { 1571 if ((eventInfo[event].count > 0) && (eventInfo[event].time > 0.0)) flopr = eventInfo[event].flops; 1572 else flopr = 0.0; 1573 1574 ierr = MPI_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1575 ierr = MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1576 ierr = MPI_Allreduce(&eventInfo[event].flops, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1577 ierr = MPI_Allreduce(&eventInfo[event].time, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1578 ierr = MPI_Allreduce(&eventInfo[event].time, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1579 ierr = MPI_Allreduce(&eventInfo[event].time, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1580 ierr = MPI_Allreduce(&eventInfo[event].numMessages, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1581 ierr = MPI_Allreduce(&eventInfo[event].messageLength, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1582 ierr = MPI_Allreduce(&eventInfo[event].numReductions, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1583 ierr = MPI_Allreduce(&eventInfo[event].count, &minCt, 1, MPI_INT, MPI_MIN, comm);CHKERRQ(ierr); 1584 ierr = MPI_Allreduce(&eventInfo[event].count, &maxCt, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 1585 name = stageLog->eventLog->eventInfo[event].name; 1586 } else { 1587 flopr = 0.0; 1588 ierr = MPI_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1589 ierr = MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1590 ierr = MPI_Allreduce(&zero, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1591 ierr = MPI_Allreduce(&zero, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm);CHKERRQ(ierr); 1592 ierr = MPI_Allreduce(&zero, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm);CHKERRQ(ierr); 1593 ierr = MPI_Allreduce(&zero, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1594 ierr = MPI_Allreduce(&zero, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1595 ierr = MPI_Allreduce(&zero, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1596 ierr = MPI_Allreduce(&zero, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 1597 ierr = MPI_Allreduce(&ierr, &minCt, 1, MPI_INT, MPI_MIN, comm);CHKERRQ(ierr); 1598 ierr = MPI_Allreduce(&ierr, &maxCt, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 1599 name = ""; 1600 } 1601 if (mint < 0.0) { 1602 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); 1603 mint = 0; 1604 } 1605 if (minf < 0.0) SETERRQ2(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Minimum flops %g over all processors for %s is negative! Not possible!",minf,name); 1606 totm *= 0.5; totml *= 0.5; totr /= size; 1607 1608 if (maxCt != 0) { 1609 if (minCt != 0) ratCt = ((PetscLogDouble) maxCt)/minCt; else ratCt = 0.0; 1610 if (mint != 0.0) ratt = maxt/mint; else ratt = 0.0; 1611 if (minf != 0.0) ratf = maxf/minf; else ratf = 0.0; 1612 if (TotalTime != 0.0) fracTime = tott/TotalTime; else fracTime = 0.0; 1613 if (TotalFlops != 0.0) fracFlops = totf/TotalFlops; else fracFlops = 0.0; 1614 if (stageTime != 0.0) fracStageTime = tott/stageTime; else fracStageTime = 0.0; 1615 if (flops != 0.0) fracStageFlops = totf/flops; else fracStageFlops = 0.0; 1616 if (numMessages != 0.0) fracMess = totm/numMessages; else fracMess = 0.0; 1617 if (messageLength != 0.0) fracMessLen = totml/messageLength; else fracMessLen = 0.0; 1618 if (numReductions != 0.0) fracRed = totr/numReductions; else fracRed = 0.0; 1619 if (mess != 0.0) fracStageMess = totm/mess; else fracStageMess = 0.0; 1620 if (messLen != 0.0) fracStageMessLen = totml/messLen; else fracStageMessLen = 0.0; 1621 if (red != 0.0) fracStageRed = totr/red; else fracStageRed = 0.0; 1622 if (totm != 0.0) totml /= totm; else totml = 0.0; 1623 if (maxt != 0.0) flopr = totf/maxt; else flopr = 0.0; 1624 if (fracStageTime > 1.00) ierr = PetscFPrintf(comm, fd,"Warning -- total time of even greater than time of entire stage -- something is wrong with the timer\n");CHKERRQ(ierr); 1625 ierr = PetscFPrintf(comm, fd, 1626 "%-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\n", 1627 name, maxCt, ratCt, maxt, ratt, maxf, ratf, totm, totml, totr, 1628 100.0*fracTime, 100.0*fracFlops, 100.0*fracMess, 100.0*fracMessLen, 100.0*fracRed, 1629 100.0*fracStageTime, 100.0*fracStageFlops, 100.0*fracStageMess, 100.0*fracStageMessLen, 100.0*fracStageRed, 1630 PetscAbsReal(flopr/1.0e6));CHKERRQ(ierr); 1631 } 1632 } 1633 } 1634 1635 /* Memory usage and object creation */ 1636 ierr = PetscFPrintf(comm, fd, "------------------------------------------------------------------------------------------------------------------------\n");CHKERRQ(ierr); 1637 ierr = PetscFPrintf(comm, fd, "\n");CHKERRQ(ierr); 1638 ierr = PetscFPrintf(comm, fd, "Memory usage is given in bytes:\n\n");CHKERRQ(ierr); 1639 1640 /* Right now, only stages on the first processor are reported here, meaning only objects associated with 1641 the global communicator, or MPI_COMM_SELF for proc 1. We really should report global stats and then 1642 stats for stages local to processor sets. 1643 */ 1644 /* We should figure out the longest object name here (now 20 characters) */ 1645 ierr = PetscFPrintf(comm, fd, "Object Type Creations Destructions Memory Descendants' Mem.\n");CHKERRQ(ierr); 1646 ierr = PetscFPrintf(comm, fd, "Reports information only for process 0.\n");CHKERRQ(ierr); 1647 for (stage = 0; stage < numStages; stage++) { 1648 if (localStageUsed[stage]) { 1649 classInfo = stageLog->stageInfo[stage].classLog->classInfo; 1650 ierr = PetscFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name);CHKERRQ(ierr); 1651 for (oclass = 0; oclass < stageLog->stageInfo[stage].classLog->numClasses; oclass++) { 1652 if ((classInfo[oclass].creations > 0) || (classInfo[oclass].destructions > 0)) { 1653 ierr = PetscFPrintf(comm, fd, "%20s %5d %5d %11.0f %g\n", stageLog->classLog->classInfo[oclass].name, 1654 classInfo[oclass].creations, classInfo[oclass].destructions, classInfo[oclass].mem, 1655 classInfo[oclass].descMem);CHKERRQ(ierr); 1656 } 1657 } 1658 } else { 1659 ierr = PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage);CHKERRQ(ierr); 1660 } 1661 } 1662 1663 ierr = PetscFree(localStageUsed);CHKERRQ(ierr); 1664 ierr = PetscFree(stageUsed);CHKERRQ(ierr); 1665 ierr = PetscFree(localStageVisible);CHKERRQ(ierr); 1666 ierr = PetscFree(stageVisible);CHKERRQ(ierr); 1667 1668 /* Information unrelated to this particular run */ 1669 ierr = PetscFPrintf(comm, fd, "========================================================================================================================\n");CHKERRQ(ierr); 1670 PetscTime(&y); 1671 PetscTime(&x); 1672 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 1673 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 1674 ierr = PetscFPrintf(comm,fd,"Average time to get PetscTime(): %g\n", (y-x)/10.0);CHKERRQ(ierr); 1675 /* MPI information */ 1676 if (size > 1) { 1677 MPI_Status status; 1678 PetscMPIInt tag; 1679 MPI_Comm newcomm; 1680 1681 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1682 PetscTime(&x); 1683 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1684 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1685 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1686 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1687 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1688 PetscTime(&y); 1689 ierr = PetscFPrintf(comm, fd, "Average time for MPI_Barrier(): %g\n", (y-x)/5.0);CHKERRQ(ierr); 1690 ierr = PetscCommDuplicate(comm,&newcomm, &tag);CHKERRQ(ierr); 1691 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 1692 if (rank) { 1693 ierr = MPI_Recv(0, 0, MPI_INT, rank-1, tag, newcomm, &status);CHKERRQ(ierr); 1694 ierr = MPI_Send(0, 0, MPI_INT, (rank+1)%size, tag, newcomm);CHKERRQ(ierr); 1695 } else { 1696 PetscTime(&x); 1697 ierr = MPI_Send(0, 0, MPI_INT, 1, tag, newcomm);CHKERRQ(ierr); 1698 ierr = MPI_Recv(0, 0, MPI_INT, size-1, tag, newcomm, &status);CHKERRQ(ierr); 1699 PetscTime(&y); 1700 ierr = PetscFPrintf(comm,fd,"Average time for zero size MPI_Send(): %g\n", (y-x)/size);CHKERRQ(ierr); 1701 } 1702 ierr = PetscCommDestroy(&newcomm);CHKERRQ(ierr); 1703 } 1704 ierr = PetscOptionsView(viewer);CHKERRQ(ierr); 1705 1706 /* Machine and compile information */ 1707 #if defined(PETSC_USE_FORTRAN_KERNELS) 1708 ierr = PetscFPrintf(comm, fd, "Compiled with FORTRAN kernels\n");CHKERRQ(ierr); 1709 #else 1710 ierr = PetscFPrintf(comm, fd, "Compiled without FORTRAN kernels\n");CHKERRQ(ierr); 1711 #endif 1712 #if defined(PETSC_USE_REAL_SINGLE) 1713 ierr = PetscFPrintf(comm, fd, "Compiled with single precision PetscScalar and PetscReal\n");CHKERRQ(ierr); 1714 #elif defined(PETSC_USE_LONGDOUBLE) 1715 ierr = PetscFPrintf(comm, fd, "Compiled with long double precision PetscScalar and PetscReal\n");CHKERRQ(ierr); 1716 #endif 1717 1718 #if defined(PETSC_USE_REAL_MAT_SINGLE) 1719 ierr = PetscFPrintf(comm, fd, "Compiled with single precision matrices\n");CHKERRQ(ierr); 1720 #else 1721 ierr = PetscFPrintf(comm, fd, "Compiled with full precision matrices (default)\n");CHKERRQ(ierr); 1722 #endif 1723 ierr = PetscFPrintf(comm, fd, "sizeof(short) %d sizeof(int) %d sizeof(long) %d sizeof(void*) %d sizeof(PetscScalar) %d sizeof(PetscInt) %d\n", 1724 (int) sizeof(short), (int) sizeof(int), (int) sizeof(long), (int) sizeof(void*),(int) sizeof(PetscScalar),(int) sizeof(PetscInt));CHKERRQ(ierr); 1725 1726 ierr = PetscFPrintf(comm, fd, "Configure run at: %s\n",petscconfigureruntime);CHKERRQ(ierr); 1727 ierr = PetscFPrintf(comm, fd, "Configure options: %s",petscconfigureoptions);CHKERRQ(ierr); 1728 ierr = PetscFPrintf(comm, fd, "%s", petscmachineinfo);CHKERRQ(ierr); 1729 ierr = PetscFPrintf(comm, fd, "%s", petsccompilerinfo);CHKERRQ(ierr); 1730 ierr = PetscFPrintf(comm, fd, "%s", petsccompilerflagsinfo);CHKERRQ(ierr); 1731 ierr = PetscFPrintf(comm, fd, "%s", petsclinkerinfo);CHKERRQ(ierr); 1732 1733 /* Cleanup */ 1734 ierr = PetscFPrintf(comm, fd, "\n");CHKERRQ(ierr); 1735 ierr = PetscStageLogPush(stageLog, lastStage);CHKERRQ(ierr); 1736 PetscFunctionReturn(0); 1737 } 1738 1739 #undef __FUNCT__ 1740 #define __FUNCT__ "PetscLogPrintDetailed" 1741 /*@C 1742 PetscLogPrintDetailed - Each process prints the times for its own events 1743 1744 Collective over MPI_Comm 1745 1746 Input Parameter: 1747 + comm - The MPI communicator (only one processor prints output) 1748 - file - [Optional] The output file name 1749 1750 Options Database Keys: 1751 . -log_summary_detailed - Prints summary of log information (for code compiled with PETSC_USE_LOG) 1752 1753 Usage: 1754 .vb 1755 PetscInitialize(...); 1756 PetscLogBegin(); 1757 ... code ... 1758 PetscLogPrintDetailed(MPI_Comm,filename); 1759 PetscFinalize(...); 1760 .ve 1761 1762 Notes: 1763 By default the summary is printed to stdout. 1764 1765 Level: beginner 1766 1767 .keywords: log, dump, print 1768 .seealso: PetscLogBegin(), PetscLogDump(), PetscLogView() 1769 @*/ 1770 PetscErrorCode PetscLogPrintDetailed(MPI_Comm comm, const char filename[]) 1771 { 1772 FILE *fd = PETSC_STDOUT; 1773 PetscStageLog stageLog; 1774 PetscStageInfo *stageInfo = NULL; 1775 PetscEventPerfInfo *eventInfo = NULL; 1776 const char *name = NULL; 1777 PetscLogDouble TotalTime; 1778 PetscLogDouble stageTime, flops, flopr, mess, messLen, red; 1779 PetscLogDouble maxf, totf, maxt, tott, totm, totml, totr = 0.0; 1780 PetscMPIInt maxCt; 1781 PetscBool *stageUsed; 1782 PetscBool *stageVisible; 1783 int numStages, numEvents; 1784 int stage; 1785 PetscLogEvent event; 1786 PetscErrorCode ierr; 1787 1788 PetscFunctionBegin; 1789 if (!PetscLogBegin_PrivateCalled) SETERRQ(comm, PETSC_ERR_ORDER, "No call to PetscLogBegin() before PetscLogPrintDetailed()"); 1790 /* Pop off any stages the user forgot to remove */ 1791 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 1792 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 1793 while (stage >= 0) { 1794 ierr = PetscStageLogPop(stageLog);CHKERRQ(ierr); 1795 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 1796 } 1797 /* Get the total elapsed time */ 1798 PetscTime(&TotalTime); TotalTime -= petsc_BaseTime; 1799 /* Open the summary file */ 1800 if (filename) { 1801 ierr = PetscFOpen(comm, filename, "w", &fd);CHKERRQ(ierr); 1802 } 1803 1804 ierr = PetscFPrintf(comm, fd, "************************************************************************************************************************\n");CHKERRQ(ierr); 1805 ierr = PetscFPrintf(comm, fd, "*** WIDEN YOUR WINDOW TO 120 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***\n");CHKERRQ(ierr); 1806 ierr = PetscFPrintf(comm, fd, "************************************************************************************************************************\n");CHKERRQ(ierr); 1807 1808 1809 numStages = stageLog->numStages; 1810 ierr = PetscMalloc(numStages * sizeof(PetscBool), &stageUsed);CHKERRQ(ierr); 1811 ierr = PetscMalloc(numStages * sizeof(PetscBool), &stageVisible);CHKERRQ(ierr); 1812 if (numStages > 0) { 1813 stageInfo = stageLog->stageInfo; 1814 for (stage = 0; stage < numStages; stage++) { 1815 if (stage < stageLog->numStages) { 1816 stageUsed[stage] = stageInfo[stage].used; 1817 stageVisible[stage] = stageInfo[stage].perfInfo.visible; 1818 } else { 1819 stageUsed[stage] = PETSC_FALSE; 1820 stageVisible[stage] = PETSC_TRUE; 1821 } 1822 } 1823 } 1824 1825 /* Report events */ 1826 ierr = PetscFPrintf(comm, fd,"Event Count Time (sec) Flops/sec \n");CHKERRQ(ierr); 1827 ierr = PetscFPrintf(comm, fd," Mess Avg len Reduct \n");CHKERRQ(ierr); 1828 ierr = PetscFPrintf(comm,fd,"-----------------------------------------------------------------------------------\n");CHKERRQ(ierr); 1829 /* Problem: The stage name will not show up unless the stage executed on proc 1 */ 1830 for (stage = 0; stage < numStages; stage++) { 1831 if (!stageVisible[stage]) continue; 1832 if (stageUsed[stage]) { 1833 ierr = PetscSynchronizedFPrintf(comm, fd, "\n--- Event Stage %d: %s\n\n", stage, stageInfo[stage].name);CHKERRQ(ierr); 1834 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1835 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1836 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1837 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1838 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1839 } 1840 mess *= 0.5; messLen *= 0.5; 1841 1842 /* Get total number of events in this stage -- 1843 */ 1844 if (stageUsed[stage]) { 1845 eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo; 1846 numEvents = stageLog->stageInfo[stage].eventLog->numEvents; 1847 } else numEvents = 0; 1848 1849 for (event = 0; event < numEvents; event++) { 1850 if (stageUsed[stage] && (event < stageLog->stageInfo[stage].eventLog->numEvents)) { 1851 if ((eventInfo[event].count > 0) && (eventInfo[event].time > 0.0)) flopr = eventInfo[event].flops/eventInfo[event].time; 1852 else flopr = 0.0; 1853 1854 ierr = MPI_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, PETSC_COMM_SELF);CHKERRQ(ierr); 1855 ierr = MPI_Allreduce(&eventInfo[event].flops, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1856 ierr = MPI_Allreduce(&eventInfo[event].time, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, PETSC_COMM_SELF);CHKERRQ(ierr); 1857 ierr = MPI_Allreduce(&eventInfo[event].time, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1858 ierr = MPI_Allreduce(&eventInfo[event].numMessages, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1859 ierr = MPI_Allreduce(&eventInfo[event].messageLength, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, PETSC_COMM_SELF);CHKERRQ(ierr); 1860 totr = eventInfo[event].numReductions; 1861 ierr = MPI_Allreduce(&eventInfo[event].count, &maxCt, 1, MPI_INT, MPI_MAX, PETSC_COMM_SELF);CHKERRQ(ierr); 1862 name = stageLog->eventLog->eventInfo[event].name; 1863 totm *= 0.5; totml *= 0.5; 1864 } 1865 1866 if (maxCt != 0) { 1867 if (totm != 0.0) totml /= totm; 1868 else totml = 0.0; 1869 ierr = PetscSynchronizedFPrintf(comm, fd,"%-16s %7d %5.4e %3.2e %2.1e %2.1e %2.1e\n",name, maxCt, maxt, maxf, totm, totml, totr);CHKERRQ(ierr); 1870 } 1871 } 1872 } 1873 ierr = PetscSynchronizedFlush(comm);CHKERRQ(ierr); 1874 1875 ierr = PetscFree(stageUsed);CHKERRQ(ierr); 1876 ierr = PetscFree(stageVisible);CHKERRQ(ierr); 1877 1878 ierr = PetscFClose(comm, fd);CHKERRQ(ierr); 1879 PetscFunctionReturn(0); 1880 } 1881 1882 /*----------------------------------------------- Counter Functions -------------------------------------------------*/ 1883 #undef __FUNCT__ 1884 #define __FUNCT__ "PetscGetFlops" 1885 /*@C 1886 PetscGetFlops - Returns the number of flops used on this processor 1887 since the program began. 1888 1889 Not Collective 1890 1891 Output Parameter: 1892 flops - number of floating point operations 1893 1894 Notes: 1895 A global counter logs all PETSc flop counts. The user can use 1896 PetscLogFlops() to increment this counter to include flops for the 1897 application code. 1898 1899 PETSc automatically logs library events if the code has been 1900 compiled with -DPETSC_USE_LOG (which is the default), and -log, 1901 -log_summary, or -log_all are specified. PetscLogFlops() is 1902 intended for logging user flops to supplement this PETSc 1903 information. 1904 1905 Level: intermediate 1906 1907 .keywords: log, flops, floating point operations 1908 1909 .seealso: PetscTime(), PetscLogFlops() 1910 @*/ 1911 PetscErrorCode PetscGetFlops(PetscLogDouble *flops) 1912 { 1913 PetscFunctionBegin; 1914 *flops = petsc_TotalFlops; 1915 PetscFunctionReturn(0); 1916 } 1917 1918 #undef __FUNCT__ 1919 #define __FUNCT__ "PetscLogObjectState" 1920 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...) 1921 { 1922 PetscErrorCode ierr; 1923 size_t fullLength; 1924 va_list Argp; 1925 1926 PetscFunctionBegin; 1927 if (!petsc_logObjects) PetscFunctionReturn(0); 1928 va_start(Argp, format); 1929 ierr = PetscVSNPrintf(petsc_objects[obj->id].info, 64,format,&fullLength, Argp);CHKERRQ(ierr); 1930 va_end(Argp); 1931 PetscFunctionReturn(0); 1932 } 1933 1934 1935 /*MC 1936 PetscLogFlops - Adds floating point operations to the global counter. 1937 1938 Synopsis: 1939 #include "petsclog.h" 1940 PetscErrorCode PetscLogFlops(PetscLogDouble f) 1941 1942 Not Collective 1943 1944 Input Parameter: 1945 . f - flop counter 1946 1947 1948 Usage: 1949 .vb 1950 PetscLogEvent USER_EVENT; 1951 PetscLogEventRegister("User event",0,&USER_EVENT); 1952 PetscLogEventBegin(USER_EVENT,0,0,0,0); 1953 [code segment to monitor] 1954 PetscLogFlops(user_flops) 1955 PetscLogEventEnd(USER_EVENT,0,0,0,0); 1956 .ve 1957 1958 Notes: 1959 A global counter logs all PETSc flop counts. The user can use 1960 PetscLogFlops() to increment this counter to include flops for the 1961 application code. 1962 1963 PETSc automatically logs library events if the code has been 1964 compiled with -DPETSC_USE_LOG (which is the default), and -log, 1965 -log_summary, or -log_all are specified. PetscLogFlops() is 1966 intended for logging user flops to supplement this PETSc 1967 information. 1968 1969 Level: intermediate 1970 1971 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscGetFlops() 1972 1973 .keywords: log, flops, floating point operations 1974 M*/ 1975 1976 /*MC 1977 PetscPreLoadBegin - Begin a segment of code that may be preloaded (run twice) 1978 to get accurate timings 1979 1980 Synopsis: 1981 #include "petsclog.h" 1982 void PetscPreLoadBegin(PetscBool flag,char *name); 1983 1984 Not Collective 1985 1986 Input Parameter: 1987 + flag - PETSC_TRUE to run twice, PETSC_FALSE to run once, may be overridden 1988 with command line option -preload true or -preload false 1989 - name - name of first stage (lines of code timed separately with -log_summary) to 1990 be preloaded 1991 1992 Usage: 1993 .vb 1994 PetscPreLoadBegin(PETSC_TRUE,"first stage); 1995 lines of code 1996 PetscPreLoadStage("second stage"); 1997 lines of code 1998 PetscPreLoadEnd(); 1999 .ve 2000 2001 Notes: Only works in C/C++, not Fortran 2002 2003 Flags available within the macro. 2004 + PetscPreLoadingUsed - true if we are or have done preloading 2005 . PetscPreLoadingOn - true if it is CURRENTLY doing preload 2006 . PetscPreLoadIt - 0 for the first computation (with preloading turned off it is only 0) 1 for the second 2007 - PetscPreLoadMax - number of times it will do the computation, only one when preloading is turned on 2008 The first two variables are available throughout the program, the second two only between the PetscPreLoadBegin() 2009 and PetscPreLoadEnd() 2010 2011 Level: intermediate 2012 2013 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadEnd(), PetscPreLoadStage() 2014 2015 Concepts: preloading 2016 Concepts: timing^accurate 2017 Concepts: paging^eliminating effects of 2018 2019 2020 M*/ 2021 2022 /*MC 2023 PetscPreLoadEnd - End a segment of code that may be preloaded (run twice) 2024 to get accurate timings 2025 2026 Synopsis: 2027 #include "petsclog.h" 2028 void PetscPreLoadEnd(void); 2029 2030 Not Collective 2031 2032 Usage: 2033 .vb 2034 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2035 lines of code 2036 PetscPreLoadStage("second stage"); 2037 lines of code 2038 PetscPreLoadEnd(); 2039 .ve 2040 2041 Notes: only works in C/C++ not fortran 2042 2043 Level: intermediate 2044 2045 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadStage() 2046 2047 M*/ 2048 2049 /*MC 2050 PetscPreLoadStage - Start a new segment of code to be timed separately. 2051 to get accurate timings 2052 2053 Synopsis: 2054 #include "petsclog.h" 2055 void PetscPreLoadStage(char *name); 2056 2057 Not Collective 2058 2059 Usage: 2060 .vb 2061 PetscPreLoadBegin(PETSC_TRUE,"first stage); 2062 lines of code 2063 PetscPreLoadStage("second stage"); 2064 lines of code 2065 PetscPreLoadEnd(); 2066 .ve 2067 2068 Notes: only works in C/C++ not fortran 2069 2070 Level: intermediate 2071 2072 .seealso: PetscLogEventRegister(), PetscLogEventBegin(), PetscLogEventEnd(), PetscPreLoadBegin(), PetscPreLoadEnd() 2073 2074 M*/ 2075 2076 #undef __FUNCT__ 2077 #define __FUNCT__ "PetscLogViewPython" 2078 /*@ 2079 PetscLogViewPython - Saves logging information in a Python format. 2080 2081 Collective on PetscViewer 2082 2083 Input Paramter: 2084 . viewer - viewer to save Python data 2085 2086 Level: intermediate 2087 2088 @*/ 2089 PetscErrorCode PetscLogViewPython(PetscViewer viewer) 2090 { 2091 FILE *fd; 2092 PetscLogDouble zero = 0.0; 2093 PetscStageLog stageLog; 2094 PetscStageInfo *stageInfo = NULL; 2095 PetscEventPerfInfo *eventInfo = NULL; 2096 const char *name; 2097 char stageName[2048]; 2098 char eventName[2048]; 2099 PetscLogDouble locTotalTime, TotalTime = 0, TotalFlops = 0; 2100 PetscLogDouble numMessages = 0, messageLength = 0, avgMessLen, numReductions = 0; 2101 PetscLogDouble stageTime, flops, mem, mess, messLen, red; 2102 PetscLogDouble fracTime, fracFlops, fracMessages, fracLength; 2103 PetscLogDouble fracReductions; 2104 PetscLogDouble tot,avg,x,y,*mydata; 2105 PetscMPIInt maxCt; 2106 PetscMPIInt size, rank, *mycount; 2107 PetscBool *localStageUsed, *stageUsed; 2108 PetscBool *localStageVisible, *stageVisible; 2109 int numStages, localNumEvents, numEvents; 2110 int stage, lastStage; 2111 PetscLogEvent event; 2112 PetscErrorCode ierr; 2113 PetscInt i; 2114 MPI_Comm comm; 2115 2116 PetscFunctionBegin; 2117 ierr = PetscObjectGetComm((PetscObject)viewer,&comm);CHKERRQ(ierr); 2118 if (!PetscLogBegin_PrivateCalled) SETERRQ(comm, PETSC_ERR_ORDER, "No call to PetscLogBegin() before PetscLogViewPython()"); 2119 ierr = MPI_Comm_size(comm, &size);CHKERRQ(ierr); 2120 ierr = MPI_Comm_rank(comm, &rank);CHKERRQ(ierr); 2121 ierr = PetscMalloc(size*sizeof(PetscLogDouble), &mydata);CHKERRQ(ierr); 2122 ierr = PetscMalloc(size*sizeof(PetscMPIInt), &mycount);CHKERRQ(ierr); 2123 ierr = PetscViewerASCIIGetPointer(viewer,&fd);CHKERRQ(ierr); 2124 2125 /* Pop off any stages the user forgot to remove */ 2126 lastStage = 0; 2127 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 2128 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 2129 while (stage >= 0) { 2130 lastStage = stage; 2131 ierr = PetscStageLogPop(stageLog);CHKERRQ(ierr); 2132 ierr = PetscStageLogGetCurrent(stageLog, &stage);CHKERRQ(ierr); 2133 } 2134 /* Get the total elapsed time */ 2135 PetscTime(&locTotalTime); locTotalTime -= petsc_BaseTime; 2136 2137 ierr = PetscFPrintf(comm, fd, "\n#------ PETSc Performance Summary ----------\n\n");CHKERRQ(ierr); 2138 ierr = PetscFPrintf(comm, fd, "Nproc = %d\n",size);CHKERRQ(ierr); 2139 2140 /* Must preserve reduction count before we go on */ 2141 red = (petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct)/((PetscLogDouble) size); 2142 2143 /* Calculate summary information */ 2144 2145 /* Time */ 2146 ierr = MPI_Gather(&locTotalTime,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2147 if (!rank) { 2148 ierr = PetscFPrintf(comm, fd, "Time = [ ");CHKERRQ(ierr); 2149 tot = 0.0; 2150 for (i=0; i<size; i++) { 2151 tot += mydata[i]; 2152 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2153 } 2154 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2155 avg = (tot)/((PetscLogDouble) size); 2156 TotalTime = tot; 2157 } 2158 2159 /* Objects */ 2160 avg = (PetscLogDouble) petsc_numObjects; 2161 ierr = MPI_Gather(&avg,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2162 if (!rank) { 2163 ierr = PetscFPrintf(comm, fd, "Objects = [ ");CHKERRQ(ierr); 2164 for (i=0; i<size; i++) { 2165 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2166 } 2167 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2168 } 2169 2170 /* Flops */ 2171 ierr = MPI_Gather(&petsc_TotalFlops,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2172 if (!rank) { 2173 ierr = PetscFPrintf(comm, fd, "Flops = [ ");CHKERRQ(ierr); 2174 tot = 0.0; 2175 for (i=0; i<size; i++) { 2176 tot += mydata[i]; 2177 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2178 } 2179 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2180 TotalFlops = tot; 2181 } 2182 2183 /* Memory */ 2184 ierr = PetscMallocGetMaximumUsage(&mem);CHKERRQ(ierr); 2185 ierr = MPI_Gather(&mem,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2186 if (!rank) { 2187 ierr = PetscFPrintf(comm, fd, "Memory = [ ");CHKERRQ(ierr); 2188 for (i=0; i<size; i++) { 2189 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2190 } 2191 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2192 } 2193 2194 /* Messages */ 2195 mess = 0.5*(petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct); 2196 ierr = MPI_Gather(&mess,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2197 if (!rank) { 2198 ierr = PetscFPrintf(comm, fd, "MPIMessages = [ ");CHKERRQ(ierr); 2199 tot = 0.0; 2200 for (i=0; i<size; i++) { 2201 tot += mydata[i]; 2202 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2203 } 2204 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2205 numMessages = tot; 2206 } 2207 2208 /* Message Lengths */ 2209 mess = 0.5*(petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len); 2210 ierr = MPI_Gather(&mess,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2211 if (!rank) { 2212 ierr = PetscFPrintf(comm, fd, "MPIMessageLengths = [ ");CHKERRQ(ierr); 2213 tot = 0.0; 2214 for (i=0; i<size; i++) { 2215 tot += mydata[i]; 2216 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2217 } 2218 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2219 messageLength = tot; 2220 } 2221 2222 /* Reductions */ 2223 ierr = MPI_Gather(&red,1,MPIU_PETSCLOGDOUBLE,mydata,1,MPIU_PETSCLOGDOUBLE,0,comm);CHKERRQ(ierr); 2224 if (!rank) { 2225 ierr = PetscFPrintf(comm, fd, "MPIReductions = [ ");CHKERRQ(ierr); 2226 tot = 0.0; 2227 for (i=0; i<size; i++) { 2228 tot += mydata[i]; 2229 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2230 } 2231 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2232 numReductions = tot; 2233 } 2234 2235 /* Get total number of stages -- 2236 Currently, a single processor can register more stages than another, but stages must all be registered in order. 2237 We can removed this requirement if necessary by having a global stage numbering and indirection on the stage ID. 2238 This seems best accomplished by assoicating a communicator with each stage. 2239 */ 2240 ierr = MPI_Allreduce(&stageLog->numStages, &numStages, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 2241 ierr = PetscMalloc(numStages * sizeof(PetscBool), &localStageUsed);CHKERRQ(ierr); 2242 ierr = PetscMalloc(numStages * sizeof(PetscBool), &stageUsed);CHKERRQ(ierr); 2243 ierr = PetscMalloc(numStages * sizeof(PetscBool), &localStageVisible);CHKERRQ(ierr); 2244 ierr = PetscMalloc(numStages * sizeof(PetscBool), &stageVisible);CHKERRQ(ierr); 2245 if (numStages > 0) { 2246 stageInfo = stageLog->stageInfo; 2247 for (stage = 0; stage < numStages; stage++) { 2248 if (stage < stageLog->numStages) { 2249 localStageUsed[stage] = stageInfo[stage].used; 2250 localStageVisible[stage] = stageInfo[stage].perfInfo.visible; 2251 } else { 2252 localStageUsed[stage] = PETSC_FALSE; 2253 localStageVisible[stage] = PETSC_TRUE; 2254 } 2255 } 2256 ierr = MPI_Allreduce(localStageUsed, stageUsed, numStages, MPIU_BOOL, MPI_LOR, comm);CHKERRQ(ierr); 2257 ierr = MPI_Allreduce(localStageVisible, stageVisible, numStages, MPIU_BOOL, MPI_LAND, comm);CHKERRQ(ierr); 2258 for (stage = 0; stage < numStages; stage++) { 2259 if (stageUsed[stage]) { 2260 ierr = PetscFPrintf(comm, fd, "\n#Summary of Stages: ----- Time ------ ----- Flops ----- --- Messages --- -- Message Lengths -- -- Reductions --\n");CHKERRQ(ierr); 2261 ierr = PetscFPrintf(comm, fd, "# Avg %%Total Avg %%Total counts %%Total Avg %%Total counts %%Total \n");CHKERRQ(ierr); 2262 break; 2263 } 2264 } 2265 for (stage = 0; stage < numStages; stage++) { 2266 if (!stageUsed[stage]) continue; 2267 if (localStageUsed[stage]) { 2268 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2269 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2270 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2271 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2272 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2273 name = stageInfo[stage].name; 2274 } else { 2275 ierr = MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2276 ierr = MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2277 ierr = MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2278 ierr = MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2279 ierr = MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2280 name = ""; 2281 } 2282 mess *= 0.5; messLen *= 0.5; red /= size; 2283 if (TotalTime != 0.0) fracTime = stageTime/TotalTime; else fracTime = 0.0; 2284 if (TotalFlops != 0.0) fracFlops = flops/TotalFlops; else fracFlops = 0.0; 2285 /* Talk to Barry if (stageTime != 0.0) flops = (size*flops)/stageTime; else flops = 0.0; */ 2286 if (numMessages != 0.0) fracMessages = mess/numMessages; else fracMessages = 0.0; 2287 if (numMessages != 0.0) avgMessLen = messLen/numMessages; else avgMessLen = 0.0; 2288 if (messageLength != 0.0) fracLength = messLen/messageLength; else fracLength = 0.0; 2289 if (numReductions != 0.0) fracReductions = red/numReductions; else fracReductions = 0.0; 2290 ierr = PetscFPrintf(comm, fd, "# "); 2291 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", 2292 stage, name, stageTime/size, 100.0*fracTime, flops, 100.0*fracFlops, 2293 mess, 100.0*fracMessages, avgMessLen, 100.0*fracLength, red, 100.0*fracReductions);CHKERRQ(ierr); 2294 } 2295 } 2296 2297 /* Report events */ 2298 ierr = PetscFPrintf(comm,fd,"\n# Event\n");CHKERRQ(ierr); 2299 ierr = PetscFPrintf(comm,fd,"# ------------------------------------------------------\n");CHKERRQ(ierr); 2300 ierr = PetscFPrintf(comm,fd,"class Stage(object):\n");CHKERRQ(ierr); 2301 ierr = PetscFPrintf(comm,fd," def __init__(self, name, time, flops, numMessages, messageLength, numReductions):\n");CHKERRQ(ierr); 2302 ierr = PetscFPrintf(comm,fd," # The time and flops represent totals across processes, whereas reductions are only counted once\n");CHKERRQ(ierr); 2303 ierr = PetscFPrintf(comm,fd," self.name = name\n");CHKERRQ(ierr); 2304 ierr = PetscFPrintf(comm,fd," self.time = time\n");CHKERRQ(ierr); 2305 ierr = PetscFPrintf(comm,fd," self.flops = flops\n");CHKERRQ(ierr); 2306 ierr = PetscFPrintf(comm,fd," self.numMessages = numMessages\n");CHKERRQ(ierr); 2307 ierr = PetscFPrintf(comm,fd," self.messageLength = messageLength\n");CHKERRQ(ierr); 2308 ierr = PetscFPrintf(comm,fd," self.numReductions = numReductions\n");CHKERRQ(ierr); 2309 ierr = PetscFPrintf(comm,fd," self.event = {}\n");CHKERRQ(ierr); 2310 ierr = PetscFPrintf(comm,fd, "class Dummy(object):\n");CHKERRQ(ierr); 2311 ierr = PetscFPrintf(comm,fd, " pass\n");CHKERRQ(ierr); 2312 /* Problem: The stage name will not show up unless the stage executed on proc 1 */ 2313 for (stage = 0; stage < numStages; stage++) { 2314 if (!stageVisible[stage]) continue; 2315 if (localStageUsed[stage]) { 2316 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2317 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2318 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2319 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2320 ierr = MPI_Allreduce(&stageInfo[stage].perfInfo.numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2321 } else { 2322 ierr = PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage);CHKERRQ(ierr); 2323 ierr = MPI_Allreduce(&zero, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2324 ierr = MPI_Allreduce(&zero, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2325 ierr = MPI_Allreduce(&zero, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2326 ierr = MPI_Allreduce(&zero, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2327 ierr = MPI_Allreduce(&zero, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm);CHKERRQ(ierr); 2328 } 2329 mess *= 0.5; messLen *= 0.5; red /= size; 2330 2331 /* Get total number of events in this stage -- 2332 Currently, a single processor can register more events than another, but events must all be registered in order, 2333 just like stages. We can removed this requirement if necessary by having a global event numbering and indirection 2334 on the event ID. This seems best accomplished by assoicating a communicator with each stage. 2335 2336 Problem: If the event did not happen on proc 1, its name will not be available. 2337 Problem: Event visibility is not implemented 2338 */ 2339 2340 { 2341 size_t len, c; 2342 2343 ierr = PetscStrcpy(stageName, stageInfo[stage].name);CHKERRQ(ierr); 2344 ierr = PetscStrlen(stageName, &len);CHKERRQ(ierr); 2345 for (c = 0; c < len; ++c) { 2346 if (stageName[c] == ' ') stageName[c] = '_'; 2347 } 2348 } 2349 if (!rank) { 2350 ierr = PetscFPrintf(comm, fd, "%s = Stage('%s', %g, %g, %g, %g, %g)\n", stageName, stageName, stageTime, flops, mess, messLen, red);CHKERRQ(ierr); 2351 } 2352 2353 if (localStageUsed[stage]) { 2354 eventInfo = stageLog->stageInfo[stage].eventLog->eventInfo; 2355 localNumEvents = stageLog->stageInfo[stage].eventLog->numEvents; 2356 } else localNumEvents = 0; 2357 2358 ierr = MPI_Allreduce(&localNumEvents, &numEvents, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 2359 for (event = 0; event < numEvents; event++) { 2360 PetscBool hasEvent = PETSC_TRUE; 2361 PetscMPIInt tmpI; 2362 PetscLogDouble tmpR; 2363 2364 if (localStageUsed[stage] && (event < stageLog->stageInfo[stage].eventLog->numEvents) && (eventInfo[event].depth == 0)) { 2365 size_t len, c; 2366 2367 ierr = MPI_Allreduce(&eventInfo[event].count, &maxCt, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 2368 ierr = PetscStrcpy(eventName, stageLog->eventLog->eventInfo[event].name);CHKERRQ(ierr); 2369 ierr = PetscStrlen(eventName, &len);CHKERRQ(ierr); 2370 for (c = 0; c < len; ++c) { 2371 if (eventName[c] == ' ') eventName[c] = '_'; 2372 } 2373 } else { 2374 ierr = MPI_Allreduce(&ierr, &maxCt, 1, MPI_INT, MPI_MAX, comm);CHKERRQ(ierr); 2375 eventName[0] = 0; 2376 hasEvent = PETSC_FALSE; 2377 } 2378 2379 if (maxCt != 0) { 2380 ierr = PetscFPrintf(comm, fd,"#\n");CHKERRQ(ierr); 2381 if (!rank) { 2382 ierr = PetscFPrintf(comm, fd, "%s = Dummy()\n",eventName);CHKERRQ(ierr); 2383 ierr = PetscFPrintf(comm, fd, "%s.event['%s'] = %s\n",stageName,eventName,eventName);CHKERRQ(ierr); 2384 } 2385 /* Count */ 2386 if (hasEvent) tmpI = eventInfo[event].count; 2387 else tmpI = 0; 2388 ierr = MPI_Gather(&tmpI,1, MPI_INT, mycount, 1, MPI_INT, 0, comm);CHKERRQ(ierr); 2389 ierr = PetscFPrintf(comm, fd, "%s.Count = [ ", eventName);CHKERRQ(ierr); 2390 if (!rank) { 2391 for (i=0; i<size; i++) { 2392 ierr = PetscFPrintf(comm, fd, " %7d,",mycount[i]);CHKERRQ(ierr); 2393 } 2394 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2395 } 2396 /* Time */ 2397 if (hasEvent) tmpR = eventInfo[event].time; 2398 else tmpR = 0.0; 2399 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2400 if (!rank) { 2401 ierr = PetscFPrintf(comm, fd, "%s.Time = [ ", eventName);CHKERRQ(ierr); 2402 for (i=0; i<size; i++) { 2403 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2404 } 2405 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2406 } 2407 if (hasEvent) tmpR = eventInfo[event].time2; 2408 else tmpR = 0.0; 2409 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2410 if (!rank) { 2411 ierr = PetscFPrintf(comm, fd, "%s.Time2 = [ ", eventName);CHKERRQ(ierr); 2412 for (i=0; i<size; i++) { 2413 ierr = PetscFPrintf(comm, fd, " %5.3e,", mydata[i]);CHKERRQ(ierr); 2414 } 2415 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2416 } 2417 /* Flops */ 2418 if (hasEvent) tmpR = eventInfo[event].flops; 2419 else tmpR = 0.0; 2420 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2421 if (!rank) { 2422 ierr = PetscFPrintf(comm, fd, "%s.Flops = [ ", eventName);CHKERRQ(ierr); 2423 for (i=0; i<size; i++) { 2424 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2425 } 2426 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2427 } 2428 if (hasEvent) tmpR = eventInfo[event].flops2; 2429 else tmpR = 0.0; 2430 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2431 if (!rank) { 2432 ierr = PetscFPrintf(comm, fd, "%s.Flops2 = [ ", eventName);CHKERRQ(ierr); 2433 for (i=0; i<size; i++) { 2434 ierr = PetscFPrintf(comm, fd, " %5.3e,", mydata[i]);CHKERRQ(ierr); 2435 } 2436 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2437 } 2438 /* Num Messages */ 2439 if (hasEvent) tmpR = eventInfo[event].numMessages; 2440 else tmpR = 0.0; 2441 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2442 ierr = PetscFPrintf(comm, fd, "%s.NumMessages = [ ", eventName);CHKERRQ(ierr); 2443 if (!rank) { 2444 for (i=0; i<size; i++) { 2445 ierr = PetscFPrintf(comm, fd, " %7.1e,",mydata[i]);CHKERRQ(ierr); 2446 } 2447 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2448 } 2449 /* Message Length */ 2450 if (hasEvent) tmpR = eventInfo[event].messageLength; 2451 else tmpR = 0.0; 2452 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2453 if (!rank) { 2454 ierr = PetscFPrintf(comm, fd, "%s.MessageLength = [ ", eventName);CHKERRQ(ierr); 2455 for (i=0; i<size; i++) { 2456 ierr = PetscFPrintf(comm, fd, " %5.3e,",mydata[i]);CHKERRQ(ierr); 2457 } 2458 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2459 } 2460 /* Num Reductions */ 2461 if (hasEvent) tmpR = eventInfo[event].numReductions; 2462 else tmpR = 0.0; 2463 ierr = MPI_Gather(&tmpR, 1, MPIU_PETSCLOGDOUBLE, mydata, 1, MPIU_PETSCLOGDOUBLE, 0, comm);CHKERRQ(ierr); 2464 ierr = PetscFPrintf(comm, fd, "%s.NumReductions = [ ", eventName);CHKERRQ(ierr); 2465 if (!rank) { 2466 for (i=0; i<size; i++) { 2467 ierr = PetscFPrintf(comm, fd, " %7.1e,",mydata[i]);CHKERRQ(ierr); 2468 } 2469 ierr = PetscFPrintf(comm, fd, "]\n");CHKERRQ(ierr); 2470 } 2471 } 2472 } 2473 } 2474 2475 /* Right now, only stages on the first processor are reported here, meaning only objects associated with 2476 the global communicator, or MPI_COMM_SELF for proc 1. We really should report global stats and then 2477 stats for stages local to processor sets. 2478 */ 2479 for (stage = 0; stage < numStages; stage++) { 2480 if (!localStageUsed[stage]) { 2481 ierr = PetscFPrintf(comm, fd, "\n--- Event Stage %d: Unknown\n\n", stage);CHKERRQ(ierr); 2482 } 2483 } 2484 2485 ierr = PetscFree(localStageUsed);CHKERRQ(ierr); 2486 ierr = PetscFree(stageUsed);CHKERRQ(ierr); 2487 ierr = PetscFree(localStageVisible);CHKERRQ(ierr); 2488 ierr = PetscFree(stageVisible);CHKERRQ(ierr); 2489 ierr = PetscFree(mydata);CHKERRQ(ierr); 2490 ierr = PetscFree(mycount);CHKERRQ(ierr); 2491 2492 /* Information unrelated to this particular run */ 2493 ierr = PetscFPrintf(comm, fd, "# ========================================================================================================================\n");CHKERRQ(ierr); 2494 PetscTime(&y); 2495 PetscTime(&x); 2496 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 2497 PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); PetscTime(&y); 2498 ierr = PetscFPrintf(comm,fd,"AveragetimetogetPetscTime = %g\n", (y-x)/10.0);CHKERRQ(ierr); 2499 /* MPI information */ 2500 if (size > 1) { 2501 MPI_Status status; 2502 PetscMPIInt tag; 2503 MPI_Comm newcomm; 2504 2505 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2506 PetscTime(&x); 2507 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2508 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2509 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2510 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2511 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2512 PetscTime(&y); 2513 ierr = PetscFPrintf(comm, fd, "AveragetimeforMPI_Barrier = %g\n", (y-x)/5.0);CHKERRQ(ierr); 2514 ierr = PetscCommDuplicate(comm,&newcomm, &tag);CHKERRQ(ierr); 2515 ierr = MPI_Barrier(comm);CHKERRQ(ierr); 2516 if (rank) { 2517 ierr = MPI_Recv(0, 0, MPI_INT, rank-1, tag, newcomm, &status);CHKERRQ(ierr); 2518 ierr = MPI_Send(0, 0, MPI_INT, (rank+1)%size, tag, newcomm);CHKERRQ(ierr); 2519 } else { 2520 PetscTime(&x); 2521 ierr = MPI_Send(0, 0, MPI_INT, 1, tag, newcomm);CHKERRQ(ierr); 2522 ierr = MPI_Recv(0, 0, MPI_INT, size-1, tag, newcomm, &status);CHKERRQ(ierr); 2523 PetscTime(&y); 2524 ierr = PetscFPrintf(comm,fd,"AveragetimforzerosizeMPI_Send = %g\n", (y-x)/size);CHKERRQ(ierr); 2525 } 2526 ierr = PetscCommDestroy(&newcomm);CHKERRQ(ierr); 2527 } 2528 2529 /* Cleanup */ 2530 ierr = PetscFPrintf(comm, fd, "\n");CHKERRQ(ierr); 2531 ierr = PetscStageLogPush(stageLog, lastStage);CHKERRQ(ierr); 2532 PetscFunctionReturn(0); 2533 } 2534 2535 #else /* end of -DPETSC_USE_LOG section */ 2536 2537 #undef __FUNCT__ 2538 #define __FUNCT__ "PetscLogObjectState" 2539 PetscErrorCode PetscLogObjectState(PetscObject obj, const char format[], ...) 2540 { 2541 PetscFunctionBegin; 2542 PetscFunctionReturn(0); 2543 } 2544 2545 #endif /* PETSC_USE_LOG*/ 2546 2547 2548 PetscClassId PETSC_LARGEST_CLASSID = PETSC_SMALLEST_CLASSID; 2549 PetscClassId PETSC_OBJECT_CLASSID = 0; 2550 2551 #undef __FUNCT__ 2552 #define __FUNCT__ "PetscClassIdRegister" 2553 /*@C 2554 PetscClassIdRegister - Registers a new class name for objects and logging operations in an application code. 2555 2556 Not Collective 2557 2558 Input Parameter: 2559 . name - The class name 2560 2561 Output Parameter: 2562 . oclass - The class id or classid 2563 2564 Level: developer 2565 2566 .keywords: log, class, register 2567 2568 @*/ 2569 PetscErrorCode PetscClassIdRegister(const char name[],PetscClassId *oclass) 2570 { 2571 #if defined(PETSC_USE_LOG) 2572 PetscStageLog stageLog; 2573 PetscInt stage; 2574 PetscErrorCode ierr; 2575 #endif 2576 2577 PetscFunctionBegin; 2578 *oclass = ++PETSC_LARGEST_CLASSID; 2579 #if defined(PETSC_USE_LOG) 2580 ierr = PetscLogGetStageLog(&stageLog);CHKERRQ(ierr); 2581 ierr = PetscClassRegLogRegister(stageLog->classLog, name, *oclass);CHKERRQ(ierr); 2582 for (stage = 0; stage < stageLog->numStages; stage++) { 2583 ierr = ClassPerfLogEnsureSize(stageLog->stageInfo[stage].classLog, stageLog->classLog->numClasses);CHKERRQ(ierr); 2584 } 2585 #endif 2586 PetscFunctionReturn(0); 2587 } 2588 2589 #if defined(PETSC_USE_LOG) && defined(PETSC_HAVE_MPE) 2590 #include <mpe.h> 2591 2592 PetscBool PetscBeganMPE = PETSC_FALSE; 2593 2594 PETSC_INTERN PetscErrorCode PetscLogEventBeginMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject); 2595 PETSC_INTERN PetscErrorCode PetscLogEventEndMPE(PetscLogEvent,int,PetscObject,PetscObject,PetscObject,PetscObject); 2596 2597 #undef __FUNCT__ 2598 #define __FUNCT__ "PetscLogMPEBegin" 2599 /*@C 2600 PetscLogMPEBegin - Turns on MPE logging of events. This creates large log files 2601 and slows the program down. 2602 2603 Collective over PETSC_COMM_WORLD 2604 2605 Options Database Keys: 2606 . -log_mpe - Prints extensive log information (for code compiled with PETSC_USE_LOG) 2607 2608 Notes: 2609 A related routine is PetscLogBegin() (with the options key -log_summary), which is 2610 intended for production runs since it logs only flop rates and object 2611 creation (and should not significantly slow the programs). 2612 2613 Level: advanced 2614 2615 Concepts: logging^MPE 2616 Concepts: logging^message passing 2617 2618 .seealso: PetscLogDump(), PetscLogBegin(), PetscLogAllBegin(), PetscLogEventActivate(), 2619 PetscLogEventDeactivate() 2620 @*/ 2621 PetscErrorCode PetscLogMPEBegin(void) 2622 { 2623 PetscErrorCode ierr; 2624 2625 PetscFunctionBegin; 2626 /* Do MPE initialization */ 2627 if (!MPE_Initialized_logging()) { /* This function exists in mpich 1.1.2 and higher */ 2628 ierr = PetscInfo(0,"Initializing MPE.\n");CHKERRQ(ierr); 2629 ierr = MPE_Init_log();CHKERRQ(ierr); 2630 2631 PetscBeganMPE = PETSC_TRUE; 2632 } else { 2633 ierr = PetscInfo(0,"MPE already initialized. Not attempting to reinitialize.\n");CHKERRQ(ierr); 2634 } 2635 ierr = PetscLogSet(PetscLogEventBeginMPE, PetscLogEventEndMPE);CHKERRQ(ierr); 2636 PetscFunctionReturn(0); 2637 } 2638 2639 #undef __FUNCT__ 2640 #define __FUNCT__ "PetscLogMPEDump" 2641 /*@C 2642 PetscLogMPEDump - Dumps the MPE logging info to file for later use with Jumpshot. 2643 2644 Collective over PETSC_COMM_WORLD 2645 2646 Level: advanced 2647 2648 .seealso: PetscLogDump(), PetscLogAllBegin(), PetscLogMPEBegin() 2649 @*/ 2650 PetscErrorCode PetscLogMPEDump(const char sname[]) 2651 { 2652 char name[PETSC_MAX_PATH_LEN]; 2653 PetscErrorCode ierr; 2654 2655 PetscFunctionBegin; 2656 if (PetscBeganMPE) { 2657 ierr = PetscInfo(0,"Finalizing MPE.\n");CHKERRQ(ierr); 2658 if (sname) { 2659 ierr = PetscStrcpy(name,sname);CHKERRQ(ierr); 2660 } else { 2661 ierr = PetscGetProgramName(name,PETSC_MAX_PATH_LEN);CHKERRQ(ierr); 2662 } 2663 ierr = MPE_Finish_log(name);CHKERRQ(ierr); 2664 } else { 2665 ierr = PetscInfo(0,"Not finalizing MPE (not started by PETSc).\n");CHKERRQ(ierr); 2666 } 2667 PetscFunctionReturn(0); 2668 } 2669 2670 #define PETSC_RGB_COLORS_MAX 39 2671 static const char *PetscLogMPERGBColors[PETSC_RGB_COLORS_MAX] = { 2672 "OliveDrab: ", 2673 "BlueViolet: ", 2674 "CadetBlue: ", 2675 "CornflowerBlue: ", 2676 "DarkGoldenrod: ", 2677 "DarkGreen: ", 2678 "DarkKhaki: ", 2679 "DarkOliveGreen: ", 2680 "DarkOrange: ", 2681 "DarkOrchid: ", 2682 "DarkSeaGreen: ", 2683 "DarkSlateGray: ", 2684 "DarkTurquoise: ", 2685 "DeepPink: ", 2686 "DarkKhaki: ", 2687 "DimGray: ", 2688 "DodgerBlue: ", 2689 "GreenYellow: ", 2690 "HotPink: ", 2691 "IndianRed: ", 2692 "LavenderBlush: ", 2693 "LawnGreen: ", 2694 "LemonChiffon: ", 2695 "LightCoral: ", 2696 "LightCyan: ", 2697 "LightPink: ", 2698 "LightSalmon: ", 2699 "LightSlateGray: ", 2700 "LightYellow: ", 2701 "LimeGreen: ", 2702 "MediumPurple: ", 2703 "MediumSeaGreen: ", 2704 "MediumSlateBlue:", 2705 "MidnightBlue: ", 2706 "MintCream: ", 2707 "MistyRose: ", 2708 "NavajoWhite: ", 2709 "NavyBlue: ", 2710 "OliveDrab: " 2711 }; 2712 2713 #undef __FUNCT__ 2714 #define __FUNCT__ "PetscLogMPEGetRGBColor" 2715 /*@C 2716 PetscLogMPEGetRGBColor - This routine returns a rgb color useable with PetscLogEventRegister() 2717 2718 Not collective. Maybe it should be? 2719 2720 Output Parameter 2721 . str - character string representing the color 2722 2723 Level: developer 2724 2725 .keywords: log, mpe , color 2726 .seealso: PetscLogEventRegister 2727 @*/ 2728 PetscErrorCode PetscLogMPEGetRGBColor(const char *str[]) 2729 { 2730 static int idx = 0; 2731 2732 PetscFunctionBegin; 2733 *str = PetscLogMPERGBColors[idx]; 2734 idx = (idx + 1)% PETSC_RGB_COLORS_MAX; 2735 PetscFunctionReturn(0); 2736 } 2737 2738 #endif /* PETSC_USE_LOG && PETSC_HAVE_MPE */ 2739