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,mx,my,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,"---------------------- \n");
59: PetscPrintf(PETSC_COMM_WORLD,"Processors: %d \n",size);
61: PetscPrintf(PETSC_COMM_WORLD,"Mesh: ");
62: for (i=0;i<nlevels;i++){
63: info = DAAppGetDA(daapplication,i,&da);CHKERRQ(info);
64: info = DAGetInfo(da,PETSC_NULL,&mx,&my,PETSC_NULL,
65: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
66: PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(info);
67: PetscPrintf(PETSC_COMM_WORLD," &%4d,%4d",mx,my);
68: }
69: PetscPrintf(PETSC_COMM_WORLD," \n");
70: PetscPrintf(PETSC_COMM_WORLD,"Times: ");
71: for (i=0;i<nlevels;i++){
72: t0=dactx->tstageend[i] - dactx->tstagebegin[i];
73: PetscPrintf(PETSC_COMM_WORLD," & %7.2f ",t0);
74: }
75: PetscPrintf(PETSC_COMM_WORLD," \n");
76: PetscPrintf(PETSC_COMM_WORLD,"Function Evaluations: ");
77: for (i=0;i<nlevels;i++){
78: PetscPrintf(PETSC_COMM_WORLD," & %7d ",dactx->nfeval[i]);
79: }
80: PetscPrintf(PETSC_COMM_WORLD," \n");
81: PetscPrintf(PETSC_COMM_WORLD,"Gradient Evaluations: ");
82: for (i=0;i<nlevels;i++){
83: PetscPrintf(PETSC_COMM_WORLD," & %7d ",dactx->ngeval[i]);
84: }
85: PetscPrintf(PETSC_COMM_WORLD," \n");
86: PetscPrintf(PETSC_COMM_WORLD,"Hessian Evaluations: ");
87: for (i=0;i<nlevels;i++){
88: PetscPrintf(PETSC_COMM_WORLD," & %7d ",dactx->nheval[i]);
89: }
90: PetscPrintf(PETSC_COMM_WORLD," \n");
91: /*
92: PetscPrintf(PETSC_COMM_WORLD,"Linear Solves: %6d ",size);
93: for (i=0;i<nlevels;i++){
94: PetscPrintf(PETSC_COMM_WORLD," & %4d ",dactx->nlsolves[i]);
95: }
96: PetscPrintf(PETSC_COMM_WORLD," \n");
97: */
98: PetscPrintf(PETSC_COMM_WORLD,"---------------------- \n");
99: }
100: return(0);
101: }
105: int DAAppTimeMonitorDestroy(void *ctx){
106: int info;
108: info=TaoApplicationFreeMemory(ctx);CHKERRQ(info);
109: return(0);
110: }
114: int DAAppPrintStageTimes(TAO_APPLICATION taoapp){
115: int info;
116: DAAppTimeMonitor *dactx;
118: PetscNew(DAAppTimeMonitor,&dactx);
119: info = DAAppSetBeforeMonitor(taoapp,DAAppTimeMonitorBefore,(void*)dactx);CHKERRQ(info);
120: info = DAAppSetAfterMonitor(taoapp,DAAppTimeMonitorAfter,(void*)dactx);CHKERRQ(info);
121: info = DAAppSetAfterMonitor(taoapp,DAAppTimeMonitorAfterAll,(void*)dactx);CHKERRQ(info);
122: info = TaoAppSetDestroyRoutine(taoapp,DAAppTimeMonitorDestroy,(void*)dactx); CHKERRQ(info);
123: return(0);
124: }
125:
128: int DAAppInterpolationMonitor(TAO_APPLICATION taoapp, DA da, int level, void *ctx){
129:
130: int info,mx,my,n;
131: PetscScalar dd=-1.0;
132: Vec X,XCoarse,XX;
133: Mat Interpolate;
136: info = DAGetInfo(da,PETSC_NULL,&mx,&my,PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL,
137: PETSC_NULL,PETSC_NULL,PETSC_NULL,PETSC_NULL);CHKERRQ(info);
138: if (level>0){
139: info = DAAppGetInterpolationMatrix(taoapp,level,&Interpolate,0); CHKERRQ(info);
140: info = DAAppGetSolution(taoapp,level,&X); CHKERRQ(info);
141: info = DAAppGetSolution(taoapp,level-1,&XCoarse); CHKERRQ(info);
142: info = VecDuplicate(X,&XX); CHKERRQ(info);
143: info = MatMult(Interpolate,XCoarse,XX); CHKERRQ(info);
144: info = VecAXPY(XX, dd, X); CHKERRQ(info);
145: info = VecNorm(XX,NORM_INFINITY,&dd); CHKERRQ(info);
146: PetscPrintf(MPI_COMM_WORLD,"Maximum Interpolation Error: %4.2e\n",dd);
147: info = VecNorm(XX,NORM_1,&dd); CHKERRQ(info);
148: info = VecGetSize(XX,&n); CHKERRQ(info);
149: PetscPrintf(MPI_COMM_WORLD,"Average Interpolation Error: %4.2e\n",dd/n);
150: info = VecDestroy(XX); CHKERRQ(info);
151: }
152: return(0);
153: return 0;
154: }
159: int DAAppPrintInterpolationError(TAO_APPLICATION taoapp){
160: int info;
162: info = DAAppSetAfterMonitor(taoapp,DAAppInterpolationMonitor,0); CHKERRQ(info);
163: return(0);
164: }