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