1 #include <petscdmnetwork.h> /*I "petscdmnetwork.h" I*/ 2 #include <petscdraw.h> 3 4 /*@ 5 DMNetworkMonitorCreate - Creates a network monitor context 6 7 Collective 8 9 Input Parameters: 10 . network - network to monitor 11 12 Output Parameters: 13 . monitorptr - Location to put network monitor context 14 15 Level: intermediate 16 17 .seealso: `DMNetworkMonitorDestroy()`, `DMNetworkMonitorAdd()` 18 @*/ 19 PetscErrorCode DMNetworkMonitorCreate(DM network,DMNetworkMonitor *monitorptr) 20 { 21 DMNetworkMonitor monitor; 22 MPI_Comm comm; 23 PetscMPIInt size; 24 25 PetscFunctionBegin; 26 PetscCall(PetscObjectGetComm((PetscObject)network,&comm)); 27 PetscCallMPI(MPI_Comm_size(comm, &size)); 28 PetscCheck(size == 1,PETSC_COMM_SELF,PETSC_ERR_SUP,"Parallel DMNetworkMonitor is not supported yet"); 29 30 PetscCall(PetscMalloc1(1,&monitor)); 31 monitor->comm = comm; 32 monitor->network = network; 33 monitor->firstnode = NULL; 34 35 *monitorptr = monitor; 36 PetscFunctionReturn(0); 37 } 38 39 /*@ 40 DMNetworkMonitorDestroy - Destroys a network monitor and all associated viewers 41 42 Collective on monitor 43 44 Input Parameters: 45 . monitor - monitor to destroy 46 47 Level: intermediate 48 49 .seealso: `DMNetworkMonitorCreate`, `DMNetworkMonitorAdd` 50 @*/ 51 PetscErrorCode DMNetworkMonitorDestroy(DMNetworkMonitor *monitor) 52 { 53 PetscFunctionBegin; 54 while ((*monitor)->firstnode) { 55 PetscCall(DMNetworkMonitorPop(*monitor)); 56 } 57 58 PetscCall(PetscFree(*monitor)); 59 PetscFunctionReturn(0); 60 } 61 62 /*@ 63 DMNetworkMonitorPop - Removes the most recently added viewer 64 65 Collective on monitor 66 67 Input Parameters: 68 . monitor - the monitor 69 70 Level: intermediate 71 72 .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()` 73 @*/ 74 PetscErrorCode DMNetworkMonitorPop(DMNetworkMonitor monitor) 75 { 76 DMNetworkMonitorList node; 77 78 PetscFunctionBegin; 79 if (monitor->firstnode) { 80 /* Update links */ 81 node = monitor->firstnode; 82 monitor->firstnode = node->next; 83 84 /* Free list node */ 85 PetscCall(PetscViewerDestroy(&(node->viewer))); 86 PetscCall(VecDestroy(&(node->v))); 87 PetscCall(PetscFree(node)); 88 } 89 PetscFunctionReturn(0); 90 } 91 92 /*@C 93 DMNetworkMonitorAdd - Adds a new viewer to monitor 94 95 Collective on monitor 96 97 Input Parameters: 98 + monitor - the monitor 99 . name - name of viewer 100 . element - vertex / edge number 101 . nodes - number of nodes 102 . start - variable starting offset 103 . blocksize - variable blocksize 104 . xmin - xmin (or PETSC_DECIDE) for viewer 105 . xmax - xmax (or PETSC_DECIDE) for viewer 106 . ymin - ymin for viewer 107 . ymax - ymax for viewer 108 - hold - determines if plot limits should be held 109 110 Level: intermediate 111 112 Notes: 113 This is written to be independent of the semantics associated to the variables 114 at a given network vertex / edge. 115 116 Precisely, the parameters nodes, start and blocksize allow you to select a general 117 strided subarray of the variables to monitor. 118 119 .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()` 120 @*/ 121 PetscErrorCode DMNetworkMonitorAdd(DMNetworkMonitor monitor,const char *name,PetscInt element,PetscInt nodes,PetscInt start,PetscInt blocksize,PetscReal xmin,PetscReal xmax,PetscReal ymin,PetscReal ymax,PetscBool hold) 122 { 123 PetscDrawLG drawlg; 124 PetscDrawAxis axis; 125 PetscMPIInt rank, size; 126 DMNetworkMonitorList node; 127 char titleBuffer[64]; 128 PetscInt vStart,vEnd,eStart,eEnd; 129 130 PetscFunctionBegin; 131 PetscCallMPI(MPI_Comm_rank(monitor->comm, &rank)); 132 PetscCallMPI(MPI_Comm_size(monitor->comm, &size)); 133 134 PetscCall(DMNetworkGetVertexRange(monitor->network, &vStart, &vEnd)); 135 PetscCall(DMNetworkGetEdgeRange(monitor->network, &eStart, &eEnd)); 136 137 /* Make window title */ 138 if (vStart <= element && element < vEnd) { 139 PetscCall(PetscSNPrintf(titleBuffer, sizeof(titleBuffer), "%s @ vertex %" PetscInt_FMT " [%d / %d]", name, element - vStart, rank, size-1)); 140 } else if (eStart <= element && element < eEnd) { 141 PetscCall(PetscSNPrintf(titleBuffer, sizeof(titleBuffer), "%s @ edge %" PetscInt_FMT " [%d / %d]", name, element - eStart, rank, size-1)); 142 } else { 143 /* vertex / edge is not on local machine, so skip! */ 144 PetscFunctionReturn(0); 145 } 146 147 PetscCall(PetscMalloc1(1, &node)); 148 149 /* Setup viewer. */ 150 PetscCall(PetscViewerDrawOpen(monitor->comm, NULL, titleBuffer, PETSC_DECIDE, PETSC_DECIDE, PETSC_DRAW_QUARTER_SIZE, PETSC_DRAW_QUARTER_SIZE, &(node->viewer))); 151 PetscCall(PetscViewerPushFormat(node->viewer, PETSC_VIEWER_DRAW_LG_XRANGE)); 152 PetscCall(PetscViewerDrawGetDrawLG(node->viewer, 0, &drawlg)); 153 PetscCall(PetscDrawLGGetAxis(drawlg, &axis)); 154 if (xmin != PETSC_DECIDE && xmax != PETSC_DECIDE) { 155 PetscCall(PetscDrawAxisSetLimits(axis, xmin, xmax, ymin, ymax)); 156 } else { 157 PetscCall(PetscDrawAxisSetLimits(axis, 0, nodes-1, ymin, ymax)); 158 } 159 PetscCall(PetscDrawAxisSetHoldLimits(axis, hold)); 160 161 /* Setup vector storage for drawing. */ 162 PetscCall(VecCreateSeq(PETSC_COMM_SELF, nodes, &(node->v))); 163 164 node->element = element; 165 node->nodes = nodes; 166 node->start = start; 167 node->blocksize = blocksize; 168 169 node->next = monitor->firstnode; 170 monitor->firstnode = node; 171 PetscFunctionReturn(0); 172 } 173 174 /*@ 175 DMNetworkMonitorView - Monitor function for TSMonitorSet. 176 177 Collectiveon DMNetworkMonitor 178 179 Input Parameters: 180 + monitor - DMNetworkMonitor object 181 - x - TS solution vector 182 183 Level: intermediate 184 185 .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()`, `DMNetworkMonitorAdd()` 186 @*/ 187 188 PetscErrorCode DMNetworkMonitorView(DMNetworkMonitor monitor,Vec x) 189 { 190 PetscInt varoffset,i,start; 191 const PetscScalar *xx; 192 PetscScalar *vv; 193 DMNetworkMonitorList node; 194 195 PetscFunctionBegin; 196 PetscCall(VecGetArrayRead(x, &xx)); 197 for (node = monitor->firstnode; node; node = node->next) { 198 PetscCall(DMNetworkGetGlobalVecOffset(monitor->network, node->element, ALL_COMPONENTS, &varoffset)); 199 PetscCall(VecGetArray(node->v, &vv)); 200 start = varoffset + node->start; 201 for (i = 0; i < node->nodes; i++) { 202 vv[i] = xx[start+i*node->blocksize]; 203 } 204 PetscCall(VecRestoreArray(node->v, &vv)); 205 PetscCall(VecView(node->v, node->viewer)); 206 } 207 PetscCall(VecRestoreArrayRead(x, &xx)); 208 PetscFunctionReturn(0); 209 } 210