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