1 #include <petsc/private/dmnetworkimpl.h> /*I "petscdmplex.h" I*/ 2 #include <petscts.h> 3 #include <petscdraw.h> 4 5 /* 6 TSMonitorLGCtxDestroy - Destroys line graph contexts that where created with TSMonitorLGCtxNetworkCreate(). 7 8 Collective 9 10 Input Parameter: 11 . ctx - the monitor context 12 13 */ 14 PetscErrorCode TSMonitorLGCtxNetworkDestroy(TSMonitorLGCtxNetwork *ctx) 15 { 16 PetscInt i; 17 18 PetscFunctionBegin; 19 for (i = 0; i < (*ctx)->nlg; i++) PetscCall(PetscDrawLGDestroy(&(*ctx)->lg[i])); 20 PetscCall(PetscFree((*ctx)->lg)); 21 PetscCall(PetscFree(*ctx)); 22 PetscFunctionReturn(PETSC_SUCCESS); 23 } 24 25 PetscErrorCode TSMonitorLGCtxNetworkCreate(TS ts, const char host[], const char label[], int x, int y, int m, int n, PetscInt howoften, TSMonitorLGCtxNetwork *ctx) 26 { 27 PetscDraw draw; 28 MPI_Comm comm; 29 DM dm; 30 PetscInt i, Start, End, e, nvar; 31 32 PetscFunctionBegin; 33 PetscCall(TSGetDM(ts, &dm)); 34 PetscCall(PetscObjectGetComm((PetscObject)ts, &comm)); 35 PetscCall(PetscNew(ctx)); 36 i = 0; 37 /* loop over edges counting number of line graphs needed */ 38 PetscCall(DMNetworkGetEdgeRange(dm, &Start, &End)); 39 for (e = Start; e < End; e++) { 40 PetscCall(DMNetworkGetComponent(dm, e, ALL_COMPONENTS, NULL, NULL, &nvar)); 41 if (!nvar) continue; 42 i++; 43 } 44 /* loop over vertices */ 45 PetscCall(DMNetworkGetVertexRange(dm, &Start, &End)); 46 for (e = Start; e < End; e++) { 47 PetscCall(DMNetworkGetComponent(dm, e, ALL_COMPONENTS, NULL, NULL, &nvar)); 48 if (!nvar) continue; 49 i++; 50 } 51 (*ctx)->nlg = i; 52 PetscCall(PetscMalloc1(i, &(*ctx)->lg)); 53 54 i = 0; 55 /* loop over edges creating all needed line graphs*/ 56 PetscCall(DMNetworkGetEdgeRange(dm, &Start, &End)); 57 for (e = Start; e < End; e++) { 58 PetscCall(DMNetworkGetComponent(dm, e, ALL_COMPONENTS, NULL, NULL, &nvar)); 59 if (!nvar) continue; 60 PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw)); 61 PetscCall(PetscDrawSetFromOptions(draw)); 62 PetscCall(PetscDrawLGCreate(draw, nvar, &(*ctx)->lg[i])); 63 PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg[i])); 64 PetscCall(PetscDrawDestroy(&draw)); 65 i++; 66 } 67 /* loop over vertices */ 68 PetscCall(DMNetworkGetVertexRange(dm, &Start, &End)); 69 for (e = Start; e < End; e++) { 70 PetscCall(DMNetworkGetComponent(dm, e, ALL_COMPONENTS, NULL, NULL, &nvar)); 71 if (!nvar) continue; 72 PetscCall(PetscDrawCreate(comm, host, label, x, y, m, n, &draw)); 73 PetscCall(PetscDrawSetFromOptions(draw)); 74 PetscCall(PetscDrawLGCreate(draw, nvar, &(*ctx)->lg[i])); 75 PetscCall(PetscDrawLGSetFromOptions((*ctx)->lg[i])); 76 PetscCall(PetscDrawDestroy(&draw)); 77 i++; 78 } 79 PetscCall(PetscDrawDestroy(&draw)); 80 (*ctx)->howoften = howoften; 81 PetscFunctionReturn(PETSC_SUCCESS); 82 } 83 84 /* 85 TSMonitorLGCtxNetworkSolution - Monitors progress of the `TS` solvers for a `DMNETWORK` solution with one window for each vertex and each edge 86 87 Collective 88 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()` 95 96 Options Database Key: 97 . -ts_monitor_lg_solution_variables 98 99 Level: intermediate 100 101 Note: 102 Each process in a parallel run displays its component solutions in a separate window 103 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; 114 115 PetscFunctionBegin; 116 if (step < 0) PetscFunctionReturn(PETSC_SUCCESS); /* -1 indicates interpolated solution */ 117 if (!step) { 118 PetscDrawAxis axis; 119 120 for (i = 0; i < ctx->nlg; i++) { 121 PetscCall(PetscDrawLGGetAxis(ctx->lg[i], &axis)); 122 PetscCall(PetscDrawAxisSetLabels(axis, "Solution as function of time", "Time", "Solution")); 123 PetscCall(PetscDrawLGReset(ctx->lg[i])); 124 } 125 } 126 127 if (ctx->semilogy) { 128 PetscInt n, j; 129 130 PetscCall(VecDuplicate(u, &uv)); 131 PetscCall(VecCopy(u, uv)); 132 PetscCall(VecGetArray(uv, &yv)); 133 PetscCall(VecGetLocalSize(uv, &n)); 134 for (j = 0; j < n; j++) { 135 if (PetscRealPart(yv[j]) <= 0.0) yv[j] = -12; 136 else yv[j] = PetscLog10Real(PetscRealPart(yv[j])); 137 } 138 xv = yv; 139 } else { 140 PetscCall(VecGetArrayRead(u, &xv)); 141 } 142 /* iterate over edges */ 143 PetscCall(TSGetDM(ts, &dm)); 144 i = 0; 145 PetscCall(DMNetworkGetEdgeRange(dm, &Start, &End)); 146 for (e = Start; e < End; e++) { 147 PetscCall(DMNetworkGetComponent(dm, e, ALL_COMPONENTS, NULL, NULL, &nvar)); 148 if (!nvar) continue; 149 150 PetscCall(DMNetworkGetLocalVecOffset(dm, e, ALL_COMPONENTS, &offset)); 151 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg[i], ptime, (const PetscReal *)(xv + offset))); 152 i++; 153 } 154 155 /* iterate over vertices */ 156 PetscCall(DMNetworkGetVertexRange(dm, &Start, &End)); 157 for (v = Start; v < End; v++) { 158 PetscCall(DMNetworkGetComponent(dm, v, ALL_COMPONENTS, NULL, NULL, &nvar)); 159 if (!nvar) continue; 160 161 PetscCall(DMNetworkGetLocalVecOffset(dm, v, ALL_COMPONENTS, &offset)); 162 PetscCall(PetscDrawLGAddCommonPoint(ctx->lg[i], ptime, (const PetscReal *)(xv + offset))); 163 i++; 164 } 165 if (ctx->semilogy) { 166 PetscCall(VecRestoreArray(uv, &yv)); 167 PetscCall(VecDestroy(&uv)); 168 } else { 169 PetscCall(VecRestoreArrayRead(u, &xv)); 170 } 171 172 PetscCall(TSGetConvergedReason(ts, &reason)); 173 if (((ctx->howoften > 0) && (!(step % ctx->howoften))) || ((ctx->howoften == -1) && reason)) { 174 for (i = 0; i < ctx->nlg; i++) { 175 PetscCall(PetscDrawLGDraw(ctx->lg[i])); 176 PetscCall(PetscDrawLGSave(ctx->lg[i])); 177 } 178 } 179 PetscFunctionReturn(PETSC_SUCCESS); 180 } 181