Actual source code: eventlog.c
2: /*
3: This defines part of the private API for logging performance information. It is intended to be used only by the
4: PETSc PetscLog...() interface and not elsewhere, nor by users. Hence the prototypes for these functions are NOT
5: in the public PETSc include files.
7: */
8: #include <petsc/private/logimpl.h>
10: PetscBool PetscLogSyncOn = PETSC_FALSE;
11: PetscBool PetscLogMemory = PETSC_FALSE;
12: #if defined(PETSC_HAVE_DEVICE)
13: PetscBool PetscLogGpuTraffic = PETSC_FALSE;
14: #endif
16: /*----------------------------------------------- Creation Functions -------------------------------------------------*/
17: /* Note: these functions do not have prototypes in a public directory, so they are considered "internal" and not exported. */
19: /*@C
20: PetscEventRegLogCreate - This creates a PetscEventRegLog object.
22: Not collective
24: Input Parameter:
25: . eventLog - The PetscEventRegLog
27: Level: developer
29: .seealso: PetscEventRegLogDestroy(), PetscStageLogCreate()
30: @*/
31: PetscErrorCode PetscEventRegLogCreate(PetscEventRegLog *eventLog)
32: {
33: PetscEventRegLog l;
35: PetscNew(&l);
36: l->numEvents = 0;
37: l->maxEvents = 100;
38: PetscMalloc1(l->maxEvents,&l->eventInfo);
39: *eventLog = l;
40: return 0;
41: }
43: /*@C
44: PetscEventRegLogDestroy - This destroys a PetscEventRegLog object.
46: Not collective
48: Input Parameter:
49: . eventLog - The PetscEventRegLog
51: Level: developer
53: .seealso: PetscEventRegLogCreate()
54: @*/
55: PetscErrorCode PetscEventRegLogDestroy(PetscEventRegLog eventLog)
56: {
57: int e;
59: for (e = 0; e < eventLog->numEvents; e++) {
60: PetscFree(eventLog->eventInfo[e].name);
61: }
62: PetscFree(eventLog->eventInfo);
63: PetscFree(eventLog);
64: return 0;
65: }
67: /*@C
68: PetscEventPerfLogCreate - This creates a PetscEventPerfLog object.
70: Not collective
72: Input Parameter:
73: . eventLog - The PetscEventPerfLog
75: Level: developer
77: .seealso: PetscEventPerfLogDestroy(), PetscStageLogCreate()
78: @*/
79: PetscErrorCode PetscEventPerfLogCreate(PetscEventPerfLog *eventLog)
80: {
81: PetscEventPerfLog l;
83: PetscNew(&l);
84: l->numEvents = 0;
85: l->maxEvents = 100;
86: PetscCalloc1(l->maxEvents,&l->eventInfo);
87: *eventLog = l;
88: return 0;
89: }
91: /*@C
92: PetscEventPerfLogDestroy - This destroys a PetscEventPerfLog object.
94: Not collective
96: Input Parameter:
97: . eventLog - The PetscEventPerfLog
99: Level: developer
101: .seealso: PetscEventPerfLogCreate()
102: @*/
103: PetscErrorCode PetscEventPerfLogDestroy(PetscEventPerfLog eventLog)
104: {
105: PetscFree(eventLog->eventInfo);
106: PetscFree(eventLog);
107: return 0;
108: }
110: /*------------------------------------------------ General Functions -------------------------------------------------*/
111: /*@C
112: PetscEventPerfInfoClear - This clears a PetscEventPerfInfo object.
114: Not collective
116: Input Parameter:
117: . eventInfo - The PetscEventPerfInfo
119: Level: developer
121: .seealso: PetscEventPerfLogCreate()
122: @*/
123: PetscErrorCode PetscEventPerfInfoClear(PetscEventPerfInfo *eventInfo)
124: {
125: eventInfo->id = -1;
126: eventInfo->active = PETSC_TRUE;
127: eventInfo->visible = PETSC_TRUE;
128: eventInfo->depth = 0;
129: eventInfo->count = 0;
130: eventInfo->flops = 0.0;
131: eventInfo->flops2 = 0.0;
132: eventInfo->flopsTmp = 0.0;
133: eventInfo->time = 0.0;
134: eventInfo->time2 = 0.0;
135: eventInfo->timeTmp = 0.0;
136: eventInfo->syncTime = 0.0;
137: eventInfo->dof[0] = -1.0;
138: eventInfo->dof[1] = -1.0;
139: eventInfo->dof[2] = -1.0;
140: eventInfo->dof[3] = -1.0;
141: eventInfo->dof[4] = -1.0;
142: eventInfo->dof[5] = -1.0;
143: eventInfo->dof[6] = -1.0;
144: eventInfo->dof[7] = -1.0;
145: eventInfo->errors[0] = -1.0;
146: eventInfo->errors[1] = -1.0;
147: eventInfo->errors[2] = -1.0;
148: eventInfo->errors[3] = -1.0;
149: eventInfo->errors[4] = -1.0;
150: eventInfo->errors[5] = -1.0;
151: eventInfo->errors[6] = -1.0;
152: eventInfo->errors[7] = -1.0;
153: eventInfo->numMessages = 0.0;
154: eventInfo->messageLength = 0.0;
155: eventInfo->numReductions = 0.0;
156: #if defined(PETSC_HAVE_DEVICE)
157: eventInfo->CpuToGpuCount = 0.0;
158: eventInfo->GpuToCpuCount = 0.0;
159: eventInfo->CpuToGpuSize = 0.0;
160: eventInfo->GpuToCpuSize = 0.0;
161: eventInfo->GpuFlops = 0.0;
162: eventInfo->GpuTime = 0.0;
163: #endif
164: return 0;
165: }
167: /*@C
168: PetscEventPerfInfoCopy - Copy the activity and visibility data in eventInfo to outInfo
170: Not collective
172: Input Parameter:
173: . eventInfo - The input PetscEventPerfInfo
175: Output Parameter:
176: . outInfo - The output PetscEventPerfInfo
178: Level: developer
180: .seealso: PetscEventPerfInfoClear()
181: @*/
182: PetscErrorCode PetscEventPerfInfoCopy(PetscEventPerfInfo *eventInfo,PetscEventPerfInfo *outInfo)
183: {
184: outInfo->id = eventInfo->id;
185: outInfo->active = eventInfo->active;
186: outInfo->visible = eventInfo->visible;
187: return 0;
188: }
190: /*@C
191: PetscEventPerfLogEnsureSize - This ensures that a PetscEventPerfLog is at least of a certain size.
193: Not collective
195: Input Parameters:
196: + eventLog - The PetscEventPerfLog
197: - size - The size
199: Level: developer
201: .seealso: PetscEventPerfLogCreate()
202: @*/
203: PetscErrorCode PetscEventPerfLogEnsureSize(PetscEventPerfLog eventLog,int size)
204: {
205: PetscEventPerfInfo *eventInfo;
207: while (size > eventLog->maxEvents) {
208: PetscCalloc1(eventLog->maxEvents*2,&eventInfo);
209: PetscArraycpy(eventInfo,eventLog->eventInfo,eventLog->maxEvents);
210: PetscFree(eventLog->eventInfo);
211: eventLog->eventInfo = eventInfo;
212: eventLog->maxEvents *= 2;
213: }
214: while (eventLog->numEvents < size) {
215: PetscEventPerfInfoClear(&eventLog->eventInfo[eventLog->numEvents++]);
216: }
217: return 0;
218: }
220: #if defined(PETSC_HAVE_MPE)
221: #include <mpe.h>
222: PETSC_INTERN PetscErrorCode PetscLogMPEGetRGBColor(const char*[]);
223: PetscErrorCode PetscLogEventBeginMPE(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
224: {
225: MPE_Log_event(petsc_stageLog->eventLog->eventInfo[event].mpe_id_begin,0,NULL);
226: return 0;
227: }
229: PetscErrorCode PetscLogEventEndMPE(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
230: {
231: MPE_Log_event(petsc_stageLog->eventLog->eventInfo[event].mpe_id_end,0,NULL);
232: return 0;
233: }
234: #endif
236: /*--------------------------------------------- Registration Functions ----------------------------------------------*/
237: /*@C
238: PetscEventRegLogRegister - Registers an event for logging operations in an application code.
240: Not Collective
242: Input Parameters:
243: + eventLog - The PetscEventLog
244: . ename - The name associated with the event
245: - classid - The classid associated to the class for this event
247: Output Parameter:
248: . event - The event
250: Example of Usage:
251: .vb
252: int USER_EVENT;
253: PetscLogDouble user_event_flops;
254: PetscLogEventRegister("User event name",0,&USER_EVENT);
255: PetscLogEventBegin(USER_EVENT,0,0,0,0);
256: [code segment to monitor]
257: PetscLogFlops(user_event_flops);
258: PetscLogEventEnd(USER_EVENT,0,0,0,0);
259: .ve
261: Notes:
263: PETSc can gather data for use with the utilities Jumpshot
264: (part of the MPICH distribution). If PETSc has been compiled
265: with flag -DPETSC_HAVE_MPE (MPE is an additional utility within
266: MPICH), the user can employ another command line option, -log_mpe,
267: to create a logfile, "mpe.log", which can be visualized
268: Jumpshot.
270: Level: developer
272: .seealso: PetscLogEventBegin(), PetscLogEventEnd(), PetscLogFlops(),
273: PetscEventLogActivate(), PetscEventLogDeactivate()
274: @*/
275: PetscErrorCode PetscEventRegLogRegister(PetscEventRegLog eventLog,const char ename[],PetscClassId classid,PetscLogEvent *event)
276: {
277: PetscEventRegInfo *eventInfo;
278: char *str;
279: int e;
283: /* Should check classid I think */
284: e = eventLog->numEvents++;
285: if (eventLog->numEvents > eventLog->maxEvents) {
286: PetscCalloc1(eventLog->maxEvents*2,&eventInfo);
287: PetscArraycpy(eventInfo,eventLog->eventInfo,eventLog->maxEvents);
288: PetscFree(eventLog->eventInfo);
289: eventLog->eventInfo = eventInfo;
290: eventLog->maxEvents *= 2;
291: }
292: PetscStrallocpy(ename,&str);
294: eventLog->eventInfo[e].name = str;
295: eventLog->eventInfo[e].classid = classid;
296: eventLog->eventInfo[e].collective = PETSC_TRUE;
297: #if defined(PETSC_HAVE_MPE)
298: if (PetscLogPLB == PetscLogEventBeginMPE) {
299: const char *color;
300: PetscMPIInt rank;
301: int beginID,endID;
303: beginID = MPE_Log_get_event_number();
304: endID = MPE_Log_get_event_number();
306: eventLog->eventInfo[e].mpe_id_begin = beginID;
307: eventLog->eventInfo[e].mpe_id_end = endID;
309: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
310: if (rank == 0) {
311: PetscLogMPEGetRGBColor(&color);
312: MPE_Describe_state(beginID,endID,str,(char*)color);
313: }
314: }
315: #endif
316: *event = e;
317: return 0;
318: }
320: /*---------------------------------------------- Activation Functions -----------------------------------------------*/
321: /*@C
322: PetscEventPerfLogActivate - Indicates that a particular event should be logged.
324: Not Collective
326: Input Parameters:
327: + eventLog - The PetscEventPerfLog
328: - event - The event
330: Usage:
331: .vb
332: PetscEventPerfLogDeactivate(log, VEC_SetValues);
333: [code where you do not want to log VecSetValues()]
334: PetscEventPerfLogActivate(log, VEC_SetValues);
335: [code where you do want to log VecSetValues()]
336: .ve
338: Note:
339: The event may be either a pre-defined PETSc event (found in
340: include/petsclog.h) or an event number obtained with PetscEventRegLogRegister().
342: Level: developer
344: .seealso: PetscEventPerfLogDeactivate(), PetscEventPerfLogDeactivatePop(), PetscEventPerfLogDeactivatePush()
345: @*/
346: PetscErrorCode PetscEventPerfLogActivate(PetscEventPerfLog eventLog,PetscLogEvent event)
347: {
348: eventLog->eventInfo[event].active = PETSC_TRUE;
349: return 0;
350: }
352: /*@C
353: PetscEventPerfLogDeactivate - Indicates that a particular event should not be logged.
355: Not Collective
357: Input Parameters:
358: + eventLog - The PetscEventPerfLog
359: - event - The event
361: Usage:
362: .vb
363: PetscEventPerfLogDeactivate(log, VEC_SetValues);
364: [code where you do not want to log VecSetValues()]
365: PetscEventPerfLogActivate(log, VEC_SetValues);
366: [code where you do want to log VecSetValues()]
367: .ve
369: Note:
370: The event may be either a pre-defined PETSc event (found in
371: include/petsclog.h) or an event number obtained with PetscEventRegLogRegister().
373: Level: developer
375: .seealso: PetscEventPerfLogActivate(), PetscEventPerfLogDeactivatePop(), PetscEventPerfLogDeactivatePush()
376: @*/
377: PetscErrorCode PetscEventPerfLogDeactivate(PetscEventPerfLog eventLog,PetscLogEvent event)
378: {
379: eventLog->eventInfo[event].active = PETSC_FALSE;
380: return 0;
381: }
383: /*@C
384: PetscEventPerfLogDeactivatePush - Indicates that a particular event should not be logged.
386: Not Collective
388: Input Parameters:
389: + eventLog - The PetscEventPerfLog
390: - event - The event
392: Usage:
393: .vb
394: PetscEventPerfLogDeactivatePush(log, VEC_SetValues);
395: [code where you do not want to log VecSetValues()]
396: PetscEventPerfLogDeactivatePop(log, VEC_SetValues);
397: [code where you do want to log VecSetValues()]
398: .ve
400: Note:
401: The event may be either a pre-defined PETSc event (found in
402: include/petsclog.h) or an event number obtained with PetscEventRegLogRegister().
404: Level: developer
406: .seealso: PetscEventPerfLogDeactivate(), PetscEventPerfLogActivate(), PetscEventPerfLogDeactivatePop()
407: @*/
408: PetscErrorCode PetscEventPerfLogDeactivatePush(PetscEventPerfLog eventLog,PetscLogEvent event)
409: {
410: eventLog->eventInfo[event].depth++;
411: return 0;
412: }
414: /*@C
415: PetscEventPerfLogDeactivatePop - Indicates that a particular event should be logged.
417: Not Collective
419: Input Parameters:
420: + eventLog - The PetscEventPerfLog
421: - event - The event
423: Usage:
424: .vb
425: PetscEventPerfLogDeactivatePush(log, VEC_SetValues);
426: [code where you do not want to log VecSetValues()]
427: PetscEventPerfLogDeactivatePop(log, VEC_SetValues);
428: [code where you do want to log VecSetValues()]
429: .ve
431: Note:
432: The event may be either a pre-defined PETSc event (found in
433: include/petsclog.h) or an event number obtained with PetscEventRegLogRegister().
435: Level: developer
437: .seealso: PetscEventPerfLogDeactivate(), PetscEventPerfLogActivate(), PetscEventPerfLogDeactivatePush()
438: @*/
439: PetscErrorCode PetscEventPerfLogDeactivatePop(PetscEventPerfLog eventLog,PetscLogEvent event)
440: {
441: eventLog->eventInfo[event].depth--;
442: return 0;
443: }
445: /*@C
446: PetscEventPerfLogActivateClass - Activates event logging for a PETSc object class.
448: Not Collective
450: Input Parameters:
451: + eventLog - The PetscEventPerfLog
452: . eventRegLog - The PetscEventRegLog
453: - classid - The class id, for example MAT_CLASSID, SNES_CLASSID,
455: Level: developer
457: .seealso: PetscEventPerfLogDeactivateClass(), PetscEventPerfLogActivate(), PetscEventPerfLogDeactivate()
458: @*/
459: PetscErrorCode PetscEventPerfLogActivateClass(PetscEventPerfLog eventLog,PetscEventRegLog eventRegLog,PetscClassId classid)
460: {
461: int e;
463: for (e = 0; e < eventLog->numEvents; e++) {
464: int c = eventRegLog->eventInfo[e].classid;
465: if (c == classid) eventLog->eventInfo[e].active = PETSC_TRUE;
466: }
467: return 0;
468: }
470: /*@C
471: PetscEventPerfLogDeactivateClass - Deactivates event logging for a PETSc object class.
473: Not Collective
475: Input Parameters:
476: + eventLog - The PetscEventPerfLog
477: . eventRegLog - The PetscEventRegLog
478: - classid - The class id, for example MAT_CLASSID, SNES_CLASSID,
480: Level: developer
482: .seealso: PetscEventPerfLogDeactivateClass(), PetscEventPerfLogDeactivate(), PetscEventPerfLogActivate()
483: @*/
484: PetscErrorCode PetscEventPerfLogDeactivateClass(PetscEventPerfLog eventLog,PetscEventRegLog eventRegLog,PetscClassId classid)
485: {
486: int e;
488: for (e = 0; e < eventLog->numEvents; e++) {
489: int c = eventRegLog->eventInfo[e].classid;
490: if (c == classid) eventLog->eventInfo[e].active = PETSC_FALSE;
491: }
492: return 0;
493: }
495: /*------------------------------------------------ Query Functions --------------------------------------------------*/
496: /*@C
497: PetscEventRegLogGetEvent - This function returns the event id given the event name.
499: Not Collective
501: Input Parameters:
502: + eventLog - The PetscEventRegLog
503: - name - The stage name
505: Output Parameter:
506: . event - The event id, or -1 if not found
508: Level: developer
510: .seealso: PetscEventRegLogRegister()
511: @*/
512: PetscErrorCode PetscEventRegLogGetEvent(PetscEventRegLog eventLog,const char name[],PetscLogEvent *event)
513: {
514: PetscBool match;
515: int e;
519: *event = -1;
520: for (e = 0; e < eventLog->numEvents; e++) {
521: PetscStrcasecmp(eventLog->eventInfo[e].name,name,&match);
522: if (match) {
523: *event = e;
524: break;
525: }
526: }
527: return 0;
528: }
530: /*@C
531: PetscEventPerfLogSetVisible - This function determines whether an event is printed during PetscLogView()
533: Not Collective
535: Input Parameters:
536: + eventLog - The PetscEventPerfLog
537: . event - The event to log
538: - isVisible - The visibility flag, PETSC_TRUE for printing, otherwise PETSC_FALSE (default is PETSC_TRUE)
540: Database Options:
541: . -log_view - Activates log summary
543: Level: developer
545: .seealso: PetscEventPerfLogGetVisible(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
546: @*/
547: PetscErrorCode PetscEventPerfLogSetVisible(PetscEventPerfLog eventLog,PetscLogEvent event,PetscBool isVisible)
548: {
549: eventLog->eventInfo[event].visible = isVisible;
550: return 0;
551: }
553: /*@C
554: PetscEventPerfLogGetVisible - This function returns whether an event is printed during PetscLogView()
556: Not Collective
558: Input Parameters:
559: + eventLog - The PetscEventPerfLog
560: - event - The event id to log
562: Output Parameter:
563: . isVisible - The visibility flag, PETSC_TRUE for printing, otherwise PETSC_FALSE (default is PETSC_TRUE)
565: Database Options:
566: . -log_view - Activates log summary
568: Level: developer
570: .seealso: PetscEventPerfLogSetVisible(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
571: @*/
572: PetscErrorCode PetscEventPerfLogGetVisible(PetscEventPerfLog eventLog,PetscLogEvent event,PetscBool *isVisible)
573: {
575: *isVisible = eventLog->eventInfo[event].visible;
576: return 0;
577: }
579: /*@C
580: PetscLogEventGetPerfInfo - Return the performance information about the given event in the given stage
582: Input Parameters:
583: + stage - The stage number or PETSC_DETERMINE for the current stage
584: - event - The event number
586: Output Parameters:
587: . info - This structure is filled with the performance information
589: Level: Intermediate
591: .seealso: PetscLogEventGetFlops()
592: @*/
593: PetscErrorCode PetscLogEventGetPerfInfo(int stage,PetscLogEvent event,PetscEventPerfInfo *info)
594: {
595: PetscStageLog stageLog;
596: PetscEventPerfLog eventLog = NULL;
600: PetscLogGetStageLog(&stageLog);
601: if (stage < 0) PetscStageLogGetCurrent(stageLog,&stage);
602: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
603: *info = eventLog->eventInfo[event];
604: return 0;
605: }
607: PetscErrorCode PetscLogEventGetFlops(PetscLogEvent event,PetscLogDouble *flops)
608: {
609: PetscStageLog stageLog;
610: PetscEventPerfLog eventLog = NULL;
611: int stage;
614: PetscLogGetStageLog(&stageLog);
615: PetscStageLogGetCurrent(stageLog,&stage);
616: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
617: *flops = eventLog->eventInfo[event].flops;
618: return 0;
619: }
621: PetscErrorCode PetscLogEventZeroFlops(PetscLogEvent event)
622: {
623: PetscStageLog stageLog;
624: PetscEventPerfLog eventLog = NULL;
625: int stage;
627: PetscLogGetStageLog(&stageLog);
628: PetscStageLogGetCurrent(stageLog,&stage);
629: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
631: eventLog->eventInfo[event].flops = 0.0;
632: eventLog->eventInfo[event].flops2 = 0.0;
633: eventLog->eventInfo[event].flopsTmp = 0.0;
634: return 0;
635: }
637: PetscErrorCode PetscLogEventSynchronize(PetscLogEvent event,MPI_Comm comm)
638: {
639: PetscStageLog stageLog;
640: PetscEventRegLog eventRegLog;
641: PetscEventPerfLog eventLog = NULL;
642: int stage;
643: PetscLogDouble time = 0.0;
645: if (!PetscLogSyncOn || comm == MPI_COMM_NULL) return 0;
646: PetscLogGetStageLog(&stageLog);
647: PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
648: if (!eventRegLog->eventInfo[event].collective) return 0;
649: PetscStageLogGetCurrent(stageLog,&stage);
650: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
651: if (eventLog->eventInfo[event].depth > 0) return 0;
653: PetscTimeSubtract(&time);
654: MPI_Barrier(comm);
655: PetscTimeAdd(&time);
656: eventLog->eventInfo[event].syncTime += time;
657: return 0;
658: }
660: PetscErrorCode PetscLogEventBeginDefault(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
661: {
662: PetscStageLog stageLog;
663: PetscEventPerfLog eventLog = NULL;
664: int stage;
666: PetscLogGetStageLog(&stageLog);
667: PetscStageLogGetCurrent(stageLog,&stage);
668: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
669: /* Synchronization */
670: PetscLogEventSynchronize(event,PetscObjectComm(o1));
671: /* Check for double counting */
672: eventLog->eventInfo[event].depth++;
673: if (eventLog->eventInfo[event].depth > 1) return 0;
674: /* Log the performance info */
675: eventLog->eventInfo[event].count++;
676: eventLog->eventInfo[event].timeTmp = 0.0;
677: PetscTimeSubtract(&eventLog->eventInfo[event].timeTmp);
678: eventLog->eventInfo[event].flopsTmp = -petsc_TotalFlops;
679: eventLog->eventInfo[event].numMessages -= petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct;
680: eventLog->eventInfo[event].messageLength -= petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
681: eventLog->eventInfo[event].numReductions -= petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
682: if (PetscLogMemory) {
683: PetscLogDouble usage;
684: PetscMemoryGetCurrentUsage(&usage);
685: eventLog->eventInfo[event].memIncrease -= usage;
686: PetscMallocGetCurrentUsage(&usage);
687: eventLog->eventInfo[event].mallocSpace -= usage;
688: PetscMallocGetMaximumUsage(&usage);
689: eventLog->eventInfo[event].mallocIncrease -= usage;
690: PetscMallocPushMaximumUsage((int)event);
691: }
692: #if defined(PETSC_HAVE_DEVICE)
693: eventLog->eventInfo[event].CpuToGpuCount -= petsc_ctog_ct;
694: eventLog->eventInfo[event].GpuToCpuCount -= petsc_gtoc_ct;
695: eventLog->eventInfo[event].CpuToGpuSize -= petsc_ctog_sz;
696: eventLog->eventInfo[event].GpuToCpuSize -= petsc_gtoc_sz;
697: eventLog->eventInfo[event].GpuFlops -= petsc_gflops;
698: eventLog->eventInfo[event].GpuTime -= petsc_gtime;
699: #endif
700: return 0;
701: }
703: PetscErrorCode PetscLogEventEndDefault(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
704: {
705: PetscStageLog stageLog;
706: PetscEventPerfLog eventLog = NULL;
707: int stage;
709: PetscLogGetStageLog(&stageLog);
710: PetscStageLogGetCurrent(stageLog,&stage);
711: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
712: /* Check for double counting */
713: eventLog->eventInfo[event].depth--;
714: if (eventLog->eventInfo[event].depth > 0) return 0;
716: /* Log performance info */
717: PetscTimeAdd(&eventLog->eventInfo[event].timeTmp);
718: eventLog->eventInfo[event].time += eventLog->eventInfo[event].timeTmp;
719: eventLog->eventInfo[event].time2 += eventLog->eventInfo[event].timeTmp*eventLog->eventInfo[event].timeTmp;
720: eventLog->eventInfo[event].flopsTmp += petsc_TotalFlops;
721: eventLog->eventInfo[event].flops += eventLog->eventInfo[event].flopsTmp;
722: eventLog->eventInfo[event].flops2 += eventLog->eventInfo[event].flopsTmp*eventLog->eventInfo[event].flopsTmp;
723: eventLog->eventInfo[event].numMessages += petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct;
724: eventLog->eventInfo[event].messageLength += petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
725: eventLog->eventInfo[event].numReductions += petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
726: if (PetscLogMemory) {
727: PetscLogDouble usage,musage;
728: PetscMemoryGetCurrentUsage(&usage);
729: eventLog->eventInfo[event].memIncrease += usage;
730: PetscMallocGetCurrentUsage(&usage);
731: eventLog->eventInfo[event].mallocSpace += usage;
732: PetscMallocPopMaximumUsage((int)event,&musage);
733: eventLog->eventInfo[event].mallocIncreaseEvent = PetscMax(musage-usage,eventLog->eventInfo[event].mallocIncreaseEvent);
734: PetscMallocGetMaximumUsage(&usage);
735: eventLog->eventInfo[event].mallocIncrease += usage;
736: }
737: #if defined(PETSC_HAVE_DEVICE)
738: eventLog->eventInfo[event].CpuToGpuCount += petsc_ctog_ct;
739: eventLog->eventInfo[event].GpuToCpuCount += petsc_gtoc_ct;
740: eventLog->eventInfo[event].CpuToGpuSize += petsc_ctog_sz;
741: eventLog->eventInfo[event].GpuToCpuSize += petsc_gtoc_sz;
742: eventLog->eventInfo[event].GpuFlops += petsc_gflops;
743: eventLog->eventInfo[event].GpuTime += petsc_gtime;
744: #endif
745: return 0;
746: }
748: PetscErrorCode PetscLogEventBeginComplete(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
749: {
750: PetscStageLog stageLog;
751: PetscEventRegLog eventRegLog;
752: PetscEventPerfLog eventPerfLog = NULL;
753: Action *tmpAction;
754: PetscLogDouble start,end;
755: PetscLogDouble curTime;
756: int stage;
758: /* Dynamically enlarge logging structures */
759: if (petsc_numActions >= petsc_maxActions) {
760: PetscTime(&start);
761: PetscCalloc1(petsc_maxActions*2,&tmpAction);
762: PetscArraycpy(tmpAction,petsc_actions,petsc_maxActions);
763: PetscFree(petsc_actions);
765: petsc_actions = tmpAction;
766: petsc_maxActions *= 2;
767: PetscTime(&end);
768: petsc_BaseTime += (end - start);
769: }
770: /* Record the event */
771: PetscLogGetStageLog(&stageLog);
772: PetscStageLogGetCurrent(stageLog,&stage);
773: PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
774: PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
775: PetscTime(&curTime);
776: if (petsc_logActions) {
777: petsc_actions[petsc_numActions].time = curTime - petsc_BaseTime;
778: petsc_actions[petsc_numActions].action = ACTIONBEGIN;
779: petsc_actions[petsc_numActions].event = event;
780: petsc_actions[petsc_numActions].classid = eventRegLog->eventInfo[event].classid;
781: if (o1) petsc_actions[petsc_numActions].id1 = o1->id;
782: else petsc_actions[petsc_numActions].id1 = -1;
783: if (o2) petsc_actions[petsc_numActions].id2 = o2->id;
784: else petsc_actions[petsc_numActions].id2 = -1;
785: if (o3) petsc_actions[petsc_numActions].id3 = o3->id;
786: else petsc_actions[petsc_numActions].id3 = -1;
787: petsc_actions[petsc_numActions].flops = petsc_TotalFlops;
789: PetscMallocGetCurrentUsage(&petsc_actions[petsc_numActions].mem);
790: PetscMallocGetMaximumUsage(&petsc_actions[petsc_numActions].maxmem);
791: petsc_numActions++;
792: }
793: /* Check for double counting */
794: eventPerfLog->eventInfo[event].depth++;
795: if (eventPerfLog->eventInfo[event].depth > 1) return 0;
796: /* Log the performance info */
797: eventPerfLog->eventInfo[event].count++;
798: eventPerfLog->eventInfo[event].time -= curTime;
799: eventPerfLog->eventInfo[event].flops -= petsc_TotalFlops;
800: eventPerfLog->eventInfo[event].numMessages -= petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct;
801: eventPerfLog->eventInfo[event].messageLength -= petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
802: eventPerfLog->eventInfo[event].numReductions -= petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
803: return 0;
804: }
806: PetscErrorCode PetscLogEventEndComplete(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
807: {
808: PetscStageLog stageLog;
809: PetscEventRegLog eventRegLog;
810: PetscEventPerfLog eventPerfLog = NULL;
811: Action *tmpAction;
812: PetscLogDouble start,end;
813: PetscLogDouble curTime;
814: int stage;
816: /* Dynamically enlarge logging structures */
817: if (petsc_numActions >= petsc_maxActions) {
818: PetscTime(&start);
819: PetscCalloc1(petsc_maxActions*2,&tmpAction);
820: PetscArraycpy(tmpAction,petsc_actions,petsc_maxActions);
821: PetscFree(petsc_actions);
823: petsc_actions = tmpAction;
824: petsc_maxActions *= 2;
825: PetscTime(&end);
826: petsc_BaseTime += (end - start);
827: }
828: /* Record the event */
829: PetscLogGetStageLog(&stageLog);
830: PetscStageLogGetCurrent(stageLog,&stage);
831: PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
832: PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
833: PetscTime(&curTime);
834: if (petsc_logActions) {
835: petsc_actions[petsc_numActions].time = curTime - petsc_BaseTime;
836: petsc_actions[petsc_numActions].action = ACTIONEND;
837: petsc_actions[petsc_numActions].event = event;
838: petsc_actions[petsc_numActions].classid = eventRegLog->eventInfo[event].classid;
839: if (o1) petsc_actions[petsc_numActions].id1 = o1->id;
840: else petsc_actions[petsc_numActions].id1 = -1;
841: if (o2) petsc_actions[petsc_numActions].id2 = o2->id;
842: else petsc_actions[petsc_numActions].id2 = -1;
843: if (o3) petsc_actions[petsc_numActions].id3 = o3->id;
844: else petsc_actions[petsc_numActions].id3 = -1;
845: petsc_actions[petsc_numActions].flops = petsc_TotalFlops;
847: PetscMallocGetCurrentUsage(&petsc_actions[petsc_numActions].mem);
848: PetscMallocGetMaximumUsage(&petsc_actions[petsc_numActions].maxmem);
849: petsc_numActions++;
850: }
851: /* Check for double counting */
852: eventPerfLog->eventInfo[event].depth--;
853: if (eventPerfLog->eventInfo[event].depth > 0) return 0;
855: /* Log the performance info */
856: eventPerfLog->eventInfo[event].count++;
857: eventPerfLog->eventInfo[event].time += curTime;
858: eventPerfLog->eventInfo[event].flops += petsc_TotalFlops;
859: eventPerfLog->eventInfo[event].numMessages += petsc_irecv_ct + petsc_isend_ct + petsc_recv_ct + petsc_send_ct;
860: eventPerfLog->eventInfo[event].messageLength += petsc_irecv_len + petsc_isend_len + petsc_recv_len + petsc_send_len;
861: eventPerfLog->eventInfo[event].numReductions += petsc_allreduce_ct + petsc_gather_ct + petsc_scatter_ct;
862: return 0;
863: }
865: PetscErrorCode PetscLogEventBeginTrace(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
866: {
867: PetscStageLog stageLog;
868: PetscEventRegLog eventRegLog;
869: PetscEventPerfLog eventPerfLog = NULL;
870: PetscLogDouble cur_time;
871: PetscMPIInt rank;
872: int stage,err;
874: if (!petsc_tracetime) PetscTime(&petsc_tracetime);
876: petsc_tracelevel++;
877: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
878: PetscLogGetStageLog(&stageLog);
879: PetscStageLogGetCurrent(stageLog,&stage);
880: PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
881: PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
882: /* Check for double counting */
883: eventPerfLog->eventInfo[event].depth++;
884: if (eventPerfLog->eventInfo[event].depth > 1) return 0;
885: /* Log performance info */
886: PetscTime(&cur_time);
887: PetscFPrintf(PETSC_COMM_SELF,petsc_tracefile,"%s[%d] %g Event begin: %s\n",petsc_tracespace,rank,cur_time-petsc_tracetime,eventRegLog->eventInfo[event].name);
888: PetscStrncpy(petsc_tracespace,petsc_traceblanks,2*petsc_tracelevel);
890: petsc_tracespace[2*petsc_tracelevel] = 0;
892: err = fflush(petsc_tracefile);
894: return 0;
895: }
897: PetscErrorCode PetscLogEventEndTrace(PetscLogEvent event,int t,PetscObject o1,PetscObject o2,PetscObject o3,PetscObject o4)
898: {
899: PetscStageLog stageLog;
900: PetscEventRegLog eventRegLog;
901: PetscEventPerfLog eventPerfLog = NULL;
902: PetscLogDouble cur_time;
903: int stage,err;
904: PetscMPIInt rank;
906: petsc_tracelevel--;
907: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
908: PetscLogGetStageLog(&stageLog);
909: PetscStageLogGetCurrent(stageLog,&stage);
910: PetscStageLogGetEventRegLog(stageLog,&eventRegLog);
911: PetscStageLogGetEventPerfLog(stageLog,stage,&eventPerfLog);
912: /* Check for double counting */
913: eventPerfLog->eventInfo[event].depth--;
914: if (eventPerfLog->eventInfo[event].depth > 0) return 0;
917: /* Log performance info */
918: if (petsc_tracelevel) {
919: PetscStrncpy(petsc_tracespace,petsc_traceblanks,2*petsc_tracelevel);
920: }
921: petsc_tracespace[2*petsc_tracelevel] = 0;
922: PetscTime(&cur_time);
923: PetscFPrintf(PETSC_COMM_SELF,petsc_tracefile,"%s[%d] %g Event end: %s\n",petsc_tracespace,rank,cur_time-petsc_tracetime,eventRegLog->eventInfo[event].name);
924: err = fflush(petsc_tracefile);
926: return 0;
927: }
929: /*@C
930: PetscLogEventSetDof - Set the nth number of degrees of freedom associated with this event
932: Not Collective
934: Input Parameters:
935: + event - The event id to log
936: . n - The dof index, in [0, 8)
937: - dof - The number of dofs
939: Database Options:
940: . -log_view - Activates log summary
942: Note: This is to enable logging of convergence
944: Level: developer
946: .seealso: PetscLogEventSetError(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
947: @*/
948: PetscErrorCode PetscLogEventSetDof(PetscLogEvent event, PetscInt n, PetscLogDouble dof)
949: {
950: PetscStageLog stageLog;
951: PetscEventPerfLog eventLog = NULL;
952: int stage;
955: PetscLogGetStageLog(&stageLog);
956: PetscStageLogGetCurrent(stageLog,&stage);
957: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
958: eventLog->eventInfo[event].dof[n] = dof;
959: return 0;
960: }
962: /*@C
963: PetscLogEventSetError - Set the nth error associated with this event
965: Not Collective
967: Input Parameters:
968: + event - The event id to log
969: . n - The error index, in [0, 8)
970: - error - The error
972: Database Options:
973: . -log_view - Activates log summary
975: Note: This is to enable logging of convergence, and enable users to interpret the errors as they wish. For example,
976: as different norms, or as errors for different fields
978: Level: developer
980: .seealso: PetscLogEventSetDof(), PetscEventRegLogRegister(), PetscStageLogGetEventLog()
981: @*/
982: PetscErrorCode PetscLogEventSetError(PetscLogEvent event, PetscInt n, PetscLogDouble error)
983: {
984: PetscStageLog stageLog;
985: PetscEventPerfLog eventLog = NULL;
986: int stage;
989: PetscLogGetStageLog(&stageLog);
990: PetscStageLogGetCurrent(stageLog,&stage);
991: PetscStageLogGetEventPerfLog(stageLog,stage,&eventLog);
992: eventLog->eventInfo[event].errors[n] = error;
993: return 0;
994: }