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