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: PetscInt nfeval[20],ngeval[20],nheval[20],nlsolves[20];
9: } DAAppTimeMonitor;
15: int DAAppTimeMonitorBefore(TAO_APPLICATION daapplication, DA da, PetscInt 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, PetscInt level, void *ctx){
28: int info;
29: PetscInt stats[4];
30: DA_APPLICATION daapp;
31: DAAppTimeMonitor *dactx=(DAAppTimeMonitor*)ctx;
32: PetscLogDouble t0,t1,t2;
34: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
35: info=PetscGetTime(&t2);CHKERRQ(info);
36: dactx->tstageend[level]=t2;
37: t0=t2 - dactx->tstagebegin[0];
38: t1=t2 - dactx->tstagebegin[level];
39: info=PetscPrintf(PETSC_COMM_WORLD," Grid %d: Time: %4.4e \n",level,t1);
40: info=PetscPrintf(PETSC_COMM_WORLD," Total Time: %4.4e \n",t0);
41: info=TaoAppCounters(daapplication,stats);CHKERRQ(info);
42: dactx->nfeval[level]=stats[0];
43: dactx->ngeval[level]=stats[1];
44: dactx->nheval[level]=stats[2];
45: dactx->nlsolves[level]=stats[3];
46: return(0);
47: }
51: int DAAppTimeMonitorAfterAll(TAO_APPLICATION daapplication, DA da, PetscInt level, void *ctx){
52: int info,size;
53: PetscInt i,nlevels;
54: PetscInt mx,my;
55: DAAppTimeMonitor *dactx=(DAAppTimeMonitor*)ctx;
56: PetscLogDouble t0;
58: info=DAAppGetNumberOfDAGrids(daapplication,&nlevels);CHKERRQ(info);
59: if (nlevels==level+1){
60: MPI_Comm_size(PETSC_COMM_WORLD,&size);CHKERRQ(info);
61: PetscPrintf(PETSC_COMM_WORLD,"---------------------- \n");
62: PetscPrintf(PETSC_COMM_WORLD,"Processors: %d \n",(PetscInt)size);
64: PetscPrintf(PETSC_COMM_WORLD,"Mesh: ");
65: for (i=0;i<nlevels;i++){
66: info = DAAppGetDA(daapplication,i,&da);CHKERRQ(info);
67: info = DAGetInfo(da,PETSC_NULL,&mx,&my,PETSC_NULL,
68: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
69: PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(info);
70: PetscPrintf(PETSC_COMM_WORLD," &%4d,%4d",mx,my);
71: }
72: PetscPrintf(PETSC_COMM_WORLD," \n");
73: PetscPrintf(PETSC_COMM_WORLD,"Times: ");
74: for (i=0;i<nlevels;i++){
75: t0=dactx->tstageend[i] - dactx->tstagebegin[i];
76: PetscPrintf(PETSC_COMM_WORLD," & %7.2f ",t0);
77: }
78: PetscPrintf(PETSC_COMM_WORLD," \n");
79: PetscPrintf(PETSC_COMM_WORLD,"Function Evaluations: ");
80: for (i=0;i<nlevels;i++){
81: PetscPrintf(PETSC_COMM_WORLD," & %7d ",dactx->nfeval[i]);
82: }
83: PetscPrintf(PETSC_COMM_WORLD," \n");
84: PetscPrintf(PETSC_COMM_WORLD,"Gradient Evaluations: ");
85: for (i=0;i<nlevels;i++){
86: PetscPrintf(PETSC_COMM_WORLD," & %7d ",dactx->ngeval[i]);
87: }
88: PetscPrintf(PETSC_COMM_WORLD," \n");
89: PetscPrintf(PETSC_COMM_WORLD,"Hessian Evaluations: ");
90: for (i=0;i<nlevels;i++){
91: PetscPrintf(PETSC_COMM_WORLD," & %7d ",dactx->nheval[i]);
92: }
93: PetscPrintf(PETSC_COMM_WORLD," \n");
94: /*
95: PetscPrintf(PETSC_COMM_WORLD,"Linear Solves: %6d ",size);
96: for (i=0;i<nlevels;i++){
97: PetscPrintf(PETSC_COMM_WORLD," & %4d ",dactx->nlsolves[i]);
98: }
99: PetscPrintf(PETSC_COMM_WORLD," \n");
100: */
101: PetscPrintf(PETSC_COMM_WORLD,"---------------------- \n");
102: }
103: return(0);
104: }
108: int DAAppTimeMonitorDestroy(void *ctx){
109: int info;
111: info=TaoApplicationFreeMemory(ctx);CHKERRQ(info);
112: return(0);
113: }
117: int DAAppPrintStageTimes(TAO_APPLICATION taoapp){
118: int info;
119: DAAppTimeMonitor *dactx;
121: PetscNew(DAAppTimeMonitor,&dactx);
122: info = DAAppSetBeforeMonitor(taoapp,DAAppTimeMonitorBefore,(void*)dactx);CHKERRQ(info);
123: info = DAAppSetAfterMonitor(taoapp,DAAppTimeMonitorAfter,(void*)dactx);CHKERRQ(info);
124: info = DAAppSetAfterMonitor(taoapp,DAAppTimeMonitorAfterAll,(void*)dactx);CHKERRQ(info);
125: info = TaoAppSetDestroyRoutine(taoapp,DAAppTimeMonitorDestroy,(void*)dactx); CHKERRQ(info);
126: return(0);
127: }
128:
131: int DAAppInterpolationMonitor(TAO_APPLICATION taoapp, DA da, PetscInt level, void *ctx){
132:
133: int info;
134: PetscInt mx,my,n;
135: PetscScalar dd=-1.0;
136: Vec X,XCoarse,XX;
137: Mat Interpolate;
140: info = DAGetInfo(da,PETSC_NULL,&mx,&my,PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
141: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(info);
142: if (level>0){
143: info = DAAppGetInterpolationMatrix(taoapp,level,&Interpolate,0); CHKERRQ(info);
144: info = DAAppGetSolution(taoapp,level,&X); CHKERRQ(info);
145: info = DAAppGetSolution(taoapp,level-1,&XCoarse); CHKERRQ(info);
146: info = VecDuplicate(X,&XX); CHKERRQ(info);
147: info = MatMult(Interpolate,XCoarse,XX); CHKERRQ(info);
148: info = VecAXPY(XX, dd, X); CHKERRQ(info);
149: info = VecNorm(XX,NORM_INFINITY,&dd); CHKERRQ(info);
150: PetscPrintf(MPI_COMM_WORLD,"Maximum Interpolation Error: %4.2e\n",dd);
151: info = VecNorm(XX,NORM_1,&dd); CHKERRQ(info);
152: info = VecGetSize(XX,&n); CHKERRQ(info);
153: PetscPrintf(MPI_COMM_WORLD,"Average Interpolation Error: %4.2e\n",dd/n);
154: info = VecDestroy(XX); CHKERRQ(info);
155: }
156: return(0);
157: return 0;
158: }
163: int DAAppPrintInterpolationError(TAO_APPLICATION taoapp){
164: int info;
166: info = DAAppSetAfterMonitor(taoapp,DAAppInterpolationMonitor,0); CHKERRQ(info);
167: return(0);
168: }