Actual source code: dmnetworkts.c

  1: #include <petsc/private/dmnetworkimpl.h>
  2: #include <petscts.h>
  3: #include <petscdraw.h>

  5: /*
  6:    TSMonitorLGCtxDestroy - Destroys  line graph contexts that where created with TSMonitorLGCtxNetworkCreate().

  8:    Collective on TSMonitorLGCtx_Network

 10:    Input Parameter:
 11: .  ctx - the monitor context

 13: */
 14: PetscErrorCode  TSMonitorLGCtxNetworkDestroy(TSMonitorLGCtxNetwork *ctx)
 15: {
 16:   PetscInt       i;

 18:   for (i=0; i<(*ctx)->nlg; i++) {
 19:     PetscDrawLGDestroy(&(*ctx)->lg[i]);
 20:   }
 21:   PetscFree((*ctx)->lg);
 22:   PetscFree(*ctx);
 23:   return 0;
 24: }

 26: PetscErrorCode  TSMonitorLGCtxNetworkCreate(TS ts,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TSMonitorLGCtxNetwork *ctx)
 27: {
 28:   PetscDraw      draw;
 29:   MPI_Comm       comm;
 30:   DM             dm;
 31:   PetscInt       i,Start,End,e,nvar;

 33:   TSGetDM(ts,&dm);
 34:   PetscObjectGetComm((PetscObject)ts,&comm);
 35:   PetscNew(ctx);
 36:   i = 0;
 37:   /* loop over edges counting number of line graphs needed */
 38:   DMNetworkGetEdgeRange(dm,&Start,&End);
 39:   for (e=Start; e<End; e++) {
 40:     DMNetworkGetComponent(dm,e,ALL_COMPONENTS,NULL,NULL,&nvar);
 41:     if (!nvar) continue;
 42:     i++;
 43:   }
 44:   /* loop over vertices */
 45:   DMNetworkGetVertexRange(dm,&Start,&End);
 46:   for (e=Start; e<End; e++) {
 47:     DMNetworkGetComponent(dm,e,ALL_COMPONENTS,NULL,NULL,&nvar);
 48:     if (!nvar) continue;
 49:     i++;
 50:   }
 51:   (*ctx)->nlg = i;
 52:   PetscMalloc1(i,&(*ctx)->lg);

 54:   i = 0;
 55:   /* loop over edges creating all needed line graphs*/
 56:   DMNetworkGetEdgeRange(dm,&Start,&End);
 57:   for (e=Start; e<End; e++) {
 58:     DMNetworkGetComponent(dm,e,ALL_COMPONENTS,NULL,NULL,&nvar);
 59:     if (!nvar) continue;
 60:     PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
 61:     PetscDrawSetFromOptions(draw);
 62:     PetscDrawLGCreate(draw,nvar,&(*ctx)->lg[i]);
 63:     PetscDrawLGSetFromOptions((*ctx)->lg[i]);
 64:     PetscDrawDestroy(&draw);
 65:     i++;
 66:   }
 67:   /* loop over vertices */
 68:   DMNetworkGetVertexRange(dm,&Start,&End);
 69:   for (e=Start; e<End; e++) {
 70:     DMNetworkGetComponent(dm,e,ALL_COMPONENTS,NULL,NULL,&nvar);
 71:     if (!nvar) continue;
 72:     PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
 73:     PetscDrawSetFromOptions(draw);
 74:     PetscDrawLGCreate(draw,nvar,&(*ctx)->lg[i]);
 75:     PetscDrawLGSetFromOptions((*ctx)->lg[i]);
 76:     PetscDrawDestroy(&draw);
 77:     i++;
 78:   }
 79:   PetscDrawDestroy(&draw);
 80:   (*ctx)->howoften = howoften;
 81:   return 0;
 82: }

 84: /*
 85:    TSMonitorLGCtxNetworkSolution - Monitors progress of the TS solvers for a DMNetwork solution with one window for each vertex and each edge

 87:    Collective on TS

 89:    Input Parameters:
 90: +  ts - the TS context
 91: .  step - current time-step
 92: .  ptime - current time
 93: .  u - current solution
 94: -  dctx - the TSMonitorLGCtxNetwork object that contains all the options for the monitoring, this is created with TSMonitorLGCtxCreateNetwork()

 96:    Options Database:
 97: .   -ts_monitor_lg_solution_variables

 99:    Level: intermediate

101:    Notes:
102:     Each process in a parallel run displays its component solutions in a separate window

104: */
105: PetscErrorCode  TSMonitorLGCtxNetworkSolution(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dctx)
106: {
107:   TSMonitorLGCtxNetwork ctx = (TSMonitorLGCtxNetwork)dctx;
108:   const PetscScalar     *xv;
109:   PetscScalar           *yv;
110:   PetscInt              i,v,Start,End,offset,nvar,e;
111:   TSConvergedReason     reason;
112:   DM                    dm;
113:   Vec                   uv;

115:   if (step < 0) return 0; /* -1 indicates interpolated solution */
116:   if (!step) {
117:     PetscDrawAxis axis;

119:     for (i=0; i<ctx->nlg; i++) {
120:       PetscDrawLGGetAxis(ctx->lg[i],&axis);
121:       PetscDrawAxisSetLabels(axis,"Solution as function of time","Time","Solution");
122:       PetscDrawLGReset(ctx->lg[i]);
123:     }
124:   }

126:   if (ctx->semilogy) {
127:     PetscInt n,j;

129:     VecDuplicate(u,&uv);
130:     VecCopy(u,uv);
131:     VecGetArray(uv,&yv);
132:     VecGetLocalSize(uv,&n);
133:     for (j=0; j<n; j++) {
134:       if (PetscRealPart(yv[j]) <= 0.0) yv[j] = -12;
135:       else yv[j] = PetscLog10Real(PetscRealPart(yv[j]));
136:     }
137:     xv = yv;
138:   } else {
139:     VecGetArrayRead(u,&xv);
140:   }
141:   /* iterate over edges */
142:   TSGetDM(ts,&dm);
143:   i = 0;
144:   DMNetworkGetEdgeRange(dm,&Start,&End);
145:   for (e=Start; e<End; e++) {
146:     DMNetworkGetComponent(dm,e,ALL_COMPONENTS,NULL,NULL,&nvar);
147:     if (!nvar) continue;

149:     DMNetworkGetLocalVecOffset(dm,e,ALL_COMPONENTS,&offset);
150:     PetscDrawLGAddCommonPoint(ctx->lg[i],ptime,(const PetscReal*)(xv+offset));
151:     i++;
152:   }

154:   /* iterate over vertices */
155:   DMNetworkGetVertexRange(dm,&Start,&End);
156:   for (v=Start; v<End; v++) {
157:     DMNetworkGetComponent(dm,v,ALL_COMPONENTS,NULL,NULL,&nvar);
158:     if (!nvar) continue;

160:     DMNetworkGetLocalVecOffset(dm,v,ALL_COMPONENTS,&offset);
161:     PetscDrawLGAddCommonPoint(ctx->lg[i],ptime,(const PetscReal*)(xv+offset));
162:     i++;
163:   }
164:   if (ctx->semilogy) {
165:     VecRestoreArray(uv,&yv);
166:     VecDestroy(&uv);
167:   } else {
168:     VecRestoreArrayRead(u,&xv);
169:   }

171:   TSGetConvergedReason(ts,&reason);
172:   if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && reason)) {
173:     for (i=0; i<ctx->nlg; i++) {
174:       PetscDrawLGDraw(ctx->lg[i]);
175:       PetscDrawLGSave(ctx->lg[i]);
176:     }
177:   }
178:   return 0;
179: }