Actual source code: ex7.c

  1: const char help[] = "How to create a log handler using the PetscLogHandler interface";

  3: #include <petscsys.h>
  4: #include <petsc/private/hashmapi.h>
  5: #include <petsctime.h>
  6: #include <petscviewer.h>
  7: #include <petsc/private/loghandlerimpl.h>

  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.

 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. */

 19: #define PETSCLOGHANDLEREX7 "ex7"

 21: typedef struct _HandlerCtx *HandlerCtx;

 23: struct _HandlerCtx {
 24:   PetscHMapI running;
 25:   PetscInt   num_objects_created;
 26:   PetscInt   num_objects_destroyed;
 27: };

 29: static PetscErrorCode HandlerCtxCreate(HandlerCtx *ctx_p)
 30: {
 31:   HandlerCtx ctx;

 33:   PetscFunctionBegin;
 34:   PetscCall(PetscNew(ctx_p));
 35:   ctx = *ctx_p;
 36:   PetscCall(PetscHMapICreate(&ctx->running));
 37:   PetscFunctionReturn(PETSC_SUCCESS);
 38: }

 40: static PetscErrorCode HandlerCtxDestroy(HandlerCtx *ctx_p)
 41: {
 42:   HandlerCtx ctx;

 44:   PetscFunctionBegin;
 45:   ctx    = *ctx_p;
 46:   *ctx_p = NULL;
 47:   PetscCall(PetscHMapIDestroy(&ctx->running));
 48:   PetscCall(PetscFree(ctx));
 49:   PetscFunctionReturn(PETSC_SUCCESS);
 50: }

 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:%-33s] " format_string, _rank, _time, PETSC_FUNCTION_NAME, __VA_ARGS__)); \
 59:   } while (0)

 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;

 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");

 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: }

 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;

 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: }

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;

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: }

118: static PetscErrorCode PetscLogHandlerObjectCreate_Ex7(PetscLogHandler h, PetscObject obj)
119: {
120:   HandlerCtx ctx;

122:   PetscFunctionBegin;
123:   ctx = (HandlerCtx)h->data;
124:   ctx->num_objects_created++;
125:   PetscFunctionReturn(PETSC_SUCCESS);
126: }

128: static PetscErrorCode PetscLogHandlerObjectDestroy_Ex7(PetscLogHandler h, PetscObject obj)
129: {
130:   HandlerCtx ctx;

132:   PetscFunctionBegin;
133:   ctx = (HandlerCtx)h->data;
134:   ctx->num_objects_destroyed++;
135:   PetscFunctionReturn(PETSC_SUCCESS);
136: }

138: static PetscErrorCode PetscLogHandlerStagePush_Ex7(PetscLogHandler h, PetscLogStage new_stage)
139: {
140:   PetscLogStage     old_stage;
141:   PetscLogStageInfo new_info;
142:   PetscLogState     state;

144:   PetscFunctionBegin;
145:   PetscCall(PetscLogHandlerGetState(h, &state));
146:   PetscCall(PetscLogStateStageGetInfo(state, new_stage, &new_info));
147:   PetscCall(PetscLogStateGetCurrentStage(state, &old_stage));
148:   if (old_stage >= 0) {
149:     PetscLogStageInfo old_info;
150:     PetscCall(PetscLogStateStageGetInfo(state, old_stage, &old_info));
151:     PrintData("Pushing stage stage \"%s\" (replacing \"%s\")\n", new_info.name, old_info.name);
152:   } else {
153:     PrintData("Pushing initial stage \"%s\"\n", new_info.name);
154:   }
155:   PetscFunctionReturn(PETSC_SUCCESS);
156: }

158: static PetscErrorCode PetscLogHandlerStagePop_Ex7(PetscLogHandler h, PetscLogStage old_stage)
159: {
160:   PetscLogStage     new_stage;
161:   PetscLogStageInfo old_info;
162:   PetscLogState     state;

164:   PetscFunctionBegin;
165:   PetscCall(PetscLogHandlerGetState(h, &state));
166:   PetscCall(PetscLogStateStageGetInfo(state, old_stage, &old_info));
167:   PetscCall(PetscLogStateGetCurrentStage(state, &new_stage));
168:   if (new_stage >= 0) {
169:     PetscLogStageInfo new_info;

171:     PetscCall(PetscLogStateStageGetInfo(state, new_stage, &new_info));
172:     PrintData("Popping stage \"%s\" (back to \"%s\")\n", old_info.name, new_info.name);
173:   } else {
174:     PrintData("Popping initial stage \"%s\"\n", old_info.name);
175:   }
176:   PetscFunctionReturn(PETSC_SUCCESS);
177: }

179: static PetscErrorCode PetscLogHandlerView_Ex7(PetscLogHandler h, PetscViewer viewer)
180: {
181:   PetscBool is_ascii;

183:   PetscFunctionBegin;
184:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &is_ascii));
185:   if (is_ascii) {
186:     HandlerCtx ctx;
187:     PetscInt   num_entries;

189:     ctx = (HandlerCtx)h->data;
190:     PetscCall(PetscHMapIGetSize(ctx->running, &num_entries));
191:     PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " events were seen by the handler\n", num_entries));
192:     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));
193:   }
194:   PetscFunctionReturn(PETSC_SUCCESS);
195: }

197: // An example of overloading one of the methods defined using PetscObjectComposeFunction()
198: static PetscErrorCode PetscLogHandlerLogObjectState_Ex7(PetscLogHandler h, PetscObject obj, const char format[], va_list argp)
199: {
200:   const char *name;

202:   PetscFunctionBegin;
203:   PetscCall(PetscObjectGetName(obj, &name));
204:   PrintData("Logged state for \"%s\": ", name);
205:   PetscCall(PetscVFPrintf(PETSC_STDOUT, format, argp));
206:   PetscCall(PetscPrintf(PETSC_COMM_SELF, "\n"));
207:   PetscFunctionReturn(PETSC_SUCCESS);
208: }

210: static PetscErrorCode PetscLogHandlerDestroy_Ex7(PetscLogHandler h)
211: {
212:   HandlerCtx ctx;

214:   PetscFunctionBegin;
215:   ctx = (HandlerCtx)h->data;
216:   PetscCall(HandlerCtxDestroy(&ctx));
217:   PetscCall(PetscObjectComposeFunction((PetscObject)h, "PetscLogHandlerLogObjectState_C", NULL));
218:   PetscFunctionReturn(PETSC_SUCCESS);
219: }

221: static PetscErrorCode PetscLogHandlerCreate_Ex7(PetscLogHandler handler)
222: {
223:   HandlerCtx ctx;

225:   PetscFunctionBegin;
226:   PetscCall(HandlerCtxCreate(&ctx));
227:   handler->data               = (void *)ctx;
228:   handler->ops->destroy       = PetscLogHandlerDestroy_Ex7;
229:   handler->ops->view          = PetscLogHandlerView_Ex7;
230:   handler->ops->eventbegin    = PetscLogHandlerEventBegin_Ex7;
231:   handler->ops->eventend      = PetscLogHandlerEventEnd_Ex7;
232:   handler->ops->eventsync     = PetscLogHandlerEventSync_Ex7;
233:   handler->ops->objectcreate  = PetscLogHandlerObjectCreate_Ex7;
234:   handler->ops->objectdestroy = PetscLogHandlerObjectDestroy_Ex7;
235:   handler->ops->stagepush     = PetscLogHandlerStagePush_Ex7;
236:   handler->ops->stagepop      = PetscLogHandlerStagePop_Ex7;
237:   PetscCall(PetscObjectComposeFunction((PetscObject)handler, "PetscLogHandlerLogObjectState_C", PetscLogHandlerLogObjectState_Ex7));
238:   PetscFunctionReturn(PETSC_SUCCESS);
239: }

241: int main(int argc, char **argv)
242: {
243:   PetscClassId        user_classid;
244:   PetscLogEvent       event_1, event_2;
245:   PetscLogStage       stage_1;
246:   PetscContainer      user_object;
247:   PetscLogHandler     h;
248:   PetscLogDouble      time;
249:   PetscLogHandlerType type;

251:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
252:   PetscCall(PetscLogHandlerRegister(PETSCLOGHANDLEREX7, PetscLogHandlerCreate_Ex7));

254:   PetscCall(PetscLogHandlerCreate(PETSC_COMM_WORLD, &h));
255:   PetscCall(PetscLogHandlerSetType(h, PETSCLOGHANDLEREX7));
256:   PetscCall(PetscLogHandlerGetType(h, &type));
257:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Log handler type is: %s\n", type));
258:   PetscCall(PetscLogHandlerStart(h));

260:   PetscCall(PetscClassIdRegister("User class", &user_classid));
261:   PetscCall(PetscLogEventRegister("Event 1", user_classid, &event_1));
262:   PetscCall(PetscLogEventRegister("Event 2", user_classid, &event_2));
263:   PetscCall(PetscLogStageRegister("Stage 1", &stage_1));

265:   PetscCall(PetscLogEventBegin(event_1, NULL, NULL, NULL, NULL));
266:   PetscCall(PetscLogStagePush(stage_1));
267:   PetscCall(PetscLogEventBegin(event_2, NULL, NULL, NULL, NULL));
268:   PetscCall(PetscLogEventSync(event_1, PETSC_COMM_WORLD));
269:   PetscCall(PetscLogEventBegin(event_1, NULL, NULL, NULL, NULL));
270:   PetscCall(PetscTime(&time));
271:   PetscCall(PetscContainerCreate(PETSC_COMM_SELF, &user_object));
272:   PetscCall(PetscObjectSetName((PetscObject)user_object, "User Container"));
273:   PetscCall(PetscLogHandlerLogObjectState(h, (PetscObject)user_object, "Created at %e", time));
274:   PetscCall(PetscContainerDestroy(&user_object));
275:   PetscCall(PetscLogEventEnd(event_1, NULL, NULL, NULL, NULL));
276:   PetscCall(PetscLogEventEnd(event_2, NULL, NULL, NULL, NULL));
277:   PetscCall(PetscLogStagePop());
278:   PetscCall(PetscLogEventEnd(event_1, NULL, NULL, NULL, NULL));

280:   PetscCall(PetscLogHandlerStop(h));
281:   PetscCall(PetscLogHandlerView(h, PETSC_VIEWER_STDOUT_WORLD));
282:   PetscCall(PetscLogHandlerDestroy(&h));
283:   PetscCall(PetscFinalize());
284:   return 0;
285: }

287: /*TEST

289:   test:
290:     requires: defined(PETSC_USE_LOG)
291:     suffix: 0
292:     filter: sed -E "s/:[^:]+:/:time_removed:/g"

294: TEST*/