1 const char help[] = "How to create a log handler using the PetscLogHandler interface"; 2 3 #include <petscsys.h> 4 #include <petsc/private/hashmapi.h> // use PetscHMapI: a PetscInt -> PetscInt hashmap 5 #include <petsctime.h> // use PetscTimeSubtract() and PetscTimeAdd() 6 #include <petscviewer.h> 7 #include <petsc/private/loghandlerimpl.h> // use the struct _p_PetscLogHandler behind PetscLogHandler 8 9 /* Log handlers that use the PetscLogHandler interface get their information 10 from the PetscLogState available to each handler and the user-defined 11 context pointer. Compare this example to src/sys/tutorials/ex6.c. 12 13 A logging event can be started multiple times before it stops: for example, 14 a linear solve may involve a subsolver, so PetscLogEventBegin() can be 15 called for the event KSP_Solve multiple times before a call to 16 PetscLogEventEnd(). The user defined handler in this example shows how many 17 times an event is running. */ 18 19 #define PETSC_LOG_HANDLER_EX7 "ex7" 20 21 typedef struct _HandlerCtx *HandlerCtx; 22 23 struct _HandlerCtx { 24 PetscHMapI running; 25 PetscInt num_objects_created; 26 PetscInt num_objects_destroyed; 27 }; 28 29 static PetscErrorCode HandlerCtxCreate(HandlerCtx *ctx_p) 30 { 31 HandlerCtx ctx; 32 33 PetscFunctionBegin; 34 PetscCall(PetscNew(ctx_p)); 35 ctx = *ctx_p; 36 PetscCall(PetscHMapICreate(&ctx->running)); 37 PetscFunctionReturn(PETSC_SUCCESS); 38 } 39 40 static PetscErrorCode HandlerCtxDestroy(HandlerCtx *ctx_p) 41 { 42 HandlerCtx ctx; 43 44 PetscFunctionBegin; 45 ctx = *ctx_p; 46 *ctx_p = NULL; 47 PetscCall(PetscHMapIDestroy(&ctx->running)); 48 PetscCall(PetscFree(ctx)); 49 PetscFunctionReturn(PETSC_SUCCESS); 50 } 51 52 #define PrintData(format_string, ...) \ 53 do { \ 54 PetscMPIInt _rank; \ 55 PetscLogDouble _time; \ 56 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &_rank)); \ 57 PetscCall(PetscTime(&_time)); \ 58 PetscCall(PetscPrintf(PETSC_COMM_SELF, "[%d:%-7g:%-29s] " format_string, _rank, _time, PETSC_FUNCTION_NAME, __VA_ARGS__)); \ 59 } while (0) 60 61 static PetscErrorCode PetscLogHandlerEventBegin_Ex7(PetscLogHandler h, PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4) 62 { 63 HandlerCtx ctx; 64 PetscInt count; 65 PetscLogState state; 66 PetscLogEventInfo event_info; 67 PetscBool is_active; 68 69 PetscFunctionBegin; 70 // This callback will only be invoked if the event is active 71 PetscCall(PetscLogHandlerGetState(h, &state)); 72 PetscCall(PetscLogStateEventGetActive(state, PETSC_DEFAULT, e, &is_active)); 73 PetscAssert(is_active, PETSC_COMM_SELF, PETSC_ERR_PLIB, "Event handler called for inactive event"); 74 75 ctx = (HandlerCtx)h->data; 76 PetscCall(PetscHMapIGetWithDefault(ctx->running, (PetscInt)e, 0, &count)); 77 count += 1; 78 PetscCall(PetscLogStateEventGetInfo(state, e, &event_info)); 79 PrintData("Event \"%s\" started: now running %" PetscInt_FMT " times\n", event_info.name, count); 80 PetscCall(PetscHMapISet(ctx->running, (PetscInt)e, count)); 81 PetscFunctionReturn(PETSC_SUCCESS); 82 } 83 84 static PetscErrorCode PetscLogHandlerEventEnd_Ex7(PetscLogHandler h, PetscLogEvent e, PetscObject o1, PetscObject o2, PetscObject o3, PetscObject o4) 85 { 86 HandlerCtx ctx; 87 PetscInt count; 88 PetscLogState state; 89 PetscLogEventInfo event_info; 90 91 PetscFunctionBegin; 92 ctx = (HandlerCtx)h->data; 93 PetscCall(PetscLogHandlerGetState(h, &state)); 94 PetscCall(PetscHMapIGetWithDefault(ctx->running, (PetscInt)e, 0, &count)); 95 count -= 1; 96 PetscCall(PetscLogStateEventGetInfo(state, e, &event_info)); 97 PrintData("Event \"%s\" stopped: now running %" PetscInt_FMT " times\n", event_info.name, count); 98 PetscCall(PetscHMapISet(ctx->running, (PetscInt)e, count)); 99 PetscFunctionReturn(PETSC_SUCCESS); 100 } 101 102 static PetscErrorCode PetscLogHandlerEventSync_Ex7(PetscLogHandler h, PetscLogEvent e, MPI_Comm comm) 103 { 104 PetscLogState state; 105 PetscLogEventInfo event_info; 106 PetscLogDouble time = 0.0; 107 108 PetscFunctionBegin; 109 PetscCall(PetscTimeSubtract(&time)); 110 PetscCallMPI(MPI_Barrier(comm)); 111 PetscCall(PetscTimeAdd(&time)); 112 PetscCall(PetscLogHandlerGetState(h, &state)); 113 PetscCall(PetscLogStateEventGetInfo(state, e, &event_info)); 114 PrintData("Event \"%s\" synced: took %g seconds\n", event_info.name, (double)time); 115 PetscFunctionReturn(PETSC_SUCCESS); 116 } 117 118 static PetscErrorCode PetscLogHandlerObjectCreate_Ex7(PetscLogHandler h, PetscObject obj) 119 { 120 HandlerCtx ctx; 121 122 PetscFunctionBegin; 123 ctx = (HandlerCtx)h->data; 124 ctx->num_objects_created++; 125 PetscFunctionReturn(PETSC_SUCCESS); 126 } 127 128 static PetscErrorCode PetscLogHandlerObjectDestroy_Ex7(PetscLogHandler h, PetscObject obj) 129 { 130 HandlerCtx ctx; 131 132 PetscFunctionBegin; 133 ctx = (HandlerCtx)h->data; 134 ctx->num_objects_destroyed++; 135 PetscFunctionReturn(PETSC_SUCCESS); 136 } 137 138 static PetscErrorCode PetscLogHandlerDestroy_Ex7(PetscLogHandler h) 139 { 140 HandlerCtx ctx; 141 142 PetscFunctionBegin; 143 ctx = (HandlerCtx)h->data; 144 PetscCall(HandlerCtxDestroy(&ctx)); 145 PetscFunctionReturn(PETSC_SUCCESS); 146 } 147 148 static PetscErrorCode PetscLogHandlerStagePush_Ex7(PetscLogHandler h, PetscLogStage new_stage) 149 { 150 PetscLogStage old_stage; 151 PetscLogStageInfo new_info; 152 PetscLogState state; 153 154 PetscFunctionBegin; 155 PetscCall(PetscLogHandlerGetState(h, &state)); 156 PetscCall(PetscLogStateStageGetInfo(state, new_stage, &new_info)); 157 PetscCall(PetscLogStateGetCurrentStage(state, &old_stage)); 158 if (old_stage >= 0) { 159 PetscLogStageInfo old_info; 160 PetscCall(PetscLogStateStageGetInfo(state, old_stage, &old_info)); 161 PrintData("Pushing stage stage \"%s\" (replacing \"%s\")\n", new_info.name, old_info.name); 162 } else { 163 PrintData("Pushing initial stage \"%s\"\n", new_info.name); 164 } 165 PetscFunctionReturn(PETSC_SUCCESS); 166 } 167 168 static PetscErrorCode PetscLogHandlerStagePop_Ex7(PetscLogHandler h, PetscLogStage old_stage) 169 { 170 PetscLogStage new_stage; 171 PetscLogStageInfo old_info; 172 PetscLogState state; 173 174 PetscFunctionBegin; 175 PetscCall(PetscLogHandlerGetState(h, &state)); 176 PetscCall(PetscLogStateStageGetInfo(state, old_stage, &old_info)); 177 PetscCall(PetscLogStateGetCurrentStage(state, &new_stage)); 178 if (new_stage >= 0) { 179 PetscLogStageInfo new_info; 180 181 PetscCall(PetscLogStateStageGetInfo(state, new_stage, &new_info)); 182 PrintData("Popping stage \"%s\" (back to \"%s\")\n", old_info.name, new_info.name); 183 } else { 184 PrintData("Popping initial stage \"%s\"\n", old_info.name); 185 } 186 PetscFunctionReturn(PETSC_SUCCESS); 187 } 188 189 static PetscErrorCode PetscLogHandlerView_Ex7(PetscLogHandler h, PetscViewer viewer) 190 { 191 PetscBool is_ascii; 192 193 PetscFunctionBegin; 194 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &is_ascii)); 195 if (is_ascii) { 196 HandlerCtx ctx; 197 PetscInt num_entries; 198 199 ctx = (HandlerCtx)h->data; 200 PetscCall(PetscHMapIGetSize(ctx->running, &num_entries)); 201 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " events were seen by the handler\n", num_entries)); 202 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " object(s) were created and %" PetscInt_FMT " object(s) were destroyed\n", ctx->num_objects_created, ctx->num_objects_created)); 203 } 204 PetscFunctionReturn(PETSC_SUCCESS); 205 } 206 207 static PetscErrorCode PetscLogHandlerCreate_Ex7(PetscLogHandler handler) 208 { 209 HandlerCtx ctx; 210 211 PetscFunctionBegin; 212 PetscCall(HandlerCtxCreate(&ctx)); 213 handler->data = (void *)ctx; 214 handler->ops->destroy = PetscLogHandlerDestroy_Ex7; 215 handler->ops->view = PetscLogHandlerView_Ex7; 216 handler->ops->eventbegin = PetscLogHandlerEventBegin_Ex7; 217 handler->ops->eventend = PetscLogHandlerEventEnd_Ex7; 218 handler->ops->eventsync = PetscLogHandlerEventSync_Ex7; 219 handler->ops->objectcreate = PetscLogHandlerObjectCreate_Ex7; 220 handler->ops->objectdestroy = PetscLogHandlerObjectDestroy_Ex7; 221 handler->ops->stagepush = PetscLogHandlerStagePush_Ex7; 222 handler->ops->stagepop = PetscLogHandlerStagePop_Ex7; 223 PetscFunctionReturn(PETSC_SUCCESS); 224 } 225 226 int main(int argc, char **argv) 227 { 228 PetscClassId user_classid; 229 PetscLogEvent event_1, event_2; 230 PetscLogStage stage_1; 231 PetscContainer user_object; 232 PetscLogHandler h; 233 PetscLogHandlerType type; 234 235 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 236 PetscCall(PetscLogHandlerRegister(PETSC_LOG_HANDLER_EX7, PetscLogHandlerCreate_Ex7)); 237 238 PetscCall(PetscLogHandlerCreate(PETSC_COMM_WORLD, &h)); 239 PetscCall(PetscLogHandlerSetType(h, PETSC_LOG_HANDLER_EX7)); 240 PetscCall(PetscLogHandlerGetType(h, &type)); 241 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Log handler type is: %s\n", type)); 242 PetscCall(PetscLogHandlerStart(h)); 243 244 PetscCall(PetscClassIdRegister("User class", &user_classid)); 245 PetscCall(PetscLogEventRegister("Event 1", user_classid, &event_1)); 246 PetscCall(PetscLogEventRegister("Event 2", user_classid, &event_2)); 247 PetscCall(PetscLogStageRegister("Stage 1", &stage_1)); 248 249 PetscCall(PetscLogEventBegin(event_1, NULL, NULL, NULL, NULL)); 250 PetscCall(PetscLogStagePush(stage_1)); 251 PetscCall(PetscLogEventBegin(event_2, NULL, NULL, NULL, NULL)); 252 PetscCall(PetscLogEventSync(event_1, PETSC_COMM_WORLD)); 253 PetscCall(PetscLogEventBegin(event_1, NULL, NULL, NULL, NULL)); 254 PetscCall(PetscContainerCreate(PETSC_COMM_SELF, &user_object)); 255 PetscCall(PetscContainerDestroy(&user_object)); 256 PetscCall(PetscLogEventEnd(event_1, NULL, NULL, NULL, NULL)); 257 PetscCall(PetscLogEventEnd(event_2, NULL, NULL, NULL, NULL)); 258 PetscCall(PetscLogStagePop()); 259 PetscCall(PetscLogEventEnd(event_1, NULL, NULL, NULL, NULL)); 260 261 PetscCall(PetscLogHandlerStop(h)); 262 PetscCall(PetscLogHandlerView(h, PETSC_VIEWER_STDOUT_WORLD)); 263 PetscCall(PetscLogHandlerDestroy(&h)); 264 PetscCall(PetscFinalize()); 265 return 0; 266 } 267 268 /*TEST 269 270 test: 271 requires: defined(PETSC_USE_LOG) 272 suffix: 0 273 filter: sed -E "s/:[^:]+:/:time_removed:/g" 274 275 TEST*/ 276