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