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