Actual source code: daapp_monitor.c
1: #include "tao.h"
2: #include "petsc.h"
3: #include "daapp_impl.h"
4: #include taodaapplication.h
6: typedef struct {
7: PetscLogDouble tstagebegin[20],tstageend[20];
8: int nfeval[20],ngeval[20],nheval[20],nlsolves[20];
9: } DAAppTimeMonitor;
15: int DAAppTimeMonitorBefore(TAO_APPLICATION daapplication, DA da, int level, void *ctx){
16: int info;
17: DAAppTimeMonitor *dactx=(DAAppTimeMonitor*)ctx;
18: PetscLogDouble t2;
20: info=PetscGetTime(&t2);CHKERRQ(info);
21: dactx->tstagebegin[level]=t2;
22: info=TaoAppResetCounters(daapplication);CHKERRQ(info);
23: return(0);
24: }
27: int DAAppTimeMonitorAfter(TAO_APPLICATION daapplication, DA da, int level, void *ctx){
28: int info,stats[4];
29: DA_APPLICATION daapp;
30: DAAppTimeMonitor *dactx=(DAAppTimeMonitor*)ctx;
31: PetscLogDouble t0,t1,t2;
33: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
34: info=PetscGetTime(&t2);CHKERRQ(info);
35: dactx->tstageend[level]=t2;
36: t0=t2 - dactx->tstagebegin[0];
37: t1=t2 - dactx->tstagebegin[level];
38: info=PetscPrintf(PETSC_COMM_WORLD," Grid %d: Time: %4.4e \n",level,t1);
39: info=PetscPrintf(PETSC_COMM_WORLD," Total Time: %4.4e \n",t0);
40: info=TaoAppCounters(daapplication,stats);CHKERRQ(info);
41: dactx->nfeval[level]=stats[0];
42: dactx->ngeval[level]=stats[1];
43: dactx->nheval[level]=stats[2];
44: dactx->nlsolves[level]=stats[3];
45: return(0);
46: }
50: int DAAppTimeMonitorAfterAll(TAO_APPLICATION daapplication, DA da, int level, void *ctx){
51: int i,info,nlevels,size;
52: DAAppTimeMonitor *dactx=(DAAppTimeMonitor*)ctx;
53: PetscLogDouble t0;
55: info=DAAppGetNumberOfDAGrids(daapplication,&nlevels);CHKERRQ(info);
56: if (nlevels==level+1){
57: MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(info);
58: PetscPrintf(PETSC_COMM_WORLD,"Stage Times: %6d ",size);
59: for (i=0;i<nlevels;i++){
60: t0=dactx->tstageend[i] - dactx->tstagebegin[i];
61: PetscPrintf(PETSC_COMM_WORLD," & %8.2f ",t0);
62: }
63: PetscPrintf(PETSC_COMM_WORLD," \n");
64: PetscPrintf(PETSC_COMM_WORLD,"Function Evaluations: %6d ",size);
65: for (i=0;i<nlevels;i++){
66: PetscPrintf(PETSC_COMM_WORLD," & %4d ",dactx->nfeval[i]);
67: }
68: PetscPrintf(PETSC_COMM_WORLD," \n");
69: PetscPrintf(PETSC_COMM_WORLD,"Gradient Evaluations: %6d ",size);
70: for (i=0;i<nlevels;i++){
71: PetscPrintf(PETSC_COMM_WORLD," & %4d ",dactx->ngeval[i]);
72: }
73: PetscPrintf(PETSC_COMM_WORLD," \n");
74: PetscPrintf(PETSC_COMM_WORLD,"Hessian Evaluations: %6d ",size);
75: for (i=0;i<nlevels;i++){
76: PetscPrintf(PETSC_COMM_WORLD," & %4d ",dactx->nheval[i]);
77: }
78: PetscPrintf(PETSC_COMM_WORLD," \n");
79: /*
80: PetscPrintf(PETSC_COMM_WORLD,"Linear Solves: %6d ",size);
81: for (i=0;i<nlevels;i++){
82: PetscPrintf(PETSC_COMM_WORLD," & %4d ",dactx->nlsolves[i]);
83: }
84: PetscPrintf(PETSC_COMM_WORLD," \n");
85: */
86: }
87: return(0);
88: }
92: int DAAppTimeMonitorDestroy(void *ctx){
93: int info;
95: info=TaoApplicationFreeMemory(ctx);CHKERRQ(info);
96: return(0);
97: }
101: int DAAppPrintStageTimes(TAO_APPLICATION taoapp){
102: int info;
103: DAAppTimeMonitor *dactx;
105: PetscNew(DAAppTimeMonitor,&dactx);
106: info = DAAppSetBeforeMonitor(taoapp,DAAppTimeMonitorBefore,(void*)dactx);CHKERRQ(info);
107: info = DAAppSetAfterMonitor(taoapp,DAAppTimeMonitorAfter,(void*)dactx);CHKERRQ(info);
108: info = DAAppSetAfterMonitor(taoapp,DAAppTimeMonitorAfterAll,(void*)dactx);CHKERRQ(info);
109: info = TaoAppSetDestroyRoutine(taoapp,DAAppTimeMonitorDestroy,(void*)dactx); CHKERRQ(info);
110: return(0);
111: }
112:
115: int DAAppInterpolationMonitor(TAO_APPLICATION taoapp, DA da, int level, void *ctx){
116:
117: int info,mx,my,n;
118: PetscScalar dd=-1.0;
119: Vec X,XCoarse,XX;
120: Mat Interpolate;
123: info = DAGetInfo(da,PETSC_NULL,&mx,&my,PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
124: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(info);
125: if (level>0){
126: info = DAAppGetInterpolationMatrix(taoapp,level,&Interpolate,0); CHKERRQ(info);
127: info = DAAppGetSolution(taoapp,level,&X); CHKERRQ(info);
128: info = DAAppGetSolution(taoapp,level-1,&XCoarse); CHKERRQ(info);
129: info = VecDuplicate(X,&XX); CHKERRQ(info);
130: info = MatMult(Interpolate,XCoarse,XX); CHKERRQ(info);
131: info = VecAXPY(XX, dd, X); CHKERRQ(info);
132: info = VecNorm(XX,NORM_INFINITY,&dd); CHKERRQ(info);
133: PetscPrintf(MPI_COMM_WORLD,"Maximum Interpolation Error: %4.2e\n",dd);
134: info = VecNorm(XX,NORM_1,&dd); CHKERRQ(info);
135: info = VecGetSize(XX,&n); CHKERRQ(info);
136: PetscPrintf(MPI_COMM_WORLD,"Average Interpolation Error: %4.2e\n",dd/n);
137: info = VecDestroy(XX); CHKERRQ(info);
138: }
139: return(0);
140: return 0;
141: }
146: int DAAppPrintInterpolationError(TAO_APPLICATION taoapp){
147: int info;
149: info = DAAppSetAfterMonitor(taoapp,DAAppInterpolationMonitor,0); CHKERRQ(info);
150: return(0);
151: }