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