Actual source code: dmnetworkts.c
petsc-3.14.6 2021-03-30
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: {
17: PetscInt i;
20: for (i=0; i<(*ctx)->nlg; i++) {
21: PetscDrawLGDestroy(&(*ctx)->lg[i]);
22: }
23: PetscFree((*ctx)->lg);
24: PetscFree(*ctx);
25: return(0);
26: }
28: PetscErrorCode TSMonitorLGCtxNetworkCreate(TS ts,const char host[],const char label[],int x,int y,int m,int n,PetscInt howoften,TSMonitorLGCtxNetwork *ctx)
29: {
30: PetscDraw draw;
32: MPI_Comm comm;
33: DM dm;
34: PetscInt i,Start,End,e,nvar;
37: TSGetDM(ts,&dm);
38: PetscObjectGetComm((PetscObject)ts,&comm);
39: PetscNew(ctx);
40: i = 0;
41: /* loop over edges counting number of line graphs needed */
42: DMNetworkGetEdgeRange(dm,&Start,&End);
43: for (e=Start; e<End; e++) {
44: DMNetworkGetNumVariables(dm,e,&nvar);
45: if (!nvar) continue;
46: i++;
47: }
48: /* loop over vertices */
49: DMNetworkGetVertexRange(dm,&Start,&End);
50: for (e=Start; e<End; e++) {
51: DMNetworkGetNumVariables(dm,e,&nvar);
52: if (!nvar) continue;
53: i++;
54: }
55: (*ctx)->nlg = i;
56: PetscMalloc1(i,&(*ctx)->lg);
58: i = 0;
59: /* loop over edges creating all needed line graphs*/
60: DMNetworkGetEdgeRange(dm,&Start,&End);
61: for (e=Start; e<End; e++) {
62: DMNetworkGetNumVariables(dm,e,&nvar);
63: if (!nvar) continue;
64: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
65: PetscDrawSetFromOptions(draw);
66: PetscDrawLGCreate(draw,nvar,&(*ctx)->lg[i]);
67: PetscDrawLGSetFromOptions((*ctx)->lg[i]);
68: PetscDrawDestroy(&draw);
69: i++;
70: }
71: /* loop over vertices */
72: DMNetworkGetVertexRange(dm,&Start,&End);
73: for (e=Start; e<End; e++) {
74: DMNetworkGetNumVariables(dm,e,&nvar);
75: if (!nvar) continue;
76: PetscDrawCreate(comm,host,label,x,y,m,n,&draw);
77: PetscDrawSetFromOptions(draw);
78: PetscDrawLGCreate(draw,nvar,&(*ctx)->lg[i]);
79: PetscDrawLGSetFromOptions((*ctx)->lg[i]);
80: PetscDrawDestroy(&draw);
81: i++;
82: }
83: PetscDrawDestroy(&draw);
84: (*ctx)->howoften = howoften;
85: return(0);
86: }
88: /*
89: TSMonitorLGCtxNetworkSolution - Monitors progress of the TS solvers for a DMNetwork solution with one window for each vertex and each edge
91: Collective on TS
93: Input Parameters:
94: + ts - the TS context
95: . step - current time-step
96: . ptime - current time
97: . u - current solution
98: - dctx - the TSMonitorLGCtxNetwork object that contains all the options for the monitoring, this is created with TSMonitorLGCtxCreateNetwork()
100: Options Database:
101: . -ts_monitor_lg_solution_variables
103: Level: intermediate
105: Notes:
106: Each process in a parallel run displays its component solutions in a separate window
108: */
109: PetscErrorCode TSMonitorLGCtxNetworkSolution(TS ts,PetscInt step,PetscReal ptime,Vec u,void *dctx)
110: {
111: PetscErrorCode ierr;
112: TSMonitorLGCtxNetwork ctx = (TSMonitorLGCtxNetwork)dctx;
113: const PetscScalar *xv;
114: PetscScalar *yv;
115: PetscInt i,v,Start,End,offset,nvar,e;
116: TSConvergedReason reason;
117: DM dm;
118: Vec uv;
121: if (step < 0) return(0); /* -1 indicates interpolated solution */
122: if (!step) {
123: PetscDrawAxis axis;
125: for (i=0; i<ctx->nlg; i++) {
126: PetscDrawLGGetAxis(ctx->lg[i],&axis);
127: PetscDrawAxisSetLabels(axis,"Solution as function of time","Time","Solution");
128: PetscDrawLGReset(ctx->lg[i]);
129: }
130: }
132: if (ctx->semilogy) {
133: PetscInt n,j;
135: VecDuplicate(u,&uv);
136: VecCopy(u,uv);
137: VecGetArray(uv,&yv);
138: VecGetLocalSize(uv,&n);
139: for (j=0; j<n; j++) {
140: if (PetscRealPart(yv[j]) <= 0.0) yv[j] = -12;
141: else yv[j] = PetscLog10Real(PetscRealPart(yv[j]));
142: }
143: xv = yv;
144: } else {
145: VecGetArrayRead(u,&xv);
146: }
147: /* iterate over edges */
148: TSGetDM(ts,&dm);
149: i = 0;
150: DMNetworkGetEdgeRange(dm,&Start,&End);
151: for (e=Start; e<End; e++) {
152: DMNetworkGetNumVariables(dm,e,&nvar);
153: if (!nvar) continue;
155: DMNetworkGetVariableOffset(dm,e,&offset);
156: PetscDrawLGAddCommonPoint(ctx->lg[i],ptime,(const PetscReal*)(xv+offset));
157: i++;
158: }
160: /* iterate over vertices */
161: DMNetworkGetVertexRange(dm,&Start,&End);
162: for (v=Start; v<End; v++) {
163: DMNetworkGetNumVariables(dm,v,&nvar);
164: if (!nvar) continue;
166: DMNetworkGetVariableOffset(dm,v,&offset);
167: PetscDrawLGAddCommonPoint(ctx->lg[i],ptime,(const PetscReal*)(xv+offset));
168: i++;
169: }
170: if (ctx->semilogy) {
171: VecRestoreArray(uv,&yv);
172: VecDestroy(&uv);
173: } else {
174: VecRestoreArrayRead(u,&xv);
175: }
177: TSGetConvergedReason(ts,&reason);
178: if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && reason)) {
179: for (i=0; i<ctx->nlg; i++) {
180: PetscDrawLGDraw(ctx->lg[i]);
181: PetscDrawLGSave(ctx->lg[i]);
182: }
183: }
184: return(0);
185: }