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