1 #include <petscviewer.h> 2 #include <petscdevice.h> 3 #include <petsc/private/logimpl.h> /*I "petscsys.h" I*/ 4 #include <petsc/private/loghandlerimpl.h> 5 #include <petsc/private/deviceimpl.h> 6 #include <petscconfiginfo.h> 7 #include <petscmachineinfo.h> 8 #include "logdefault.h" 9 10 static PetscErrorCode PetscEventPerfInfoInit(PetscEventPerfInfo *eventInfo) 11 { 12 PetscFunctionBegin; 13 PetscCall(PetscMemzero(eventInfo, sizeof(*eventInfo))); 14 eventInfo->visible = PETSC_TRUE; 15 eventInfo->id = -1; 16 eventInfo->dof[0] = -1.0; 17 eventInfo->dof[1] = -1.0; 18 eventInfo->dof[2] = -1.0; 19 eventInfo->dof[3] = -1.0; 20 eventInfo->dof[4] = -1.0; 21 eventInfo->dof[5] = -1.0; 22 eventInfo->dof[6] = -1.0; 23 eventInfo->dof[7] = -1.0; 24 eventInfo->errors[0] = -1.0; 25 eventInfo->errors[1] = -1.0; 26 eventInfo->errors[2] = -1.0; 27 eventInfo->errors[3] = -1.0; 28 eventInfo->errors[4] = -1.0; 29 eventInfo->errors[5] = -1.0; 30 eventInfo->errors[6] = -1.0; 31 eventInfo->errors[7] = -1.0; 32 PetscFunctionReturn(PETSC_SUCCESS); 33 } 34 35 static PetscErrorCode PetscEventPerfInfoTic_Internal(PetscEventPerfInfo *eventInfo, PetscLogDouble time, PetscBool logMemory, int event, PetscBool resume) 36 { 37 PetscFunctionBegin; 38 if (resume) { 39 eventInfo->timeTmp -= time; 40 eventInfo->flopsTmp -= petsc_TotalFlops_th; 41 } else { 42 eventInfo->timeTmp = -time; 43 eventInfo->flopsTmp = -petsc_TotalFlops_th; 44 } 45 eventInfo->numMessages -= petsc_irecv_ct_th + petsc_isend_ct_th + petsc_recv_ct_th + petsc_send_ct_th; 46 eventInfo->messageLength -= petsc_irecv_len_th + petsc_isend_len_th + petsc_recv_len_th + petsc_send_len_th; 47 eventInfo->numReductions -= petsc_allreduce_ct_th + petsc_gather_ct_th + petsc_scatter_ct_th; 48 #if defined(PETSC_HAVE_DEVICE) 49 eventInfo->CpuToGpuCount -= petsc_ctog_ct_th; 50 eventInfo->GpuToCpuCount -= petsc_gtoc_ct_th; 51 eventInfo->CpuToGpuSize -= petsc_ctog_sz_th; 52 eventInfo->GpuToCpuSize -= petsc_gtoc_sz_th; 53 eventInfo->GpuFlops -= petsc_gflops_th; 54 eventInfo->GpuTime -= petsc_gtime; 55 #endif 56 if (logMemory) { 57 PetscLogDouble usage; 58 PetscCall(PetscMemoryGetCurrentUsage(&usage)); 59 eventInfo->memIncrease -= usage; 60 PetscCall(PetscMallocGetCurrentUsage(&usage)); 61 eventInfo->mallocSpace -= usage; 62 PetscCall(PetscMallocGetMaximumUsage(&usage)); 63 eventInfo->mallocIncrease -= usage; 64 PetscCall(PetscMallocPushMaximumUsage(event)); 65 } 66 PetscFunctionReturn(PETSC_SUCCESS); 67 } 68 69 static PetscErrorCode PetscEventPerfInfoTic(PetscEventPerfInfo *eventInfo, PetscLogDouble time, PetscBool logMemory, int event) 70 { 71 PetscFunctionBegin; 72 PetscCall(PetscEventPerfInfoTic_Internal(eventInfo, time, logMemory, event, PETSC_FALSE)); 73 PetscFunctionReturn(PETSC_SUCCESS); 74 } 75 76 static PetscErrorCode PetscEventPerfInfoResume(PetscEventPerfInfo *eventInfo, PetscLogDouble time, PetscBool logMemory, int event) 77 { 78 PetscFunctionBegin; 79 PetscCall(PetscEventPerfInfoTic_Internal(eventInfo, time, logMemory, event, PETSC_TRUE)); 80 PetscFunctionReturn(PETSC_SUCCESS); 81 } 82 83 static PetscErrorCode PetscEventPerfInfoToc_Internal(PetscEventPerfInfo *eventInfo, PetscLogDouble time, PetscBool logMemory, int event, PetscBool pause) 84 { 85 PetscFunctionBegin; 86 eventInfo->timeTmp += time; 87 eventInfo->flopsTmp += petsc_TotalFlops_th; 88 if (!pause) { 89 eventInfo->time += eventInfo->timeTmp; 90 eventInfo->time2 += eventInfo->timeTmp * eventInfo->timeTmp; 91 eventInfo->flops += eventInfo->flopsTmp; 92 eventInfo->flops2 += eventInfo->flopsTmp * eventInfo->flopsTmp; 93 } 94 eventInfo->numMessages += petsc_irecv_ct_th + petsc_isend_ct_th + petsc_recv_ct_th + petsc_send_ct_th; 95 eventInfo->messageLength += petsc_irecv_len_th + petsc_isend_len_th + petsc_recv_len + petsc_send_len_th; 96 eventInfo->numReductions += petsc_allreduce_ct_th + petsc_gather_ct_th + petsc_scatter_ct_th; 97 #if defined(PETSC_HAVE_DEVICE) 98 eventInfo->CpuToGpuCount += petsc_ctog_ct_th; 99 eventInfo->GpuToCpuCount += petsc_gtoc_ct_th; 100 eventInfo->CpuToGpuSize += petsc_ctog_sz_th; 101 eventInfo->GpuToCpuSize += petsc_gtoc_sz_th; 102 eventInfo->GpuFlops += petsc_gflops_th; 103 eventInfo->GpuTime += petsc_gtime; 104 #endif 105 if (logMemory) { 106 PetscLogDouble usage, musage; 107 PetscCall(PetscMemoryGetCurrentUsage(&usage)); /* the comments below match the column labels printed in PetscLogView_Default() */ 108 eventInfo->memIncrease += usage; /* RMI */ 109 PetscCall(PetscMallocGetCurrentUsage(&usage)); 110 eventInfo->mallocSpace += usage; /* Malloc */ 111 PetscCall(PetscMallocPopMaximumUsage(event, &musage)); 112 eventInfo->mallocIncreaseEvent = PetscMax(musage - usage, eventInfo->mallocIncreaseEvent); /* EMalloc */ 113 PetscCall(PetscMallocGetMaximumUsage(&usage)); 114 eventInfo->mallocIncrease += usage; /* MMalloc */ 115 } 116 PetscFunctionReturn(PETSC_SUCCESS); 117 } 118 119 static PetscErrorCode PetscEventPerfInfoToc(PetscEventPerfInfo *eventInfo, PetscLogDouble time, PetscBool logMemory, int event) 120 { 121 PetscFunctionBegin; 122 PetscCall(PetscEventPerfInfoToc_Internal(eventInfo, time, logMemory, event, PETSC_FALSE)); 123 PetscFunctionReturn(PETSC_SUCCESS); 124 } 125 126 static PetscErrorCode PetscEventPerfInfoPause(PetscEventPerfInfo *eventInfo, PetscLogDouble time, PetscBool logMemory, int event) 127 { 128 PetscFunctionBegin; 129 PetscCall(PetscEventPerfInfoToc_Internal(eventInfo, time, logMemory, event, PETSC_TRUE)); 130 PetscFunctionReturn(PETSC_SUCCESS); 131 } 132 133 static PetscErrorCode PetscEventPerfInfoAdd_Internal(const PetscEventPerfInfo *eventInfo, PetscEventPerfInfo *outInfo) 134 { 135 PetscFunctionBegin; 136 outInfo->count += eventInfo->count; 137 outInfo->time += eventInfo->time; 138 outInfo->time2 += eventInfo->time2; 139 outInfo->flops += eventInfo->flops; 140 outInfo->flops2 += eventInfo->flops2; 141 outInfo->numMessages += eventInfo->numMessages; 142 outInfo->messageLength += eventInfo->messageLength; 143 outInfo->numReductions += eventInfo->numReductions; 144 #if defined(PETSC_HAVE_DEVICE) 145 outInfo->CpuToGpuCount += eventInfo->CpuToGpuCount; 146 outInfo->GpuToCpuCount += eventInfo->GpuToCpuCount; 147 outInfo->CpuToGpuSize += eventInfo->CpuToGpuSize; 148 outInfo->GpuToCpuSize += eventInfo->GpuToCpuSize; 149 outInfo->GpuFlops += eventInfo->GpuFlops; 150 outInfo->GpuTime += eventInfo->GpuTime; 151 #endif 152 outInfo->memIncrease += eventInfo->memIncrease; 153 outInfo->mallocSpace += eventInfo->mallocSpace; 154 outInfo->mallocIncreaseEvent += eventInfo->mallocIncreaseEvent; 155 outInfo->mallocIncrease += eventInfo->mallocIncrease; 156 PetscFunctionReturn(PETSC_SUCCESS); 157 } 158 159 PETSC_LOG_RESIZABLE_ARRAY(EventPerfArray, PetscEventPerfInfo, PetscLogEvent, PetscEventPerfInfoInit, NULL, NULL) 160 161 /* --- PetscClassPerf --- */ 162 163 typedef struct { 164 int creations; /* The number of objects of this class created */ 165 int destructions; /* The number of objects of this class destroyed */ 166 PetscLogDouble mem; /* The total memory allocated by objects of this class; this is completely wrong and should possibly be removed */ 167 PetscLogDouble descMem; /* The total memory allocated by descendents of these objects; this is completely wrong and should possibly be removed */ 168 } PetscClassPerf; 169 170 static PetscErrorCode PetscClassPerfInit(PetscClassPerf *classInfo) 171 { 172 PetscFunctionBegin; 173 PetscCall(PetscMemzero(classInfo, sizeof(*classInfo))); 174 PetscFunctionReturn(PETSC_SUCCESS); 175 } 176 177 PETSC_LOG_RESIZABLE_ARRAY(ClassPerfArray, PetscClassPerf, PetscLogClass, PetscClassPerfInit, NULL, NULL) 178 179 /* --- PetscStagePerf --- */ 180 181 typedef struct _PetscStagePerf { 182 PetscBool used; /* The stage was pushed on this processor */ 183 PetscEventPerfInfo perfInfo; /* The stage performance information */ 184 PetscLogEventPerfArray eventLog; /* The event information for this stage */ 185 PetscLogClassPerfArray classLog; /* The class information for this stage */ 186 } PetscStagePerf; 187 188 static PetscErrorCode PetscStageInfoInit(PetscStagePerf *stageInfo) 189 { 190 PetscFunctionBegin; 191 PetscCall(PetscMemzero(stageInfo, sizeof(*stageInfo))); 192 PetscCall(PetscLogEventPerfArrayCreate(128, &stageInfo->eventLog)); 193 PetscCall(PetscLogClassPerfArrayCreate(128, &stageInfo->classLog)); 194 PetscCall(PetscEventPerfInfoInit(&stageInfo->perfInfo)); 195 PetscFunctionReturn(PETSC_SUCCESS); 196 } 197 198 static PetscErrorCode PetscStageInfoReset(PetscStagePerf *stageInfo) 199 { 200 PetscFunctionBegin; 201 PetscCall(PetscLogEventPerfArrayDestroy(&stageInfo->eventLog)); 202 PetscCall(PetscLogClassPerfArrayDestroy(&stageInfo->classLog)); 203 PetscFunctionReturn(PETSC_SUCCESS); 204 } 205 206 PETSC_LOG_RESIZABLE_ARRAY(StageInfoArray, PetscStagePerf, PetscLogStage, PetscStageInfoInit, PetscStageInfoReset, NULL) 207 208 /* --- Action --- */ 209 210 /* The structure for action logging */ 211 typedef enum { 212 PETSC_LOG_ACTION_CREATE, 213 PETSC_LOG_ACTION_DESTROY, 214 PETSC_LOG_ACTION_BEGIN, 215 PETSC_LOG_ACTION_END, 216 } PetscLogActionType; 217 218 typedef struct _Action { 219 PetscLogActionType action; /* The type of execution */ 220 PetscLogEvent event; /* The event number */ 221 PetscClassId classid; /* The event class id */ 222 PetscLogDouble time; /* The time of occurrence */ 223 PetscLogDouble flops; /* The cumulative flops */ 224 PetscLogDouble mem; /* The current memory usage */ 225 PetscLogDouble maxmem; /* The maximum memory usage */ 226 PetscObjectId id1, id2, id3; /* The ids of associated objects */ 227 } Action; 228 229 PETSC_LOG_RESIZABLE_ARRAY(ActionArray, Action, PetscLogEvent, NULL, NULL, NULL) 230 231 /* --- Object --- */ 232 233 /* The structure for object logging */ 234 typedef struct _Object { 235 PetscObject obj; /* The associated PetscObject */ 236 int parent; /* The parent id */ 237 PetscLogDouble mem; /* The memory associated with the object */ 238 char name[64]; /* The object name */ 239 char info[64]; /* The information string */ 240 } Object; 241 242 PETSC_LOG_RESIZABLE_ARRAY(ObjectArray, Object, PetscObject, NULL, NULL, NULL) 243 244 /* Map from (threadid,stage,event) to perfInfo data struct */ 245 #include <petsc/private/hashmapijk.h> 246 247 PETSC_HASH_MAP(HMapEvent, PetscHashIJKKey, PetscEventPerfInfo *, PetscHashIJKKeyHash, PetscHashIJKKeyEqual, NULL) 248 249 typedef struct _n_PetscLogHandler_Default *PetscLogHandler_Default; 250 struct _n_PetscLogHandler_Default { 251 PetscLogStageInfoArray stages; 252 PetscSpinlock lock; 253 PetscLogActionArray petsc_actions; 254 PetscLogObjectArray petsc_objects; 255 PetscBool petsc_logActions; 256 PetscBool petsc_logObjects; 257 int petsc_numObjectsCreated; 258 int petsc_numObjectsDestroyed; 259 PetscHMapEvent eventInfoMap_th; 260 int pause_depth; 261 PetscBool use_threadsafe; 262 }; 263 264 /* --- PetscLogHandler_Default --- */ 265 266 static PetscErrorCode PetscLogHandlerContextCreate_Default(PetscLogHandler_Default *def_p) 267 { 268 PetscLogHandler_Default def; 269 270 PetscFunctionBegin; 271 PetscCall(PetscNew(def_p)); 272 def = *def_p; 273 PetscCall(PetscLogStageInfoArrayCreate(8, &def->stages)); 274 PetscCall(PetscLogActionArrayCreate(64, &def->petsc_actions)); 275 PetscCall(PetscLogObjectArrayCreate(64, &def->petsc_objects)); 276 PetscCall(PetscSpinlockCreate(&def->lock)); 277 278 PetscCall(PetscOptionsGetBool(NULL, NULL, "-log_include_actions", &def->petsc_logActions, NULL)); 279 PetscCall(PetscOptionsGetBool(NULL, NULL, "-log_include_objects", &def->petsc_logObjects, NULL)); 280 PetscCall(PetscOptionsGetBool(NULL, NULL, "-log_handler_default_use_threadsafe_events", &def->use_threadsafe, NULL)); 281 if (PetscDefined(HAVE_THREADSAFETY) || def->use_threadsafe) { PetscCall(PetscHMapEventCreate(&def->eventInfoMap_th)); } 282 PetscFunctionReturn(PETSC_SUCCESS); 283 } 284 285 static PetscErrorCode PetscLogHandlerDestroy_Default(PetscLogHandler h) 286 { 287 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 288 289 PetscFunctionBegin; 290 PetscCall(PetscLogStageInfoArrayDestroy(&def->stages)); 291 PetscCall(PetscLogActionArrayDestroy(&def->petsc_actions)); 292 PetscCall(PetscLogObjectArrayDestroy(&def->petsc_objects)); 293 PetscCall(PetscSpinlockDestroy(&def->lock)); 294 if (def->eventInfoMap_th) { 295 PetscEventPerfInfo **array; 296 PetscInt n, off = 0; 297 298 PetscCall(PetscHMapEventGetSize(def->eventInfoMap_th, &n)); 299 PetscCall(PetscMalloc1(n, &array)); 300 PetscCall(PetscHMapEventGetVals(def->eventInfoMap_th, &off, array)); 301 for (PetscInt i = 0; i < n; i++) PetscCall(PetscFree(array[i])); 302 PetscCall(PetscFree(array)); 303 PetscCall(PetscHMapEventDestroy(&def->eventInfoMap_th)); 304 } 305 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerGetEventPerfInfo_C", NULL)); 306 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerGetStagePerfInfo_C", NULL)); 307 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerSetLogActions_C", NULL)); 308 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerSetLogObjects_C", NULL)); 309 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerLogObjectState_C", NULL)); 310 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerGetNumObjects_C", NULL)); 311 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerEventDeactivatePush_C", NULL)); 312 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerEventDeactivatePop_C", NULL)); 313 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerEventsPause_C", NULL)); 314 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerEventsResume_C", NULL)); 315 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerDump_C", NULL)); 316 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerStageSetVisible_C", NULL)); 317 PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerStageGetVisible_C", NULL)); 318 PetscCall(PetscFree(def)); 319 PetscFunctionReturn(PETSC_SUCCESS); 320 } 321 322 static PetscErrorCode PetscLogHandlerDefaultGetStageInfo(PetscLogHandler handler, PetscLogStage stage, PetscStagePerf **stage_info_p) 323 { 324 PetscStagePerf *stage_info = NULL; 325 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 326 327 PetscFunctionBegin; 328 PetscCall(PetscLogStageInfoArrayResize(def->stages, stage + 1)); 329 PetscCall(PetscLogStageInfoArrayGetRef(def->stages, stage, &stage_info)); 330 *stage_info_p = stage_info; 331 PetscFunctionReturn(PETSC_SUCCESS); 332 } 333 334 static PetscErrorCode PetscLogHandlerGetEventPerfInfo_Default(PetscLogHandler handler, PetscLogStage stage, PetscLogEvent event, PetscEventPerfInfo **event_info_p) 335 { 336 PetscEventPerfInfo *event_info = NULL; 337 PetscStagePerf *stage_info = NULL; 338 PetscLogEventPerfArray event_log; 339 340 PetscFunctionBegin; 341 if (stage < 0) PetscCall(PetscLogStateGetCurrentStage(handler->state, &stage)); 342 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage, &stage_info)); 343 event_log = stage_info->eventLog; 344 PetscCall(PetscLogEventPerfArrayResize(event_log, event + 1)); 345 PetscCall(PetscLogEventPerfArrayGetRef(event_log, event, &event_info)); 346 event_info->id = event; 347 *event_info_p = event_info; 348 PetscFunctionReturn(PETSC_SUCCESS); 349 } 350 351 static PetscErrorCode PetscLogHandlerGetStagePerfInfo_Default(PetscLogHandler handler, PetscLogStage stage, PetscEventPerfInfo **stage_info_p) 352 { 353 PetscStagePerf *stage_perf_info = NULL; 354 355 PetscFunctionBegin; 356 if (stage < 0) PetscCall(PetscLogStateGetCurrentStage(handler->state, &stage)); 357 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage, &stage_perf_info)); 358 *stage_info_p = &stage_perf_info->perfInfo; 359 PetscFunctionReturn(PETSC_SUCCESS); 360 } 361 362 static PetscErrorCode PetscLogHandlerDefaultGetClassPerf(PetscLogHandler handler, PetscLogStage stage, PetscLogClass clss, PetscClassPerf **class_info) 363 { 364 PetscLogClassPerfArray class_log; 365 PetscStagePerf *stage_info; 366 367 PetscFunctionBegin; 368 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage, &stage_info)); 369 class_log = stage_info->classLog; 370 PetscCall(PetscLogClassPerfArrayResize(class_log, clss + 1)); 371 PetscCall(PetscLogClassPerfArrayGetRef(class_log, clss, class_info)); 372 PetscFunctionReturn(PETSC_SUCCESS); 373 } 374 375 static PetscErrorCode PetscLogHandlerObjectCreate_Default(PetscLogHandler h, PetscObject obj) 376 { 377 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 378 PetscLogState state; 379 PetscLogStage stage; 380 PetscClassPerf *classInfo; 381 int oclass = 0; 382 383 PetscFunctionBegin; 384 PetscCall(PetscLogHandlerGetState(h, &state)); 385 PetscCall(PetscSpinlockLock(&def->lock)); 386 /* Record stage info */ 387 PetscCall(PetscLogStateGetCurrentStage(state, &stage)); 388 PetscCall(PetscLogStateGetClassFromClassId(state, obj->classid, &oclass)); 389 PetscCall(PetscLogHandlerDefaultGetClassPerf(h, stage, oclass, &classInfo)); 390 classInfo->creations++; 391 /* Record the creation action */ 392 if (def->petsc_logActions) { 393 Action new_action; 394 395 PetscCall(PetscTime(&new_action.time)); 396 new_action.time -= petsc_BaseTime; 397 new_action.action = PETSC_LOG_ACTION_CREATE; 398 new_action.event = -1; 399 new_action.classid = obj->classid; 400 new_action.id1 = obj->id; 401 new_action.id2 = -1; 402 new_action.id3 = -1; 403 new_action.flops = petsc_TotalFlops; 404 PetscCall(PetscMallocGetCurrentUsage(&new_action.mem)); 405 PetscCall(PetscMallocGetMaximumUsage(&new_action.maxmem)); 406 PetscCall(PetscLogActionArrayPush(def->petsc_actions, new_action)); 407 } 408 /* We don't just use obj->id to count all objects that are created 409 because PetscLogHandlers are objects and PetscLogObjectDestroy() will not 410 be called for them: the number of objects created and destroyed as counted 411 here and below would have an imbalance */ 412 def->petsc_numObjectsCreated++; 413 /* Record the object */ 414 if (def->petsc_logObjects) { 415 Object new_object; 416 PetscInt objid; 417 418 new_object.parent = -1; 419 new_object.obj = obj; 420 new_object.mem = 0; 421 422 PetscCall(PetscMemzero(new_object.name, sizeof(new_object.name))); 423 PetscCall(PetscMemzero(new_object.info, sizeof(new_object.info))); 424 PetscAssert(obj->id >= 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Object ids from PetscObjectNewId_Internal() start at 1"); 425 PetscCall(PetscIntCast(obj->id, &objid)); 426 PetscCall(PetscLogObjectArrayResize(def->petsc_objects, objid)); 427 PetscCall(PetscLogObjectArraySet(def->petsc_objects, objid - 1, new_object)); 428 } 429 PetscCall(PetscSpinlockUnlock(&def->lock)); 430 PetscFunctionReturn(PETSC_SUCCESS); 431 } 432 433 static PetscErrorCode PetscLogHandlerObjectDestroy_Default(PetscLogHandler h, PetscObject obj) 434 { 435 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 436 PetscLogState state; 437 PetscLogStage stage; 438 PetscClassPerf *classInfo; 439 int oclass = 0; 440 441 PetscFunctionBegin; 442 PetscCall(PetscLogHandlerGetState(h, &state)); 443 /* Record stage info */ 444 PetscCall(PetscSpinlockLock(&def->lock)); 445 PetscCall(PetscLogStateGetCurrentStage(state, &stage)); 446 if (stage >= 0) { 447 /* stage < 0 can happen if the log summary is output before some things are destroyed */ 448 PetscCall(PetscLogStateGetClassFromClassId(state, obj->classid, &oclass)); 449 PetscCall(PetscLogHandlerDefaultGetClassPerf(h, stage, oclass, &classInfo)); 450 classInfo->destructions++; 451 } 452 /* Cannot Credit all ancestors with your memory because they may have already been destroyed*/ 453 def->petsc_numObjectsDestroyed++; 454 /* Dynamically enlarge logging structures */ 455 /* Record the destruction action */ 456 if (def->petsc_logActions) { 457 Action new_action; 458 459 PetscCall(PetscTime(&new_action.time)); 460 new_action.time -= petsc_BaseTime; 461 new_action.event = -1; 462 new_action.action = PETSC_LOG_ACTION_DESTROY; 463 new_action.classid = obj->classid; 464 new_action.id1 = obj->id; 465 new_action.id2 = -1; 466 new_action.id3 = -1; 467 new_action.flops = petsc_TotalFlops; 468 PetscCall(PetscMallocGetCurrentUsage(&new_action.mem)); 469 PetscCall(PetscMallocGetMaximumUsage(&new_action.maxmem)); 470 PetscCall(PetscLogActionArrayPush(def->petsc_actions, new_action)); 471 } 472 if (def->petsc_logObjects) { 473 Object *obj_entry = NULL; 474 PetscInt objid; 475 476 PetscAssert(obj->id >= 1, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Object ids from PetscObjectNewId_Internal() start at 1"); 477 PetscCall(PetscIntCast(obj->id - 1, &objid)); 478 PetscCall(PetscLogObjectArrayGetRef(def->petsc_objects, objid, &obj_entry)); 479 if (obj->name) PetscCall(PetscStrncpy(obj_entry->name, obj->name, 64)); 480 obj_entry->obj = NULL; 481 } 482 PetscCall(PetscSpinlockUnlock(&def->lock)); 483 PetscFunctionReturn(PETSC_SUCCESS); 484 } 485 486 static PetscErrorCode PetscLogHandlerEventSync_Default(PetscLogHandler h, PetscLogEvent event, MPI_Comm comm) 487 { 488 PetscLogState state; 489 PetscLogEventInfo event_info; 490 PetscEventPerfInfo *event_perf_info; 491 int stage; 492 PetscLogDouble time = 0.0; 493 494 PetscFunctionBegin; 495 if (!PetscLogSyncOn || comm == MPI_COMM_NULL) PetscFunctionReturn(PETSC_SUCCESS); 496 PetscCall(PetscLogHandlerGetState(h, &state)); 497 PetscCall(PetscLogStateEventGetInfo(state, event, &event_info)); 498 if (!event_info.collective) PetscFunctionReturn(PETSC_SUCCESS); 499 PetscCall(PetscLogStateGetCurrentStage(state, &stage)); 500 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(h, stage, event, &event_perf_info)); 501 if (event_perf_info->depth > 0) PetscFunctionReturn(PETSC_SUCCESS); 502 503 PetscCall(PetscTimeSubtract(&time)); 504 PetscCallMPI(MPI_Barrier(comm)); 505 PetscCall(PetscTimeAdd(&time)); 506 event_perf_info->syncTime += time; 507 PetscFunctionReturn(PETSC_SUCCESS); 508 } 509 510 static PetscErrorCode PetscLogGetStageEventPerfInfo_threaded(PetscLogHandler_Default def, PetscLogStage stage, PetscLogEvent event, PetscEventPerfInfo **eventInfo) 511 { 512 PetscEventPerfInfo *leventInfo = NULL; 513 PetscHashIJKKey key; 514 515 PetscFunctionBegin; 516 #if PetscDefined(HAVE_THREADSAFETY) 517 key.i = PetscLogGetTid(); 518 #else 519 key.i = 0; 520 #endif 521 key.j = stage; 522 key.k = event; 523 PetscCall(PetscSpinlockLock(&def->lock)); 524 PetscCall(PetscHMapEventGet(def->eventInfoMap_th, key, &leventInfo)); 525 if (!leventInfo) { 526 PetscCall(PetscNew(&leventInfo)); 527 leventInfo->id = event; 528 PetscCall(PetscHMapEventSet(def->eventInfoMap_th, key, leventInfo)); 529 } 530 PetscCall(PetscSpinlockUnlock(&def->lock)); 531 *eventInfo = leventInfo; 532 PetscFunctionReturn(PETSC_SUCCESS); 533 } 534 535 static PetscErrorCode PetscLogHandlerEventBegin_Default(PetscLogHandler h, PetscLogEvent event, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4) 536 { 537 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 538 PetscEventPerfInfo *event_perf_info = NULL; 539 PetscLogEventInfo event_info; 540 PetscLogDouble time; 541 PetscLogState state; 542 PetscLogStage stage; 543 544 PetscFunctionBegin; 545 PetscCall(PetscLogHandlerGetState(h, &state)); 546 if (PetscDefined(USE_DEBUG)) { 547 PetscCall(PetscLogStateEventGetInfo(state, event, &event_info)); 548 if (PetscUnlikely(o1)) PetscValidHeader(o1, 3); 549 if (PetscUnlikely(o2)) PetscValidHeader(o2, 4); 550 if (PetscUnlikely(o3)) PetscValidHeader(o3, 5); 551 if (PetscUnlikely(o4)) PetscValidHeader(o4, 6); 552 if (event_info.collective && o1) { 553 PetscInt64 b1[2], b2[2]; 554 555 b1[0] = -o1->cidx; 556 b1[1] = o1->cidx; 557 PetscCallMPI(MPIU_Allreduce(b1, b2, 2, MPIU_INT64, MPI_MAX, PetscObjectComm(o1))); 558 PetscCheck(-b2[0] == b2[1], PETSC_COMM_SELF, PETSC_ERR_PLIB, "Collective event %s not called collectively %" PetscInt64_FMT " != %" PetscInt64_FMT, event_info.name, -b2[0], b2[1]); 559 } 560 } 561 /* Synchronization */ 562 PetscCall(PetscLogHandlerEventSync_Default(h, event, PetscObjectComm(o1))); 563 PetscCall(PetscLogStateGetCurrentStage(state, &stage)); 564 if (def->pause_depth > 0) stage = 0; // in pause-mode, all events run on the main stage 565 if (PetscDefined(HAVE_THREADSAFETY) || def->use_threadsafe) { 566 PetscCall(PetscLogGetStageEventPerfInfo_threaded(def, stage, event, &event_perf_info)); 567 if (event_perf_info->depth == 0) { PetscCall(PetscEventPerfInfoInit(event_perf_info)); } 568 } else { 569 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(h, stage, event, &event_perf_info)); 570 } 571 PetscCheck(event_perf_info->depth >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Trying to begin a paused event, this is not allowed"); 572 event_perf_info->depth++; 573 /* Check for double counting */ 574 if (event_perf_info->depth > 1) PetscFunctionReturn(PETSC_SUCCESS); 575 PetscCall(PetscLogStateEventGetInfo(state, event, &event_info)); 576 /* Log the performance info */ 577 event_perf_info->count++; 578 PetscCall(PetscTime(&time)); 579 PetscCall(PetscEventPerfInfoTic(event_perf_info, time, PetscLogMemory, event)); 580 if (def->petsc_logActions) { 581 PetscLogDouble curTime; 582 Action new_action; 583 584 PetscCall(PetscTime(&curTime)); 585 new_action.time = curTime - petsc_BaseTime; 586 new_action.action = PETSC_LOG_ACTION_BEGIN; 587 new_action.event = event; 588 new_action.classid = event_info.classid; 589 new_action.id1 = o1 ? o1->id : -1; 590 new_action.id2 = o2 ? o2->id : -1; 591 new_action.id3 = o3 ? o3->id : -1; 592 new_action.flops = petsc_TotalFlops; 593 PetscCall(PetscMallocGetCurrentUsage(&new_action.mem)); 594 PetscCall(PetscMallocGetMaximumUsage(&new_action.maxmem)); 595 PetscCall(PetscLogActionArrayPush(def->petsc_actions, new_action)); 596 } 597 PetscFunctionReturn(PETSC_SUCCESS); 598 } 599 600 static PetscErrorCode PetscLogHandlerEventEnd_Default(PetscLogHandler h, PetscLogEvent event, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4) 601 { 602 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 603 PetscEventPerfInfo *event_perf_info = NULL; 604 PetscLogDouble time; 605 PetscLogState state; 606 int stage; 607 PetscLogEventInfo event_info; 608 609 PetscFunctionBegin; 610 PetscCall(PetscLogHandlerGetState(h, &state)); 611 if (PetscDefined(USE_DEBUG)) { 612 PetscCall(PetscLogStateEventGetInfo(state, event, &event_info)); 613 if (PetscUnlikely(o1)) PetscValidHeader(o1, 3); 614 if (PetscUnlikely(o2)) PetscValidHeader(o2, 4); 615 if (PetscUnlikely(o3)) PetscValidHeader(o3, 5); 616 if (PetscUnlikely(o4)) PetscValidHeader(o4, 6); 617 if (event_info.collective && o1) { 618 PetscInt64 b1[2], b2[2]; 619 620 b1[0] = -o1->cidx; 621 b1[1] = o1->cidx; 622 PetscCallMPI(MPIU_Allreduce(b1, b2, 2, MPIU_INT64, MPI_MAX, PetscObjectComm(o1))); 623 PetscCheck(-b2[0] == b2[1], PETSC_COMM_SELF, PETSC_ERR_PLIB, "Collective event %s not called collectively %" PetscInt64_FMT " != %" PetscInt64_FMT, event_info.name, -b2[0], b2[1]); 624 } 625 } 626 if (def->petsc_logActions) { 627 PetscLogDouble curTime; 628 Action new_action; 629 630 PetscCall(PetscLogStateEventGetInfo(state, event, &event_info)); 631 PetscCall(PetscTime(&curTime)); 632 new_action.time = curTime - petsc_BaseTime; 633 new_action.action = PETSC_LOG_ACTION_END; 634 new_action.event = event; 635 new_action.classid = event_info.classid; 636 new_action.id1 = o1 ? o1->id : -1; 637 new_action.id2 = o2 ? o2->id : -2; 638 new_action.id3 = o3 ? o3->id : -3; 639 new_action.flops = petsc_TotalFlops; 640 PetscCall(PetscMallocGetCurrentUsage(&new_action.mem)); 641 PetscCall(PetscMallocGetMaximumUsage(&new_action.maxmem)); 642 PetscCall(PetscLogActionArrayPush(def->petsc_actions, new_action)); 643 } 644 PetscCall(PetscLogStateGetCurrentStage(state, &stage)); 645 if (def->pause_depth > 0) stage = 0; // all events run on the main stage in pause-mode 646 if (PetscDefined(HAVE_THREADSAFETY) || def->use_threadsafe) { 647 PetscCall(PetscLogGetStageEventPerfInfo_threaded(def, stage, event, &event_perf_info)); 648 } else { 649 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(h, stage, event, &event_perf_info)); 650 } 651 PetscCheck(event_perf_info->depth > 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Trying to end paused event, not allowed"); 652 event_perf_info->depth--; 653 /* Check for double counting */ 654 if (event_perf_info->depth > 0) PetscFunctionReturn(PETSC_SUCCESS); 655 else PetscCheck(event_perf_info->depth == 0, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Logging event had unbalanced begin/end pairs"); 656 657 /* Log performance info */ 658 PetscCall(PetscTime(&time)); 659 PetscCall(PetscEventPerfInfoToc(event_perf_info, time, PetscLogMemory, event)); 660 if (PetscDefined(HAVE_THREADSAFETY) || def->use_threadsafe) { 661 PetscEventPerfInfo *event_perf_info_global; 662 PetscCall(PetscSpinlockLock(&def->lock)); 663 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(h, stage, event, &event_perf_info_global)); 664 PetscCall(PetscEventPerfInfoAdd_Internal(event_perf_info, event_perf_info_global)); 665 PetscCall(PetscSpinlockUnlock(&def->lock)); 666 } 667 PetscFunctionReturn(PETSC_SUCCESS); 668 } 669 670 static PetscErrorCode PetscLogHandlerEventDeactivatePush_Default(PetscLogHandler h, PetscLogStage stage, PetscLogEvent event) 671 { 672 PetscEventPerfInfo *event_perf_info; 673 674 PetscFunctionBegin; 675 if (stage < 0) PetscCall(PetscLogStateGetCurrentStage(h->state, &stage)); 676 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(h, stage, event, &event_perf_info)); 677 event_perf_info->depth++; 678 PetscFunctionReturn(PETSC_SUCCESS); 679 } 680 681 static PetscErrorCode PetscLogHandlerEventDeactivatePop_Default(PetscLogHandler h, PetscLogStage stage, PetscLogEvent event) 682 { 683 PetscEventPerfInfo *event_perf_info; 684 685 PetscFunctionBegin; 686 if (stage < 0) PetscCall(PetscLogStateGetCurrentStage(h->state, &stage)); 687 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(h, stage, event, &event_perf_info)); 688 event_perf_info->depth--; 689 PetscFunctionReturn(PETSC_SUCCESS); 690 } 691 692 static PetscErrorCode PetscLogHandlerEventsPause_Default(PetscLogHandler h) 693 { 694 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 695 PetscLogDouble time; 696 PetscInt num_stages; 697 698 PetscFunctionBegin; 699 if (def->pause_depth++ > 0) PetscFunctionReturn(PETSC_SUCCESS); 700 PetscCall(PetscLogStageInfoArrayGetSize(def->stages, &num_stages, NULL)); 701 PetscCall(PetscTime(&time)); 702 for (PetscInt stage = 0; stage < num_stages; stage++) { 703 PetscStagePerf *stage_info = NULL; 704 PetscInt num_events; 705 706 PetscCall(PetscLogStageInfoArrayGetRef(def->stages, stage, &stage_info)); 707 PetscCall(PetscLogEventPerfArrayGetSize(stage_info->eventLog, &num_events, NULL)); 708 for (PetscInt event = 0; event < num_events; event++) { 709 PetscEventPerfInfo *event_info = NULL; 710 PetscCall(PetscLogEventPerfArrayGetRef(stage_info->eventLog, event, &event_info)); 711 if (event_info->depth > 0) { 712 event_info->depth *= -1; 713 PetscCall(PetscEventPerfInfoPause(event_info, time, PetscLogMemory, event)); 714 } 715 } 716 if (stage > 0 && stage_info->perfInfo.depth > 0) { 717 stage_info->perfInfo.depth *= -1; 718 PetscCall(PetscEventPerfInfoPause(&stage_info->perfInfo, time, PetscLogMemory, -(stage + 2))); 719 } 720 } 721 PetscFunctionReturn(PETSC_SUCCESS); 722 } 723 724 static PetscErrorCode PetscLogHandlerEventsResume_Default(PetscLogHandler h) 725 { 726 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 727 PetscLogDouble time; 728 PetscInt num_stages; 729 730 PetscFunctionBegin; 731 if (--def->pause_depth > 0) PetscFunctionReturn(PETSC_SUCCESS); 732 PetscCall(PetscLogStageInfoArrayGetSize(def->stages, &num_stages, NULL)); 733 PetscCall(PetscTime(&time)); 734 for (PetscInt stage = 0; stage < num_stages; stage++) { 735 PetscStagePerf *stage_info = NULL; 736 PetscInt num_events; 737 738 PetscCall(PetscLogStageInfoArrayGetRef(def->stages, stage, &stage_info)); 739 PetscCall(PetscLogEventPerfArrayGetSize(stage_info->eventLog, &num_events, NULL)); 740 for (PetscInt event = 0; event < num_events; event++) { 741 PetscEventPerfInfo *event_info = NULL; 742 PetscCall(PetscLogEventPerfArrayGetRef(stage_info->eventLog, event, &event_info)); 743 if (event_info->depth < 0) { 744 event_info->depth *= -1; 745 PetscCall(PetscEventPerfInfoResume(event_info, time, PetscLogMemory, event)); 746 } 747 } 748 if (stage > 0 && stage_info->perfInfo.depth < 0) { 749 stage_info->perfInfo.depth *= -1; 750 PetscCall(PetscEventPerfInfoResume(&stage_info->perfInfo, time, PetscLogMemory, -(stage + 2))); 751 } 752 } 753 PetscFunctionReturn(PETSC_SUCCESS); 754 } 755 756 static PetscErrorCode PetscLogHandlerStagePush_Default(PetscLogHandler h, PetscLogStage new_stage) 757 { 758 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 759 PetscLogDouble time; 760 PetscLogState state; 761 PetscLogStage current_stage; 762 PetscStagePerf *new_stage_info; 763 764 PetscFunctionBegin; 765 if (def->pause_depth > 0) PetscFunctionReturn(PETSC_SUCCESS); 766 PetscCall(PetscLogHandlerGetState(h, &state)); 767 current_stage = state->current_stage; 768 PetscCall(PetscLogHandlerDefaultGetStageInfo(h, new_stage, &new_stage_info)); 769 PetscCall(PetscTime(&time)); 770 771 /* Record flops/time of previous stage */ 772 if (current_stage >= 0) { 773 if (PetscBTLookup(state->active, current_stage)) { 774 PetscStagePerf *current_stage_info; 775 PetscCall(PetscLogHandlerDefaultGetStageInfo(h, current_stage, ¤t_stage_info)); 776 PetscCall(PetscEventPerfInfoToc(¤t_stage_info->perfInfo, time, PetscLogMemory, -(current_stage + 2))); 777 } 778 } 779 new_stage_info->used = PETSC_TRUE; 780 new_stage_info->perfInfo.count++; 781 new_stage_info->perfInfo.depth++; 782 /* Subtract current quantities so that we obtain the difference when we pop */ 783 if (PetscBTLookup(state->active, new_stage)) PetscCall(PetscEventPerfInfoTic(&new_stage_info->perfInfo, time, PetscLogMemory, -(new_stage + 2))); 784 PetscFunctionReturn(PETSC_SUCCESS); 785 } 786 787 static PetscErrorCode PetscLogHandlerStagePop_Default(PetscLogHandler h, PetscLogStage old_stage) 788 { 789 PetscLogHandler_Default def = (PetscLogHandler_Default)h->data; 790 PetscLogStage current_stage; 791 PetscStagePerf *old_stage_info; 792 PetscLogState state; 793 PetscLogDouble time; 794 795 PetscFunctionBegin; 796 if (def->pause_depth > 0) PetscFunctionReturn(PETSC_SUCCESS); 797 PetscCall(PetscLogHandlerGetState(h, &state)); 798 current_stage = state->current_stage; 799 PetscCall(PetscLogHandlerDefaultGetStageInfo(h, old_stage, &old_stage_info)); 800 PetscCall(PetscTime(&time)); 801 old_stage_info->perfInfo.depth--; 802 if (PetscBTLookup(state->active, old_stage)) { PetscCall(PetscEventPerfInfoToc(&old_stage_info->perfInfo, time, PetscLogMemory, -(old_stage + 2))); } 803 if (current_stage >= 0) { 804 if (PetscBTLookup(state->active, current_stage)) { 805 PetscStagePerf *current_stage_info; 806 PetscCall(PetscLogHandlerDefaultGetStageInfo(h, current_stage, ¤t_stage_info)); 807 PetscCall(PetscEventPerfInfoTic(¤t_stage_info->perfInfo, time, PetscLogMemory, -(current_stage + 2))); 808 } 809 } 810 PetscFunctionReturn(PETSC_SUCCESS); 811 } 812 813 static PetscErrorCode PetscLogHandlerStageSetVisible_Default(PetscLogHandler h, PetscLogStage stage, PetscBool is_visible) 814 { 815 PetscStagePerf *stage_info; 816 817 PetscFunctionBegin; 818 if (stage < 0) PetscCall(PetscLogStateGetCurrentStage(h->state, &stage)); 819 PetscCall(PetscLogHandlerDefaultGetStageInfo(h, stage, &stage_info)); 820 stage_info->perfInfo.visible = is_visible; 821 PetscFunctionReturn(PETSC_SUCCESS); 822 } 823 824 static PetscErrorCode PetscLogHandlerStageGetVisible_Default(PetscLogHandler h, PetscLogStage stage, PetscBool *is_visible) 825 { 826 PetscStagePerf *stage_info; 827 828 PetscFunctionBegin; 829 if (stage < 0) PetscCall(PetscLogStateGetCurrentStage(h->state, &stage)); 830 PetscCall(PetscLogHandlerDefaultGetStageInfo(h, stage, &stage_info)); 831 *is_visible = stage_info->perfInfo.visible; 832 PetscFunctionReturn(PETSC_SUCCESS); 833 } 834 835 static PetscErrorCode PetscLogHandlerSetLogActions_Default(PetscLogHandler handler, PetscBool flag) 836 { 837 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 838 839 PetscFunctionBegin; 840 def->petsc_logActions = flag; 841 PetscFunctionReturn(PETSC_SUCCESS); 842 } 843 844 static PetscErrorCode PetscLogHandlerSetLogObjects_Default(PetscLogHandler handler, PetscBool flag) 845 { 846 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 847 848 PetscFunctionBegin; 849 def->petsc_logObjects = flag; 850 PetscFunctionReturn(PETSC_SUCCESS); 851 } 852 853 static PetscErrorCode PetscLogHandlerLogObjectState_Default(PetscLogHandler handler, PetscObject obj, const char format[], va_list Argp) 854 { 855 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 856 size_t fullLength; 857 858 PetscFunctionBegin; 859 if (def->petsc_logObjects) { 860 Object *obj_entry = NULL; 861 PetscInt objid; 862 863 PetscCall(PetscIntCast(obj->id - 1, &objid)); 864 PetscCall(PetscLogObjectArrayGetRef(def->petsc_objects, objid, &obj_entry)); 865 PetscCall(PetscVSNPrintf(obj_entry->info, 64, format, &fullLength, Argp)); 866 } 867 PetscFunctionReturn(PETSC_SUCCESS); 868 } 869 870 static PetscErrorCode PetscLogHandlerGetNumObjects_Default(PetscLogHandler handler, PetscInt *num_objects) 871 { 872 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 873 874 PetscFunctionBegin; 875 PetscCall(PetscLogObjectArrayGetSize(def->petsc_objects, num_objects, NULL)); 876 PetscFunctionReturn(PETSC_SUCCESS); 877 } 878 879 static PetscErrorCode PetscLogHandlerDump_Default(PetscLogHandler handler, const char sname[]) 880 { 881 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 882 FILE *fd; 883 char file[PETSC_MAX_PATH_LEN], fname[PETSC_MAX_PATH_LEN]; 884 PetscLogDouble flops, _TotalTime; 885 PetscMPIInt rank; 886 int curStage; 887 PetscLogState state; 888 PetscInt num_events; 889 PetscLogEvent event; 890 891 PetscFunctionBegin; 892 /* Calculate the total elapsed time */ 893 PetscCall(PetscTime(&_TotalTime)); 894 _TotalTime -= petsc_BaseTime; 895 /* Open log file */ 896 PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)handler), &rank)); 897 PetscCall(PetscSNPrintf(file, PETSC_STATIC_ARRAY_LENGTH(file), "%s.%d", sname && sname[0] ? sname : "Log", rank)); 898 PetscCall(PetscFixFilename(file, fname)); 899 PetscCall(PetscFOpen(PETSC_COMM_SELF, fname, "w", &fd)); 900 PetscCheck(!(rank == 0) || !(!fd), PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Cannot open file: %s", fname); 901 /* Output totals */ 902 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Total Flop %14e %16.8e\n", petsc_TotalFlops, _TotalTime)); 903 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Clock Resolution %g\n", 0.0)); 904 /* Output actions */ 905 if (def->petsc_logActions) { 906 PetscInt num_actions; 907 PetscCall(PetscLogActionArrayGetSize(def->petsc_actions, &num_actions, NULL)); 908 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Actions accomplished %" PetscInt_FMT "\n", num_actions)); 909 for (int a = 0; a < num_actions; a++) { 910 Action *action; 911 912 PetscCall(PetscLogActionArrayGetRef(def->petsc_actions, a, &action)); 913 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "%g %d %d %d %" PetscInt64_FMT " %" PetscInt64_FMT " %" PetscInt64_FMT " %g %g %g\n", action->time, action->action, action->event, action->classid, action->id1, action->id2, action->id3, action->flops, 914 action->mem, action->maxmem)); 915 } 916 } 917 /* Output objects */ 918 if (def->petsc_logObjects) { 919 PetscInt num_objects; 920 921 PetscCall(PetscLogObjectArrayGetSize(def->petsc_objects, &num_objects, NULL)); 922 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Objects created %d destroyed %d\n", def->petsc_numObjectsCreated, def->petsc_numObjectsDestroyed)); 923 for (int o = 0; o < num_objects; o++) { 924 Object *object = NULL; 925 926 PetscCall(PetscLogObjectArrayGetRef(def->petsc_objects, o, &object)); 927 if (object->parent != -1) continue; // object with this id wasn't logged, probably a PetscLogHandler 928 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Parent ID: %d Memory: %d\n", object->parent, (int)object->mem)); 929 if (!object->name[0]) { 930 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "No Name\n")); 931 } else { 932 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Name: %s\n", object->name)); 933 } 934 if (!object->info[0]) { 935 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "No Info\n")); 936 } else { 937 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Info: %s\n", object->info)); 938 } 939 } 940 } 941 /* Output events */ 942 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "Event log:\n")); 943 PetscCall(PetscLogHandlerGetState(handler, &state)); 944 PetscCall(PetscLogStateGetNumEvents(state, &num_events)); 945 PetscCall(PetscLogStateGetCurrentStage(state, &curStage)); 946 for (event = 0; event < num_events; event++) { 947 PetscEventPerfInfo *event_info; 948 949 PetscCall(PetscLogHandlerGetEventPerfInfo_Default(handler, curStage, event, &event_info)); 950 if (event_info->time != 0.0) flops = event_info->flops / event_info->time; 951 else flops = 0.0; 952 PetscCall(PetscFPrintf(PETSC_COMM_SELF, fd, "%d %16d %16g %16g %16g\n", event, event_info->count, event_info->flops, event_info->time, flops)); 953 } 954 PetscCall(PetscFClose(PETSC_COMM_SELF, fd)); 955 PetscFunctionReturn(PETSC_SUCCESS); 956 } 957 958 /* 959 PetscLogView_Detailed - Each process prints the times for its own events 960 961 */ 962 static PetscErrorCode PetscLogHandlerView_Default_Detailed(PetscLogHandler handler, PetscViewer viewer) 963 { 964 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 965 PetscLogDouble locTotalTime, numRed, maxMem; 966 PetscInt numStages, numEvents; 967 MPI_Comm comm = PetscObjectComm((PetscObject)viewer); 968 PetscMPIInt rank, size; 969 PetscLogGlobalNames global_stages, global_events; 970 PetscLogState state; 971 PetscEventPerfInfo zero_info; 972 973 PetscFunctionBegin; 974 PetscCall(PetscLogHandlerGetState(handler, &state)); 975 PetscCallMPI(MPI_Comm_size(comm, &size)); 976 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 977 /* Must preserve reduction count before we go on */ 978 numRed = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 979 /* Get the total elapsed time */ 980 PetscCall(PetscTime(&locTotalTime)); 981 locTotalTime -= petsc_BaseTime; 982 PetscCall(PetscViewerASCIIPrintf(viewer, "size = %d\n", size)); 983 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalTimes = {}\n")); 984 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalMessages = {}\n")); 985 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalMessageLens = {}\n")); 986 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalReductions = {}\n")); 987 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalFlop = {}\n")); 988 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalObjects = {}\n")); 989 PetscCall(PetscViewerASCIIPrintf(viewer, "LocalMemory = {}\n")); 990 PetscCall(PetscLogRegistryCreateGlobalStageNames(comm, state->registry, &global_stages)); 991 PetscCall(PetscLogRegistryCreateGlobalEventNames(comm, state->registry, &global_events)); 992 PetscCall(PetscLogGlobalNamesGetSize(global_stages, NULL, &numStages)); 993 PetscCall(PetscLogGlobalNamesGetSize(global_events, NULL, &numEvents)); 994 PetscCall(PetscMemzero(&zero_info, sizeof(zero_info))); 995 PetscCall(PetscViewerASCIIPrintf(viewer, "Stages = {}\n")); 996 for (PetscInt stage = 0; stage < numStages; stage++) { 997 PetscInt stage_id; 998 const char *stage_name; 999 1000 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_stages, stage, &stage_id)); 1001 PetscCall(PetscLogGlobalNamesGlobalGetName(global_stages, stage, &stage_name)); 1002 PetscCall(PetscViewerASCIIPrintf(viewer, "Stages[\"%s\"] = {}\n", stage_name)); 1003 PetscCall(PetscViewerASCIIPrintf(viewer, "Stages[\"%s\"][\"summary\"] = {}\n", stage_name)); 1004 for (PetscInt event = 0; event < numEvents; event++) { 1005 PetscEventPerfInfo *eventInfo = &zero_info; 1006 PetscBool is_zero = PETSC_FALSE; 1007 PetscInt event_id; 1008 const char *event_name; 1009 1010 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_events, event, &event_id)); 1011 PetscCall(PetscLogGlobalNamesGlobalGetName(global_events, event, &event_name)); 1012 if (event_id >= 0 && stage_id >= 0) PetscCall(PetscLogHandlerGetEventPerfInfo_Default(handler, stage_id, event_id, &eventInfo)); 1013 is_zero = eventInfo->count == 0 ? PETSC_TRUE : PETSC_FALSE; 1014 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &is_zero, 1, MPIU_BOOL, MPI_LAND, comm)); 1015 if (!is_zero) { PetscCall(PetscViewerASCIIPrintf(viewer, "Stages[\"%s\"][\"%s\"] = {}\n", stage_name, event_name)); } 1016 } 1017 } 1018 PetscCall(PetscMallocGetMaximumUsage(&maxMem)); 1019 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 1020 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalTimes[%d] = %g\n", rank, locTotalTime)); 1021 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalMessages[%d] = %g\n", rank, petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct)); 1022 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalMessageLens[%d] = %g\n", rank, petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len)); 1023 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalReductions[%d] = %g\n", rank, numRed)); 1024 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalFlop[%d] = %g\n", rank, petsc_TotalFlops)); 1025 { 1026 PetscInt num_objects; 1027 1028 PetscCall(PetscLogObjectArrayGetSize(def->petsc_objects, &num_objects, NULL)); 1029 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalObjects[%d] = %" PetscInt_FMT "\n", rank, num_objects)); 1030 } 1031 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "LocalMemory[%d] = %g\n", rank, maxMem)); 1032 PetscCall(PetscViewerFlush(viewer)); 1033 for (PetscInt stage = 0; stage < numStages; stage++) { 1034 PetscEventPerfInfo *stage_perf_info = &zero_info; 1035 PetscInt stage_id; 1036 const char *stage_name; 1037 1038 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_stages, stage, &stage_id)); 1039 PetscCall(PetscLogGlobalNamesGlobalGetName(global_stages, stage, &stage_name)); 1040 if (stage_id >= 0) { 1041 PetscStagePerf *stage_info; 1042 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage_id, &stage_info)); 1043 stage_perf_info = &stage_info->perfInfo; 1044 } 1045 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "Stages[\"%s\"][\"summary\"][%d] = {\"time\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g}\n", stage_name, rank, stage_perf_info->time, 1046 stage_perf_info->numMessages, stage_perf_info->messageLength, stage_perf_info->numReductions, stage_perf_info->flops)); 1047 for (PetscInt event = 0; event < numEvents; event++) { 1048 PetscEventPerfInfo *eventInfo = &zero_info; 1049 PetscBool is_zero = PETSC_FALSE; 1050 PetscInt event_id; 1051 const char *event_name; 1052 1053 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_events, event, &event_id)); 1054 PetscCall(PetscLogGlobalNamesGlobalGetName(global_events, event, &event_name)); 1055 if (event_id >= 0 && stage_id >= 0) PetscCall(PetscLogHandlerGetEventPerfInfo_Default(handler, stage_id, event_id, &eventInfo)); 1056 is_zero = eventInfo->count == 0 ? PETSC_TRUE : PETSC_FALSE; 1057 PetscCall(PetscMemcmp(eventInfo, &zero_info, sizeof(zero_info), &is_zero)); 1058 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &is_zero, 1, MPIU_BOOL, MPI_LAND, comm)); 1059 if (!is_zero) { 1060 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "Stages[\"%s\"][\"%s\"][%d] = {\"count\" : %d, \"time\" : %g, \"syncTime\" : %g, \"numMessages\" : %g, \"messageLength\" : %g, \"numReductions\" : %g, \"flop\" : %g", stage_name, event_name, rank, 1061 eventInfo->count, eventInfo->time, eventInfo->syncTime, eventInfo->numMessages, eventInfo->messageLength, eventInfo->numReductions, eventInfo->flops)); 1062 if (eventInfo->dof[0] >= 0.) { 1063 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", \"dof\" : [")); 1064 for (PetscInt d = 0; d < 8; ++d) { 1065 if (d > 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", ")); 1066 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->dof[d])); 1067 } 1068 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "]")); 1069 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", \"error\" : [")); 1070 for (PetscInt e = 0; e < 8; ++e) { 1071 if (e > 0) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ", ")); 1072 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%g", eventInfo->errors[e])); 1073 } 1074 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "]")); 1075 } 1076 } 1077 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "}\n")); 1078 } 1079 } 1080 PetscCall(PetscViewerFlush(viewer)); 1081 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 1082 PetscCall(PetscLogGlobalNamesDestroy(&global_events)); 1083 PetscCall(PetscLogGlobalNamesDestroy(&global_stages)); 1084 PetscFunctionReturn(PETSC_SUCCESS); 1085 } 1086 1087 /* 1088 PetscLogView_CSV - Each process prints the times for its own events in Comma-Separated Value Format 1089 */ 1090 static PetscErrorCode PetscLogHandlerView_Default_CSV(PetscLogHandler handler, PetscViewer viewer) 1091 { 1092 PetscLogDouble locTotalTime, maxMem; 1093 PetscInt numStages, numEvents, stage, event; 1094 MPI_Comm comm = PetscObjectComm((PetscObject)viewer); 1095 PetscMPIInt rank, size; 1096 PetscLogGlobalNames global_stages, global_events; 1097 PetscLogState state; 1098 PetscEventPerfInfo zero_info; 1099 1100 PetscFunctionBegin; 1101 PetscCall(PetscLogHandlerGetState(handler, &state)); 1102 PetscCallMPI(MPI_Comm_size(comm, &size)); 1103 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1104 /* Must preserve reduction count before we go on */ 1105 /* Get the total elapsed time */ 1106 PetscCall(PetscTime(&locTotalTime)); 1107 locTotalTime -= petsc_BaseTime; 1108 PetscCall(PetscMallocGetMaximumUsage(&maxMem)); 1109 PetscCall(PetscLogRegistryCreateGlobalStageNames(comm, state->registry, &global_stages)); 1110 PetscCall(PetscLogRegistryCreateGlobalEventNames(comm, state->registry, &global_events)); 1111 PetscCall(PetscLogGlobalNamesGetSize(global_stages, NULL, &numStages)); 1112 PetscCall(PetscLogGlobalNamesGetSize(global_events, NULL, &numEvents)); 1113 PetscCall(PetscMemzero(&zero_info, sizeof(zero_info))); 1114 PetscCall(PetscViewerASCIIPushSynchronized(viewer)); 1115 PetscCall(PetscViewerASCIIPrintf(viewer, "Stage Name,Event Name,Rank,Count,Time,Num Messages,Message Length,Num Reductions,FLOP,dof0,dof1,dof2,dof3,dof4,dof5,dof6,dof7,e0,e1,e2,e3,e4,e5,e6,e7,%d\n", size)); 1116 PetscCall(PetscViewerFlush(viewer)); 1117 for (stage = 0; stage < numStages; stage++) { 1118 PetscEventPerfInfo *stage_perf_info; 1119 PetscInt stage_id; 1120 const char *stage_name; 1121 1122 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_stages, stage, &stage_id)); 1123 PetscCall(PetscLogGlobalNamesGlobalGetName(global_stages, stage, &stage_name)); 1124 stage_perf_info = &zero_info; 1125 if (stage_id >= 0) { 1126 PetscStagePerf *stage_info; 1127 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage_id, &stage_info)); 1128 stage_perf_info = &stage_info->perfInfo; 1129 } 1130 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%s,summary,%d,1,%g,%g,%g,%g,%g\n", stage_name, rank, stage_perf_info->time, stage_perf_info->numMessages, stage_perf_info->messageLength, stage_perf_info->numReductions, stage_perf_info->flops)); 1131 for (event = 0; event < numEvents; event++) { 1132 PetscEventPerfInfo *eventInfo = &zero_info; 1133 PetscBool is_zero = PETSC_FALSE; 1134 PetscInt event_id; 1135 const char *event_name; 1136 1137 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_events, event, &event_id)); 1138 PetscCall(PetscLogGlobalNamesGlobalGetName(global_events, event, &event_name)); 1139 if (event_id >= 0 && stage_id >= 0) PetscCall(PetscLogHandlerGetEventPerfInfo_Default(handler, stage_id, event_id, &eventInfo)); 1140 PetscCall(PetscMemcmp(eventInfo, &zero_info, sizeof(zero_info), &is_zero)); 1141 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &is_zero, 1, MPIU_BOOL, MPI_LAND, comm)); 1142 if (!is_zero) { 1143 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "%s,%s,%d,%d,%g,%g,%g,%g,%g", stage_name, event_name, rank, eventInfo->count, eventInfo->time, eventInfo->numMessages, eventInfo->messageLength, eventInfo->numReductions, eventInfo->flops)); 1144 if (eventInfo->dof[0] >= 0.) { 1145 for (PetscInt d = 0; d < 8; ++d) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->dof[d])); 1146 for (PetscInt e = 0; e < 8; ++e) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, ",%g", eventInfo->errors[e])); 1147 } 1148 } 1149 PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "\n")); 1150 } 1151 } 1152 PetscCall(PetscViewerFlush(viewer)); 1153 PetscCall(PetscViewerASCIIPopSynchronized(viewer)); 1154 PetscCall(PetscLogGlobalNamesDestroy(&global_stages)); 1155 PetscCall(PetscLogGlobalNamesDestroy(&global_events)); 1156 PetscFunctionReturn(PETSC_SUCCESS); 1157 } 1158 1159 static PetscErrorCode PetscLogViewWarnSync(PetscViewer viewer) 1160 { 1161 PetscFunctionBegin; 1162 if (!PetscLogSyncOn) PetscFunctionReturn(PETSC_SUCCESS); 1163 PetscCall(PetscViewerASCIIPrintf(viewer, "\n\n")); 1164 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n")); 1165 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1166 PetscCall(PetscViewerASCIIPrintf(viewer, " # WARNING!!! #\n")); 1167 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1168 PetscCall(PetscViewerASCIIPrintf(viewer, " # This program was run with logging synchronization. #\n")); 1169 PetscCall(PetscViewerASCIIPrintf(viewer, " # This option provides more meaningful imbalance #\n")); 1170 PetscCall(PetscViewerASCIIPrintf(viewer, " # figures at the expense of slowing things down and #\n")); 1171 PetscCall(PetscViewerASCIIPrintf(viewer, " # providing a distorted view of the overall runtime. #\n")); 1172 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1173 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n\n\n")); 1174 PetscFunctionReturn(PETSC_SUCCESS); 1175 } 1176 1177 static PetscErrorCode PetscLogViewWarnDebugging(PetscViewer viewer) 1178 { 1179 PetscFunctionBegin; 1180 if (PetscDefined(USE_DEBUG)) { 1181 PetscCall(PetscViewerASCIIPrintf(viewer, "\n\n")); 1182 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n")); 1183 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1184 PetscCall(PetscViewerASCIIPrintf(viewer, " # WARNING!!! #\n")); 1185 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1186 PetscCall(PetscViewerASCIIPrintf(viewer, " # This code was compiled with a debugging option. #\n")); 1187 PetscCall(PetscViewerASCIIPrintf(viewer, " # To get timing results run ./configure #\n")); 1188 PetscCall(PetscViewerASCIIPrintf(viewer, " # using --with-debugging=no, the performance will #\n")); 1189 PetscCall(PetscViewerASCIIPrintf(viewer, " # be generally two or three times faster. #\n")); 1190 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1191 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n\n\n")); 1192 } 1193 PetscFunctionReturn(PETSC_SUCCESS); 1194 } 1195 1196 static PetscErrorCode PetscLogViewWarnNoGpuAwareMpi(PetscViewer viewer) 1197 { 1198 #if defined(PETSC_HAVE_DEVICE) 1199 PetscMPIInt size; 1200 PetscBool deviceInitialized = PETSC_FALSE; 1201 1202 PetscFunctionBegin; 1203 PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)viewer), &size)); 1204 for (int i = PETSC_DEVICE_HOST + 1; i < PETSC_DEVICE_MAX; ++i) { 1205 const PetscDeviceType dtype = PetscDeviceTypeCast(i); 1206 if (PetscDeviceInitialized(dtype)) { /* a non-host device was initialized */ 1207 deviceInitialized = PETSC_TRUE; 1208 break; 1209 } 1210 } 1211 /* the last condition says PETSc is configured with device but it is a pure CPU run, so don't print misleading warnings */ 1212 if (use_gpu_aware_mpi || size == 1 || !deviceInitialized) PetscFunctionReturn(PETSC_SUCCESS); 1213 PetscCall(PetscViewerASCIIPrintf(viewer, "\n\n")); 1214 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n")); 1215 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1216 PetscCall(PetscViewerASCIIPrintf(viewer, " # WARNING!!! #\n")); 1217 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1218 PetscCall(PetscViewerASCIIPrintf(viewer, " # This code was compiled with GPU support and you've #\n")); 1219 PetscCall(PetscViewerASCIIPrintf(viewer, " # created PETSc/GPU objects, but you intentionally #\n")); 1220 PetscCall(PetscViewerASCIIPrintf(viewer, " # used -use_gpu_aware_mpi 0, requiring PETSc to copy #\n")); 1221 PetscCall(PetscViewerASCIIPrintf(viewer, " # additional data between the GPU and CPU. To obtain #\n")); 1222 PetscCall(PetscViewerASCIIPrintf(viewer, " # meaningful timing results on multi-rank runs, use #\n")); 1223 PetscCall(PetscViewerASCIIPrintf(viewer, " # GPU-aware MPI instead. #\n")); 1224 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1225 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n\n\n")); 1226 PetscFunctionReturn(PETSC_SUCCESS); 1227 #else 1228 return PETSC_SUCCESS; 1229 #endif 1230 } 1231 1232 static PetscErrorCode PetscLogViewWarnGpuTime(PetscViewer viewer) 1233 { 1234 #if defined(PETSC_HAVE_DEVICE) && !defined(PETSC_HAVE_KOKKOS_WITHOUT_GPU) 1235 1236 PetscFunctionBegin; 1237 if (!PetscLogGpuTimeFlag || petsc_gflops == 0) PetscFunctionReturn(PETSC_SUCCESS); 1238 PetscCall(PetscViewerASCIIPrintf(viewer, "\n\n")); 1239 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n")); 1240 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1241 PetscCall(PetscViewerASCIIPrintf(viewer, " # WARNING!!! #\n")); 1242 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1243 PetscCall(PetscViewerASCIIPrintf(viewer, " # This code was run with -log_view_gpu_time #\n")); 1244 PetscCall(PetscViewerASCIIPrintf(viewer, " # This provides accurate timing within the GPU kernels #\n")); 1245 PetscCall(PetscViewerASCIIPrintf(viewer, " # but can slow down the entire computation by a #\n")); 1246 PetscCall(PetscViewerASCIIPrintf(viewer, " # measurable amount. For fastest runs we recommend #\n")); 1247 PetscCall(PetscViewerASCIIPrintf(viewer, " # not using this option. #\n")); 1248 PetscCall(PetscViewerASCIIPrintf(viewer, " # #\n")); 1249 PetscCall(PetscViewerASCIIPrintf(viewer, " ##########################################################\n\n\n")); 1250 PetscFunctionReturn(PETSC_SUCCESS); 1251 #else 1252 return PETSC_SUCCESS; 1253 #endif 1254 } 1255 1256 PETSC_INTERN int PetscGlobalArgc; 1257 PETSC_INTERN char **PetscGlobalArgs; 1258 1259 static PetscErrorCode PetscLogHandlerView_Default_Info(PetscLogHandler handler, PetscViewer viewer) 1260 { 1261 PetscLogHandler_Default def = (PetscLogHandler_Default)handler->data; 1262 char arch[128], hostname[128], username[128], pname[PETSC_MAX_PATH_LEN], date[128]; 1263 PetscLogDouble locTotalTime, TotalTime, TotalFlops; 1264 PetscLogDouble numMessages, messageLength, avgMessLen, numReductions; 1265 PetscLogDouble stageTime, flops, flopr, mem, mess, messLen, red; 1266 PetscLogDouble fracTime, fracFlops, fracMessages, fracLength, fracReductions, fracMess, fracMessLen, fracRed; 1267 PetscLogDouble fracStageTime, fracStageFlops, fracStageMess, fracStageMessLen, fracStageRed; 1268 PetscLogDouble min, max, tot, ratio, avg, x, y; 1269 PetscLogDouble minf, maxf, totf, ratf, mint, maxt, tott, ratt, ratC, totm, totml, totr, mal, malmax, emalmax; 1270 #if defined(PETSC_HAVE_DEVICE) 1271 PetscLogEvent KSP_Solve, SNES_Solve, TS_Step, TAO_Solve; /* These need to be fixed to be some events registered with certain objects */ 1272 PetscLogDouble cct, gct, csz, gsz, gmaxt, gflops, gflopr, fracgflops; 1273 #endif 1274 PetscMPIInt minC, maxC; 1275 PetscMPIInt size, rank; 1276 PetscBool *localStageUsed, *stageUsed; 1277 PetscBool *localStageVisible, *stageVisible; 1278 PetscInt numStages, numEvents; 1279 int stage, oclass; 1280 PetscLogEvent event; 1281 char version[256]; 1282 MPI_Comm comm; 1283 #if defined(PETSC_HAVE_DEVICE) && !defined(PETSC_HAVE_KOKKOS_WITHOUT_GPU) 1284 PetscInt64 nas = 0x7FF0000000000002; 1285 #endif 1286 PetscLogGlobalNames global_stages, global_events; 1287 PetscEventPerfInfo zero_info; 1288 PetscLogState state; 1289 1290 PetscFunctionBegin; 1291 PetscCall(PetscLogHandlerGetState(handler, &state)); 1292 PetscCall(PetscFPTrapPush(PETSC_FP_TRAP_OFF)); 1293 PetscCall(PetscObjectGetComm((PetscObject)viewer, &comm)); 1294 PetscCallMPI(MPI_Comm_size(comm, &size)); 1295 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 1296 /* Get the total elapsed time */ 1297 PetscCall(PetscTime(&locTotalTime)); 1298 locTotalTime -= petsc_BaseTime; 1299 1300 PetscCall(PetscViewerASCIIPrintf(viewer, "****************************************************************************************************************************************************************\n")); 1301 PetscCall(PetscViewerASCIIPrintf(viewer, "*** WIDEN YOUR WINDOW TO 160 CHARACTERS. Use 'enscript -r -fCourier9' to print this document ***\n")); 1302 PetscCall(PetscViewerASCIIPrintf(viewer, "****************************************************************************************************************************************************************\n")); 1303 PetscCall(PetscViewerASCIIPrintf(viewer, "\n------------------------------------------------------------------ PETSc Performance Summary: ------------------------------------------------------------------\n\n")); 1304 PetscCall(PetscLogViewWarnSync(viewer)); 1305 PetscCall(PetscLogViewWarnDebugging(viewer)); 1306 PetscCall(PetscLogViewWarnNoGpuAwareMpi(viewer)); 1307 PetscCall(PetscLogViewWarnGpuTime(viewer)); 1308 PetscCall(PetscGetArchType(arch, sizeof(arch))); 1309 PetscCall(PetscGetHostName(hostname, sizeof(hostname))); 1310 PetscCall(PetscGetUserName(username, sizeof(username))); 1311 PetscCall(PetscGetProgramName(pname, sizeof(pname))); 1312 PetscCall(PetscGetDate(date, sizeof(date))); 1313 PetscCall(PetscGetVersion(version, sizeof(version))); 1314 1315 #if defined(PETSC_HAVE_CUPM) 1316 const char *cupm = PetscDefined(HAVE_CUDA) ? "CUDA" : "HIP"; 1317 if (PetscDeviceCUPMRuntimeArch) 1318 PetscCall(PetscViewerASCIIPrintf(viewer, "%s on a %s named %s with %d process%s and %s architecture %d, by %s on %s\n", pname, arch, hostname, size, size > 1 ? "es" : "", cupm, PetscDeviceCUPMRuntimeArch, username, date)); 1319 else 1320 #endif 1321 PetscCall(PetscViewerASCIIPrintf(viewer, "%s on a %s named %s with %d process%s, by %s on %s\n", pname, arch, hostname, size, size > 1 ? "es" : "", username, date)); 1322 1323 #if defined(PETSC_HAVE_OPENMP) 1324 PetscCall(PetscViewerASCIIPrintf(viewer, "Using %" PetscInt_FMT " OpenMP threads\n", PetscNumOMPThreads)); 1325 #endif 1326 PetscCall(PetscViewerASCIIPrintf(viewer, "Using %s\n", version)); 1327 1328 /* Must preserve reduction count before we go on */ 1329 red = petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct; 1330 1331 /* Calculate summary information */ 1332 PetscCall(PetscViewerASCIIPrintf(viewer, "\n Max Max/Min Avg Total\n")); 1333 /* Time */ 1334 PetscCallMPI(MPIU_Allreduce(&locTotalTime, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1335 PetscCallMPI(MPIU_Allreduce(&locTotalTime, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1336 PetscCallMPI(MPIU_Allreduce(&locTotalTime, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1337 avg = tot / ((PetscLogDouble)size); 1338 if (min != 0.0) ratio = max / min; 1339 else ratio = 0.0; 1340 PetscCall(PetscViewerASCIIPrintf(viewer, "Time (sec): %5.3e %7.3f %5.3e\n", max, ratio, avg)); 1341 TotalTime = tot; 1342 /* Objects */ 1343 { 1344 PetscInt num_objects; 1345 1346 PetscCall(PetscLogObjectArrayGetSize(def->petsc_objects, &num_objects, NULL)); 1347 avg = (PetscLogDouble)num_objects; 1348 } 1349 PetscCallMPI(MPIU_Allreduce(&avg, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1350 PetscCallMPI(MPIU_Allreduce(&avg, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1351 PetscCallMPI(MPIU_Allreduce(&avg, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1352 avg = tot / ((PetscLogDouble)size); 1353 if (min != 0.0) ratio = max / min; 1354 else ratio = 0.0; 1355 PetscCall(PetscViewerASCIIPrintf(viewer, "Objects: %5.3e %7.3f %5.3e\n", max, ratio, avg)); 1356 /* Flops */ 1357 PetscCallMPI(MPIU_Allreduce(&petsc_TotalFlops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1358 PetscCallMPI(MPIU_Allreduce(&petsc_TotalFlops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1359 PetscCallMPI(MPIU_Allreduce(&petsc_TotalFlops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1360 avg = tot / ((PetscLogDouble)size); 1361 if (min != 0.0) ratio = max / min; 1362 else ratio = 0.0; 1363 PetscCall(PetscViewerASCIIPrintf(viewer, "Flops: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1364 TotalFlops = tot; 1365 /* Flops/sec -- Must talk to Barry here */ 1366 if (locTotalTime != 0.0) flops = petsc_TotalFlops / locTotalTime; 1367 else flops = 0.0; 1368 PetscCallMPI(MPIU_Allreduce(&flops, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1369 PetscCallMPI(MPIU_Allreduce(&flops, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1370 PetscCallMPI(MPIU_Allreduce(&flops, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1371 avg = tot / ((PetscLogDouble)size); 1372 if (min != 0.0) ratio = max / min; 1373 else ratio = 0.0; 1374 PetscCall(PetscViewerASCIIPrintf(viewer, "Flops/sec: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1375 /* Memory */ 1376 PetscCall(PetscMallocGetMaximumUsage(&mem)); 1377 if (mem > 0.0) { 1378 PetscCallMPI(MPIU_Allreduce(&mem, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1379 PetscCallMPI(MPIU_Allreduce(&mem, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1380 PetscCallMPI(MPIU_Allreduce(&mem, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1381 avg = tot / ((PetscLogDouble)size); 1382 if (min != 0.0) ratio = max / min; 1383 else ratio = 0.0; 1384 PetscCall(PetscViewerASCIIPrintf(viewer, "Memory (bytes): %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1385 } 1386 /* Messages */ 1387 mess = 0.5 * (petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct); 1388 PetscCallMPI(MPIU_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1389 PetscCallMPI(MPIU_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1390 PetscCallMPI(MPIU_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1391 avg = tot / ((PetscLogDouble)size); 1392 if (min != 0.0) ratio = max / min; 1393 else ratio = 0.0; 1394 PetscCall(PetscViewerASCIIPrintf(viewer, "MPI Msg Count: %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1395 numMessages = tot; 1396 /* Message Lengths */ 1397 mess = 0.5 * (petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len); 1398 PetscCallMPI(MPIU_Allreduce(&mess, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1399 PetscCallMPI(MPIU_Allreduce(&mess, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1400 PetscCallMPI(MPIU_Allreduce(&mess, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1401 if (numMessages != 0) avg = tot / numMessages; 1402 else avg = 0.0; 1403 if (min != 0.0) ratio = max / min; 1404 else ratio = 0.0; 1405 PetscCall(PetscViewerASCIIPrintf(viewer, "MPI Msg Len (bytes): %5.3e %7.3f %5.3e %5.3e\n", max, ratio, avg, tot)); 1406 messageLength = tot; 1407 /* Reductions */ 1408 PetscCallMPI(MPIU_Allreduce(&red, &min, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1409 PetscCallMPI(MPIU_Allreduce(&red, &max, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1410 PetscCallMPI(MPIU_Allreduce(&red, &tot, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1411 if (min != 0.0) ratio = max / min; 1412 else ratio = 0.0; 1413 PetscCall(PetscViewerASCIIPrintf(viewer, "MPI Reductions: %5.3e %7.3f\n", max, ratio)); 1414 numReductions = red; /* wrong because uses count from process zero */ 1415 PetscCall(PetscViewerASCIIPrintf(viewer, "\nFlop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)\n")); 1416 PetscCall(PetscViewerASCIIPrintf(viewer, " e.g., VecAXPY() for real vectors of length N --> 2N flops\n")); 1417 PetscCall(PetscViewerASCIIPrintf(viewer, " and VecAXPY() for complex vectors of length N --> 8N flops\n")); 1418 1419 PetscCall(PetscLogRegistryCreateGlobalStageNames(comm, state->registry, &global_stages)); 1420 PetscCall(PetscLogRegistryCreateGlobalEventNames(comm, state->registry, &global_events)); 1421 PetscCall(PetscLogGlobalNamesGetSize(global_stages, NULL, &numStages)); 1422 PetscCall(PetscLogGlobalNamesGetSize(global_events, NULL, &numEvents)); 1423 PetscCall(PetscMemzero(&zero_info, sizeof(zero_info))); 1424 PetscCall(PetscMalloc1(numStages, &localStageUsed)); 1425 PetscCall(PetscMalloc1(numStages, &stageUsed)); 1426 PetscCall(PetscMalloc1(numStages, &localStageVisible)); 1427 PetscCall(PetscMalloc1(numStages, &stageVisible)); 1428 if (numStages > 0) { 1429 for (stage = 0; stage < numStages; stage++) { 1430 PetscInt stage_id; 1431 1432 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_stages, stage, &stage_id)); 1433 if (stage_id >= 0) { 1434 PetscStagePerf *stage_info; 1435 1436 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage, &stage_info)); 1437 localStageUsed[stage] = stage_info->used; 1438 localStageVisible[stage] = stage_info->perfInfo.visible; 1439 } else { 1440 localStageUsed[stage] = PETSC_FALSE; 1441 localStageVisible[stage] = PETSC_TRUE; 1442 } 1443 } 1444 PetscCallMPI(MPIU_Allreduce(localStageUsed, stageUsed, numStages, MPIU_BOOL, MPI_LOR, comm)); 1445 PetscCallMPI(MPIU_Allreduce(localStageVisible, stageVisible, numStages, MPIU_BOOL, MPI_LAND, comm)); 1446 for (stage = 0; stage < numStages; stage++) { 1447 if (stageUsed[stage] && stageVisible[stage]) { 1448 PetscCall(PetscViewerASCIIPrintf(viewer, "\nSummary of Stages: ----- Time ------ ----- Flop ------ --- Messages --- -- Message Lengths -- -- Reductions --\n")); 1449 PetscCall(PetscViewerASCIIPrintf(viewer, " Avg %%Total Avg %%Total Count %%Total Avg %%Total Count %%Total\n")); 1450 break; 1451 } 1452 } 1453 for (stage = 0; stage < numStages; stage++) { 1454 PetscInt stage_id; 1455 PetscEventPerfInfo *stage_info; 1456 const char *stage_name; 1457 1458 if (!(stageUsed[stage] && stageVisible[stage])) continue; 1459 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_stages, stage, &stage_id)); 1460 PetscCall(PetscLogGlobalNamesGlobalGetName(global_stages, stage, &stage_name)); 1461 stage_info = &zero_info; 1462 if (localStageUsed[stage]) { 1463 PetscStagePerf *stage_perf_info; 1464 1465 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage, &stage_perf_info)); 1466 stage_info = &stage_perf_info->perfInfo; 1467 } 1468 PetscCallMPI(MPIU_Allreduce(&stage_info->time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1469 PetscCallMPI(MPIU_Allreduce(&stage_info->flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1470 PetscCallMPI(MPIU_Allreduce(&stage_info->numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1471 PetscCallMPI(MPIU_Allreduce(&stage_info->messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1472 PetscCallMPI(MPIU_Allreduce(&stage_info->numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1473 mess *= 0.5; 1474 messLen *= 0.5; 1475 red /= size; 1476 if (TotalTime != 0.0) fracTime = stageTime / TotalTime; 1477 else fracTime = 0.0; 1478 if (TotalFlops != 0.0) fracFlops = flops / TotalFlops; 1479 else fracFlops = 0.0; 1480 /* Talk to Barry if (stageTime != 0.0) flops = (size*flops)/stageTime; else flops = 0.0; */ 1481 if (numMessages != 0.0) fracMessages = mess / numMessages; 1482 else fracMessages = 0.0; 1483 if (mess != 0.0) avgMessLen = messLen / mess; 1484 else avgMessLen = 0.0; 1485 if (messageLength != 0.0) fracLength = messLen / messageLength; 1486 else fracLength = 0.0; 1487 if (numReductions != 0.0) fracReductions = red / numReductions; 1488 else fracReductions = 0.0; 1489 PetscCall(PetscViewerASCIIPrintf(viewer, "%2d: %15s: %6.4e %5.1f%% %6.4e %5.1f%% %5.3e %5.1f%% %5.3e %5.1f%% %5.3e %5.1f%%\n", stage, stage_name, stageTime / size, 100.0 * fracTime, flops, 100.0 * fracFlops, mess, 100.0 * fracMessages, avgMessLen, 100.0 * fracLength, red, 100.0 * fracReductions)); 1490 } 1491 } 1492 1493 PetscCall(PetscViewerASCIIPrintf(viewer, "\n------------------------------------------------------------------------------------------------------------------------\n")); 1494 PetscCall(PetscViewerASCIIPrintf(viewer, "See the 'Profiling' chapter of the users' manual for details on interpreting output.\n")); 1495 PetscCall(PetscViewerASCIIPrintf(viewer, "Phase summary info:\n")); 1496 PetscCall(PetscViewerASCIIPrintf(viewer, " Count: number of times phase was executed\n")); 1497 PetscCall(PetscViewerASCIIPrintf(viewer, " Time and Flop: Max - maximum over all processors\n")); 1498 PetscCall(PetscViewerASCIIPrintf(viewer, " Ratio - ratio of maximum to minimum over all processors\n")); 1499 PetscCall(PetscViewerASCIIPrintf(viewer, " Mess: number of messages sent\n")); 1500 PetscCall(PetscViewerASCIIPrintf(viewer, " AvgLen: average message length (bytes)\n")); 1501 PetscCall(PetscViewerASCIIPrintf(viewer, " Reduct: number of global reductions\n")); 1502 PetscCall(PetscViewerASCIIPrintf(viewer, " Global: entire computation\n")); 1503 PetscCall(PetscViewerASCIIPrintf(viewer, " Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().\n")); 1504 PetscCall(PetscViewerASCIIPrintf(viewer, " %%T - percent time in this phase %%F - percent flop in this phase\n")); 1505 PetscCall(PetscViewerASCIIPrintf(viewer, " %%M - percent messages in this phase %%L - percent message lengths in this phase\n")); 1506 PetscCall(PetscViewerASCIIPrintf(viewer, " %%R - percent reductions in this phase\n")); 1507 PetscCall(PetscViewerASCIIPrintf(viewer, " Total Mflop/s: 1e-6 * (sum of flop over all processors)/(max time over all processors)\n")); 1508 if (PetscLogMemory) { 1509 PetscCall(PetscViewerASCIIPrintf(viewer, " Memory usage is summed over all MPI processes, it is given in mega-bytes\n")); 1510 PetscCall(PetscViewerASCIIPrintf(viewer, " Malloc Mbytes: Memory allocated and kept during event (sum over all calls to event). May be negative\n")); 1511 PetscCall(PetscViewerASCIIPrintf(viewer, " EMalloc Mbytes: extra memory allocated during event and then freed (maximum over all calls to events). Never negative\n")); 1512 PetscCall(PetscViewerASCIIPrintf(viewer, " MMalloc Mbytes: Increase in high water mark of allocated memory (sum over all calls to event). Never negative\n")); 1513 PetscCall(PetscViewerASCIIPrintf(viewer, " RMI Mbytes: Increase in resident memory (sum over all calls to event)\n")); 1514 } 1515 #if defined(PETSC_HAVE_DEVICE) 1516 PetscCall(PetscViewerASCIIPrintf(viewer, " GPU Mflop/s: 1e-6 * (sum of flop on GPU over all processors)/(max GPU time over all processors)\n")); 1517 PetscCall(PetscViewerASCIIPrintf(viewer, " CpuToGpu Count: total number of CPU to GPU copies per processor\n")); 1518 PetscCall(PetscViewerASCIIPrintf(viewer, " CpuToGpu Size (Mbytes): 1e-6 * (total size of CPU to GPU copies per processor)\n")); 1519 PetscCall(PetscViewerASCIIPrintf(viewer, " GpuToCpu Count: total number of GPU to CPU copies per processor\n")); 1520 PetscCall(PetscViewerASCIIPrintf(viewer, " GpuToCpu Size (Mbytes): 1e-6 * (total size of GPU to CPU copies per processor)\n")); 1521 PetscCall(PetscViewerASCIIPrintf(viewer, " GPU %%F: percent flops on GPU in this event\n")); 1522 #endif 1523 PetscCall(PetscViewerASCIIPrintf(viewer, "------------------------------------------------------------------------------------------------------------------------\n")); 1524 1525 PetscCall(PetscLogViewWarnDebugging(viewer)); 1526 1527 /* Report events */ 1528 PetscCall(PetscViewerASCIIPrintf(viewer, "Event Count Time (sec) Flop --- Global --- --- Stage ---- Total")); 1529 if (PetscLogMemory) PetscCall(PetscViewerASCIIPrintf(viewer, " Malloc EMalloc MMalloc RMI")); 1530 #if defined(PETSC_HAVE_DEVICE) 1531 PetscCall(PetscViewerASCIIPrintf(viewer, " GPU - CpuToGpu - - GpuToCpu - GPU")); 1532 #endif 1533 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1534 PetscCall(PetscViewerASCIIPrintf(viewer, " Max Ratio Max Ratio Max Ratio Mess AvgLen Reduct %%T %%F %%M %%L %%R %%T %%F %%M %%L %%R Mflop/s")); 1535 if (PetscLogMemory) PetscCall(PetscViewerASCIIPrintf(viewer, " Mbytes Mbytes Mbytes Mbytes")); 1536 #if defined(PETSC_HAVE_DEVICE) 1537 PetscCall(PetscViewerASCIIPrintf(viewer, " Mflop/s Count Size Count Size %%F")); 1538 #endif 1539 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1540 PetscCall(PetscViewerASCIIPrintf(viewer, "------------------------------------------------------------------------------------------------------------------------")); 1541 if (PetscLogMemory) PetscCall(PetscViewerASCIIPrintf(viewer, "-----------------------------")); 1542 #if defined(PETSC_HAVE_DEVICE) 1543 PetscCall(PetscViewerASCIIPrintf(viewer, "---------------------------------------")); 1544 #endif 1545 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1546 1547 #if defined(PETSC_HAVE_DEVICE) 1548 /* this indirect way of accessing these values is needed when PETSc is build with multiple libraries since the symbols are not in libpetscsys */ 1549 PetscCall(PetscLogStateGetEventFromName(state, "TaoSolve", &TAO_Solve)); 1550 PetscCall(PetscLogStateGetEventFromName(state, "TSStep", &TS_Step)); 1551 PetscCall(PetscLogStateGetEventFromName(state, "SNESSolve", &SNES_Solve)); 1552 PetscCall(PetscLogStateGetEventFromName(state, "KSPSolve", &KSP_Solve)); 1553 #endif 1554 1555 for (stage = 0; stage < numStages; stage++) { 1556 PetscInt stage_id; 1557 PetscEventPerfInfo *stage_info; 1558 const char *stage_name; 1559 1560 if (!(stageVisible[stage] && stageUsed[stage])) continue; 1561 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_stages, stage, &stage_id)); 1562 PetscCall(PetscLogGlobalNamesGlobalGetName(global_stages, stage, &stage_name)); 1563 PetscCall(PetscViewerASCIIPrintf(viewer, "\n--- Event Stage %d: %s\n\n", stage, stage_name)); 1564 stage_info = &zero_info; 1565 if (localStageUsed[stage]) { 1566 PetscStagePerf *stage_perf_info; 1567 1568 PetscCall(PetscLogHandlerDefaultGetStageInfo(handler, stage_id, &stage_perf_info)); 1569 stage_info = &stage_perf_info->perfInfo; 1570 } 1571 PetscCallMPI(MPIU_Allreduce(&stage_info->time, &stageTime, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1572 PetscCallMPI(MPIU_Allreduce(&stage_info->flops, &flops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1573 PetscCallMPI(MPIU_Allreduce(&stage_info->numMessages, &mess, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1574 PetscCallMPI(MPIU_Allreduce(&stage_info->messageLength, &messLen, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1575 PetscCallMPI(MPIU_Allreduce(&stage_info->numReductions, &red, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1576 mess *= 0.5; 1577 messLen *= 0.5; 1578 red /= size; 1579 1580 for (event = 0; event < numEvents; event++) { 1581 PetscInt event_id; 1582 PetscEventPerfInfo *event_info = &zero_info; 1583 PetscBool is_zero = PETSC_FALSE; 1584 const char *event_name; 1585 1586 PetscCall(PetscLogGlobalNamesGlobalGetLocal(global_events, event, &event_id)); 1587 PetscCall(PetscLogGlobalNamesGlobalGetName(global_events, event, &event_name)); 1588 if (event_id >= 0 && stage_id >= 0) { PetscCall(PetscLogHandlerGetEventPerfInfo_Default(handler, stage_id, event_id, &event_info)); } 1589 PetscCall(PetscMemcmp(event_info, &zero_info, sizeof(zero_info), &is_zero)); 1590 PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, &is_zero, 1, MPIU_BOOL, MPI_LAND, comm)); 1591 if (!is_zero) { 1592 flopr = event_info->flops; 1593 PetscCallMPI(MPIU_Allreduce(&flopr, &minf, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1594 PetscCallMPI(MPIU_Allreduce(&flopr, &maxf, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1595 PetscCallMPI(MPIU_Allreduce(&event_info->flops, &totf, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1596 PetscCallMPI(MPIU_Allreduce(&event_info->time, &mint, 1, MPIU_PETSCLOGDOUBLE, MPI_MIN, comm)); 1597 PetscCallMPI(MPIU_Allreduce(&event_info->time, &maxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1598 PetscCallMPI(MPIU_Allreduce(&event_info->time, &tott, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1599 PetscCallMPI(MPIU_Allreduce(&event_info->numMessages, &totm, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1600 PetscCallMPI(MPIU_Allreduce(&event_info->messageLength, &totml, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1601 PetscCallMPI(MPIU_Allreduce(&event_info->numReductions, &totr, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1602 PetscCallMPI(MPIU_Allreduce(&event_info->count, &minC, 1, MPI_INT, MPI_MIN, comm)); 1603 PetscCallMPI(MPIU_Allreduce(&event_info->count, &maxC, 1, MPI_INT, MPI_MAX, comm)); 1604 if (PetscLogMemory) { 1605 PetscCallMPI(MPIU_Allreduce(&event_info->memIncrease, &mem, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1606 PetscCallMPI(MPIU_Allreduce(&event_info->mallocSpace, &mal, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1607 PetscCallMPI(MPIU_Allreduce(&event_info->mallocIncrease, &malmax, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1608 PetscCallMPI(MPIU_Allreduce(&event_info->mallocIncreaseEvent, &emalmax, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1609 } 1610 #if defined(PETSC_HAVE_DEVICE) 1611 PetscCallMPI(MPIU_Allreduce(&event_info->CpuToGpuCount, &cct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1612 PetscCallMPI(MPIU_Allreduce(&event_info->GpuToCpuCount, &gct, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1613 PetscCallMPI(MPIU_Allreduce(&event_info->CpuToGpuSize, &csz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1614 PetscCallMPI(MPIU_Allreduce(&event_info->GpuToCpuSize, &gsz, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1615 PetscCallMPI(MPIU_Allreduce(&event_info->GpuFlops, &gflops, 1, MPIU_PETSCLOGDOUBLE, MPI_SUM, comm)); 1616 PetscCallMPI(MPIU_Allreduce(&event_info->GpuTime, &gmaxt, 1, MPIU_PETSCLOGDOUBLE, MPI_MAX, comm)); 1617 #endif 1618 if (mint < 0.0) { 1619 PetscCall(PetscViewerASCIIPrintf(viewer, "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, event_name)); 1620 mint = 0; 1621 } 1622 PetscCheck(minf >= 0.0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Minimum flop %g over all processors for %s is negative! Not possible!", minf, event_name); 1623 #if defined(PETSC_HAVE_DEVICE) && !defined(PETSC_HAVE_KOKKOS_WITHOUT_GPU) 1624 /* Put NaN into the time for all events that may not be time accurately since they may happen asynchronously on the GPU */ 1625 if (!PetscLogGpuTimeFlag && petsc_gflops > 0) { 1626 memcpy(&gmaxt, &nas, sizeof(PetscLogDouble)); 1627 if (event_id != SNES_Solve && event_id != KSP_Solve && event_id != TS_Step && event_id != TAO_Solve) { 1628 memcpy(&mint, &nas, sizeof(PetscLogDouble)); 1629 memcpy(&maxt, &nas, sizeof(PetscLogDouble)); 1630 } 1631 } 1632 #endif 1633 totm *= 0.5; 1634 totml *= 0.5; 1635 totr /= size; 1636 1637 if (maxC != 0) { 1638 if (minC != 0) ratC = ((PetscLogDouble)maxC) / minC; 1639 else ratC = 0.0; 1640 if (mint != 0.0) ratt = maxt / mint; 1641 else ratt = 0.0; 1642 if (minf != 0.0) ratf = maxf / minf; 1643 else ratf = 0.0; 1644 if (TotalTime != 0.0) fracTime = tott / TotalTime; 1645 else fracTime = 0.0; 1646 if (TotalFlops != 0.0) fracFlops = totf / TotalFlops; 1647 else fracFlops = 0.0; 1648 if (stageTime != 0.0) fracStageTime = tott / stageTime; 1649 else fracStageTime = 0.0; 1650 if (flops != 0.0) fracStageFlops = totf / flops; 1651 else fracStageFlops = 0.0; 1652 if (numMessages != 0.0) fracMess = totm / numMessages; 1653 else fracMess = 0.0; 1654 if (messageLength != 0.0) fracMessLen = totml / messageLength; 1655 else fracMessLen = 0.0; 1656 if (numReductions != 0.0) fracRed = totr / numReductions; 1657 else fracRed = 0.0; 1658 if (mess != 0.0) fracStageMess = totm / mess; 1659 else fracStageMess = 0.0; 1660 if (messLen != 0.0) fracStageMessLen = totml / messLen; 1661 else fracStageMessLen = 0.0; 1662 if (red != 0.0) fracStageRed = totr / red; 1663 else fracStageRed = 0.0; 1664 if (totm != 0.0) totml /= totm; 1665 else totml = 0.0; 1666 if (maxt != 0.0) flopr = totf / maxt; 1667 else flopr = 0.0; 1668 if (fracStageTime > 1.0 || fracStageFlops > 1.0 || fracStageMess > 1.0 || fracStageMessLen > 1.0 || fracStageRed > 1.0) { 1669 if (PetscIsNanReal(maxt)) 1670 PetscCall(PetscViewerASCIIPrintf(viewer, "%-16s %7d %3.1f n/a n/a %3.2e %3.1f %2.1e %2.1e %2.1e %2.0f %2.0f %2.0f %2.0f %2.0f Multiple stages n/a", event_name, maxC, ratC, maxf, ratf, totm, totml, totr, 100.0 * fracTime, 100.0 * fracFlops, 100.0 * fracMess, 100.0 * fracMessLen, 100.0 * fracRed)); 1671 else 1672 PetscCall(PetscViewerASCIIPrintf(viewer, "%-16s %7d %3.1f %5.4e %3.1f %3.2e %3.1f %2.1e %2.1e %2.1e %2.0f %2.0f %2.0f %2.0f %2.0f Multiple stages %5.0f", event_name, maxC, ratC, maxt, ratt, maxf, ratf, totm, totml, totr, 100.0 * fracTime, 100.0 * fracFlops, 100.0 * fracMess, 100.0 * fracMessLen, 100.0 * fracRed, PetscAbs(flopr) / 1.0e6)); 1673 } else { 1674 if (PetscIsNanReal(maxt)) { // when maxt, ratt, flopr are NaN (i.e., run with GPUs but without -log_view_gpu_time), replace the confusing "nan" with "n/a" 1675 PetscCall(PetscViewerASCIIPrintf(viewer, "%-16s %7d %3.1f n/a n/a %3.2e %3.1f %2.1e %2.1e %2.1e %2.0f %2.0f %2.0f %2.0f %2.0f %3.0f %2.0f %2.0f %2.0f %2.0f n/a", event_name, maxC, ratC, maxf, ratf, totm, totml, totr, 100.0 * fracTime, 100.0 * fracFlops, 100.0 * fracMess, 100.0 * fracMessLen, 100.0 * fracRed, 100.0 * fracStageTime, 100.0 * fracStageFlops, 100.0 * fracStageMess, 100.0 * fracStageMessLen, 100.0 * fracStageRed)); 1676 } else { 1677 PetscCall(PetscViewerASCIIPrintf(viewer, "%-16s %7d %3.1f %5.4e %3.1f %3.2e %3.1f %2.1e %2.1e %2.1e %2.0f %2.0f %2.0f %2.0f %2.0f %3.0f %2.0f %2.0f %2.0f %2.0f %5.0f", event_name, maxC, ratC, maxt, ratt, maxf, ratf, totm, totml, totr, 100.0 * fracTime, 100.0 * fracFlops, 100.0 * fracMess, 100.0 * fracMessLen, 100.0 * fracRed, 100.0 * fracStageTime, 100.0 * fracStageFlops, 100.0 * fracStageMess, 100.0 * fracStageMessLen, 100.0 * fracStageRed, PetscAbs(flopr) / 1.0e6)); 1678 } 1679 } 1680 if (PetscLogMemory) PetscCall(PetscViewerASCIIPrintf(viewer, " %5.0f %5.0f %5.0f %5.0f", mal / 1.0e6, emalmax / 1.0e6, malmax / 1.0e6, mem / 1.0e6)); 1681 #if defined(PETSC_HAVE_DEVICE) 1682 if (totf != 0.0) fracgflops = gflops / totf; 1683 else fracgflops = 0.0; 1684 if (gmaxt != 0.0) gflopr = gflops / gmaxt; 1685 else gflopr = 0.0; 1686 if (PetscIsNanReal(gflopr)) { 1687 PetscCall(PetscViewerASCIIPrintf(viewer, " n/a %4.0f %3.2e %4.0f %3.2e % 2.0f", cct / size, csz / (1.0e6 * size), gct / size, gsz / (1.0e6 * size), 100.0 * fracgflops)); 1688 } else { 1689 PetscCall(PetscViewerASCIIPrintf(viewer, " %5.0f %4.0f %3.2e %4.0f %3.2e % 2.0f", PetscAbs(gflopr) / 1.0e6, cct / size, csz / (1.0e6 * size), gct / size, gsz / (1.0e6 * size), 100.0 * fracgflops)); 1690 } 1691 #endif 1692 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1693 } 1694 } 1695 } 1696 } 1697 1698 /* Memory usage and object creation */ 1699 PetscCall(PetscViewerASCIIPrintf(viewer, "------------------------------------------------------------------------------------------------------------------------")); 1700 if (PetscLogMemory) PetscCall(PetscViewerASCIIPrintf(viewer, "-----------------------------")); 1701 #if defined(PETSC_HAVE_DEVICE) 1702 PetscCall(PetscViewerASCIIPrintf(viewer, "---------------------------------------")); 1703 #endif 1704 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1705 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1706 1707 /* Right now, only stages on the first processor are reported here, meaning only objects associated with 1708 the global communicator, or MPI_COMM_SELF for proc 1. We really should report global stats and then 1709 stats for stages local to processor sets. 1710 */ 1711 /* We should figure out the longest object name here (now 20 characters) */ 1712 PetscCall(PetscViewerASCIIPrintf(viewer, "Object Type Creations Destructions. Reports information only for process 0.\n")); 1713 for (stage = 0; stage < numStages; stage++) { 1714 const char *stage_name; 1715 1716 PetscCall(PetscLogGlobalNamesGlobalGetName(global_stages, stage, &stage_name)); 1717 PetscCall(PetscViewerASCIIPrintf(viewer, "\n--- Event Stage %d: %s\n\n", stage, stage_name)); 1718 if (localStageUsed[stage]) { 1719 PetscInt num_classes; 1720 1721 PetscCall(PetscLogStateGetNumClasses(state, &num_classes)); 1722 for (oclass = 0; oclass < num_classes; oclass++) { 1723 PetscClassPerf *class_perf_info; 1724 1725 PetscCall(PetscLogHandlerDefaultGetClassPerf(handler, stage, oclass, &class_perf_info)); 1726 if (class_perf_info->creations > 0 || class_perf_info->destructions > 0) { 1727 PetscLogClassInfo class_reg_info; 1728 PetscBool flg = PETSC_FALSE; 1729 1730 PetscCall(PetscLogStateClassGetInfo(state, oclass, &class_reg_info)); 1731 if (stage == 0 && oclass == num_classes - 1) { 1732 if (PetscGlobalArgc == 0 && PetscGlobalArgs == NULL) { 1733 PetscCall(PetscStrcmp(class_reg_info.name, "Viewer", &flg)); 1734 PetscCheck(flg && class_perf_info->creations == PetscLogNumViewersCreated && class_perf_info->destructions == PetscLogNumViewersDestroyed, PETSC_COMM_SELF, PETSC_ERR_PLIB, "The last PetscObject type of the main PetscLogStage should be PetscViewer with only %" PetscInt_FMT " log viewers created and %" PetscInt_FMT " destroyed", PetscLogNumViewersCreated, PetscLogNumViewersDestroyed); 1735 continue; 1736 } 1737 } 1738 PetscCall(PetscViewerASCIIPrintf(viewer, "%20s %5d %5d\n", class_reg_info.name, class_perf_info->creations, class_perf_info->destructions)); 1739 } 1740 } 1741 } 1742 } 1743 1744 PetscCall(PetscFree(localStageUsed)); 1745 PetscCall(PetscFree(stageUsed)); 1746 PetscCall(PetscFree(localStageVisible)); 1747 PetscCall(PetscFree(stageVisible)); 1748 PetscCall(PetscLogGlobalNamesDestroy(&global_stages)); 1749 PetscCall(PetscLogGlobalNamesDestroy(&global_events)); 1750 1751 /* Information unrelated to this particular run */ 1752 PetscCall(PetscViewerASCIIPrintf(viewer, "========================================================================================================================\n")); 1753 PetscCall(PetscTime(&y)); 1754 PetscCall(PetscTime(&x)); 1755 PetscCall(PetscTime(&y)); 1756 PetscCall(PetscTime(&y)); 1757 PetscCall(PetscTime(&y)); 1758 PetscCall(PetscTime(&y)); 1759 PetscCall(PetscTime(&y)); 1760 PetscCall(PetscTime(&y)); 1761 PetscCall(PetscTime(&y)); 1762 PetscCall(PetscTime(&y)); 1763 PetscCall(PetscTime(&y)); 1764 PetscCall(PetscTime(&y)); 1765 PetscCall(PetscViewerASCIIPrintf(viewer, "Average time to get PetscTime(): %g\n", (y - x) / 10.0)); 1766 /* MPI information */ 1767 if (size > 1) { 1768 MPI_Status status; 1769 PetscMPIInt tag; 1770 MPI_Comm newcomm; 1771 1772 PetscCallMPI(MPI_Barrier(comm)); 1773 PetscCall(PetscTime(&x)); 1774 PetscCallMPI(MPI_Barrier(comm)); 1775 PetscCallMPI(MPI_Barrier(comm)); 1776 PetscCallMPI(MPI_Barrier(comm)); 1777 PetscCallMPI(MPI_Barrier(comm)); 1778 PetscCallMPI(MPI_Barrier(comm)); 1779 PetscCall(PetscTime(&y)); 1780 PetscCall(PetscViewerASCIIPrintf(viewer, "Average time for MPI_Barrier(): %g\n", (y - x) / 5.0)); 1781 PetscCall(PetscCommDuplicate(comm, &newcomm, &tag)); 1782 PetscCallMPI(MPI_Barrier(comm)); 1783 if (rank) { 1784 PetscCallMPI(MPI_Recv(NULL, 0, MPI_INT, rank - 1, tag, newcomm, &status)); 1785 PetscCallMPI(MPI_Send(NULL, 0, MPI_INT, (rank + 1) % size, tag, newcomm)); 1786 } else { 1787 PetscCall(PetscTime(&x)); 1788 PetscCallMPI(MPI_Send(NULL, 0, MPI_INT, 1, tag, newcomm)); 1789 PetscCallMPI(MPI_Recv(NULL, 0, MPI_INT, size - 1, tag, newcomm, &status)); 1790 PetscCall(PetscTime(&y)); 1791 PetscCall(PetscViewerASCIIPrintf(viewer, "Average time for zero size MPI_Send(): %g\n", (y - x) / size)); 1792 } 1793 PetscCall(PetscCommDestroy(&newcomm)); 1794 } 1795 PetscCall(PetscOptionsView(NULL, viewer)); 1796 1797 /* Machine and compile information */ 1798 if (PetscDefined(USE_FORTRAN_KERNELS)) { 1799 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with FORTRAN kernels\n")); 1800 } else { 1801 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled without FORTRAN kernels\n")); 1802 } 1803 if (PetscDefined(USE_64BIT_INDICES)) { 1804 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with 64-bit PetscInt\n")); 1805 } else if (PetscDefined(USE___FLOAT128)) { 1806 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with 32-bit PetscInt\n")); 1807 } 1808 if (PetscDefined(USE_REAL_SINGLE)) { 1809 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with single precision PetscScalar and PetscReal\n")); 1810 } else if (PetscDefined(USE___FLOAT128)) { 1811 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with 128 bit precision PetscScalar and PetscReal\n")); 1812 } 1813 if (PetscDefined(USE_REAL_MAT_SINGLE)) { 1814 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with single precision matrices\n")); 1815 } else { 1816 PetscCall(PetscViewerASCIIPrintf(viewer, "Compiled with full precision matrices (default)\n")); 1817 } 1818 PetscCall(PetscViewerASCIIPrintf(viewer, "sizeof(short) %d sizeof(int) %d sizeof(long) %d sizeof(void*) %d sizeof(PetscScalar) %d sizeof(PetscInt) %d\n", (int)sizeof(short), (int)sizeof(int), (int)sizeof(long), (int)sizeof(void *), (int)sizeof(PetscScalar), (int)sizeof(PetscInt))); 1819 1820 PetscCall(PetscViewerASCIIPrintf(viewer, "Configure options: %s", petscconfigureoptions)); 1821 PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petscmachineinfo)); 1822 PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petsccompilerinfo)); 1823 PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petsccompilerflagsinfo)); 1824 PetscCall(PetscViewerASCIIPrintf(viewer, "%s", petsclinkerinfo)); 1825 1826 /* Cleanup */ 1827 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1828 PetscCall(PetscLogViewWarnNoGpuAwareMpi(viewer)); 1829 PetscCall(PetscLogViewWarnDebugging(viewer)); 1830 PetscCall(PetscFPTrapPop()); 1831 PetscFunctionReturn(PETSC_SUCCESS); 1832 } 1833 1834 static PetscErrorCode PetscLogHandlerView_Default(PetscLogHandler handler, PetscViewer viewer) 1835 { 1836 PetscViewerFormat format; 1837 1838 PetscFunctionBegin; 1839 PetscCall(PetscViewerGetFormat(viewer, &format)); 1840 if (format == PETSC_VIEWER_DEFAULT || format == PETSC_VIEWER_ASCII_INFO) { 1841 PetscCall(PetscLogHandlerView_Default_Info(handler, viewer)); 1842 } else if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) { 1843 PetscCall(PetscLogHandlerView_Default_Detailed(handler, viewer)); 1844 } else if (format == PETSC_VIEWER_ASCII_CSV) { 1845 PetscCall(PetscLogHandlerView_Default_CSV(handler, viewer)); 1846 } 1847 PetscFunctionReturn(PETSC_SUCCESS); 1848 } 1849 1850 /*MC 1851 PETSCLOGHANDLERDEFAULT - PETSCLOGHANDLERDEFAULT = "default" - A `PetscLogHandler` that collects data for PETSc 1852 default profiling log viewers (`PetscLogView()` and `PetscLogDump()`). A log handler of this type is 1853 created and started (`PetscLogHandlerStart()`) by `PetscLogDefaultBegin()`. 1854 1855 Options Database Keys: 1856 + -log_include_actions - include a growing list of actions (event beginnings and endings, object creations and destructions) in `PetscLogDump()` (`PetscLogActions()`). 1857 - -log_include_objects - include a growing list of object creations and destructions in `PetscLogDump()` (`PetscLogObjects()`). 1858 1859 Level: developer 1860 1861 .seealso: [](ch_profiling), `PetscLogHandler` 1862 M*/ 1863 1864 PETSC_INTERN PetscErrorCode PetscLogHandlerCreate_Default(PetscLogHandler handler) 1865 { 1866 PetscFunctionBegin; 1867 PetscCall(PetscLogHandlerContextCreate_Default((PetscLogHandler_Default *)&handler->data)); 1868 handler->ops->destroy = PetscLogHandlerDestroy_Default; 1869 handler->ops->eventbegin = PetscLogHandlerEventBegin_Default; 1870 handler->ops->eventend = PetscLogHandlerEventEnd_Default; 1871 handler->ops->eventsync = PetscLogHandlerEventSync_Default; 1872 handler->ops->objectcreate = PetscLogHandlerObjectCreate_Default; 1873 handler->ops->objectdestroy = PetscLogHandlerObjectDestroy_Default; 1874 handler->ops->stagepush = PetscLogHandlerStagePush_Default; 1875 handler->ops->stagepop = PetscLogHandlerStagePop_Default; 1876 handler->ops->view = PetscLogHandlerView_Default; 1877 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerGetEventPerfInfo_C", PetscLogHandlerGetEventPerfInfo_Default)); 1878 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerGetStagePerfInfo_C", PetscLogHandlerGetStagePerfInfo_Default)); 1879 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerSetLogActions_C", PetscLogHandlerSetLogActions_Default)); 1880 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerSetLogObjects_C", PetscLogHandlerSetLogObjects_Default)); 1881 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerLogObjectState_C", PetscLogHandlerLogObjectState_Default)); 1882 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerGetNumObjects_C", PetscLogHandlerGetNumObjects_Default)); 1883 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerEventDeactivatePush_C", PetscLogHandlerEventDeactivatePush_Default)); 1884 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerEventDeactivatePop_C", PetscLogHandlerEventDeactivatePop_Default)); 1885 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerEventsPause_C", PetscLogHandlerEventsPause_Default)); 1886 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerEventsResume_C", PetscLogHandlerEventsResume_Default)); 1887 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerDump_C", PetscLogHandlerDump_Default)); 1888 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerStageSetVisible_C", PetscLogHandlerStageSetVisible_Default)); 1889 PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerStageGetVisible_C", PetscLogHandlerStageGetVisible_Default)); 1890 PetscFunctionReturn(PETSC_SUCCESS); 1891 } 1892