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) PetscCall(DMNetworkMonitorPop(*monitor)); 55 56 PetscCall(PetscFree(*monitor)); 57 PetscFunctionReturn(0); 58 } 59 60 /*@ 61 DMNetworkMonitorPop - Removes the most recently added viewer 62 63 Collective on monitor 64 65 Input Parameters: 66 . monitor - the monitor 67 68 Level: intermediate 69 70 .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()` 71 @*/ 72 PetscErrorCode DMNetworkMonitorPop(DMNetworkMonitor monitor) 73 { 74 DMNetworkMonitorList node; 75 76 PetscFunctionBegin; 77 if (monitor->firstnode) { 78 /* Update links */ 79 node = monitor->firstnode; 80 monitor->firstnode = node->next; 81 82 /* Free list node */ 83 PetscCall(PetscViewerDestroy(&(node->viewer))); 84 PetscCall(VecDestroy(&(node->v))); 85 PetscCall(PetscFree(node)); 86 } 87 PetscFunctionReturn(0); 88 } 89 90 /*@C 91 DMNetworkMonitorAdd - Adds a new viewer to monitor 92 93 Collective on monitor 94 95 Input Parameters: 96 + monitor - the monitor 97 . name - name of viewer 98 . element - vertex / edge number 99 . nodes - number of nodes 100 . start - variable starting offset 101 . blocksize - variable blocksize 102 . xmin - xmin (or PETSC_DECIDE) for viewer 103 . xmax - xmax (or PETSC_DECIDE) for viewer 104 . ymin - ymin for viewer 105 . ymax - ymax for viewer 106 - hold - determines if plot limits should be held 107 108 Level: intermediate 109 110 Notes: 111 This is written to be independent of the semantics associated to the variables 112 at a given network vertex / edge. 113 114 Precisely, the parameters nodes, start and blocksize allow you to select a general 115 strided subarray of the variables to monitor. 116 117 .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()` 118 @*/ 119 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) 120 { 121 PetscDrawLG drawlg; 122 PetscDrawAxis axis; 123 PetscMPIInt rank, size; 124 DMNetworkMonitorList node; 125 char titleBuffer[64]; 126 PetscInt vStart, vEnd, eStart, eEnd; 127 128 PetscFunctionBegin; 129 PetscCallMPI(MPI_Comm_rank(monitor->comm, &rank)); 130 PetscCallMPI(MPI_Comm_size(monitor->comm, &size)); 131 132 PetscCall(DMNetworkGetVertexRange(monitor->network, &vStart, &vEnd)); 133 PetscCall(DMNetworkGetEdgeRange(monitor->network, &eStart, &eEnd)); 134 135 /* Make window title */ 136 if (vStart <= element && element < vEnd) { 137 PetscCall(PetscSNPrintf(titleBuffer, sizeof(titleBuffer), "%s @ vertex %" PetscInt_FMT " [%d / %d]", name, element - vStart, rank, size - 1)); 138 } else if (eStart <= element && element < eEnd) { 139 PetscCall(PetscSNPrintf(titleBuffer, sizeof(titleBuffer), "%s @ edge %" PetscInt_FMT " [%d / %d]", name, element - eStart, rank, size - 1)); 140 } else { 141 /* vertex / edge is not on local machine, so skip! */ 142 PetscFunctionReturn(0); 143 } 144 145 PetscCall(PetscMalloc1(1, &node)); 146 147 /* Setup viewer. */ 148 PetscCall(PetscViewerDrawOpen(monitor->comm, NULL, titleBuffer, PETSC_DECIDE, PETSC_DECIDE, PETSC_DRAW_QUARTER_SIZE, PETSC_DRAW_QUARTER_SIZE, &(node->viewer))); 149 PetscCall(PetscViewerPushFormat(node->viewer, PETSC_VIEWER_DRAW_LG_XRANGE)); 150 PetscCall(PetscViewerDrawGetDrawLG(node->viewer, 0, &drawlg)); 151 PetscCall(PetscDrawLGGetAxis(drawlg, &axis)); 152 if (xmin != PETSC_DECIDE && xmax != PETSC_DECIDE) PetscCall(PetscDrawAxisSetLimits(axis, xmin, xmax, ymin, ymax)); 153 else PetscCall(PetscDrawAxisSetLimits(axis, 0, nodes - 1, ymin, ymax)); 154 PetscCall(PetscDrawAxisSetHoldLimits(axis, hold)); 155 156 /* Setup vector storage for drawing. */ 157 PetscCall(VecCreateSeq(PETSC_COMM_SELF, nodes, &(node->v))); 158 159 node->element = element; 160 node->nodes = nodes; 161 node->start = start; 162 node->blocksize = blocksize; 163 164 node->next = monitor->firstnode; 165 monitor->firstnode = node; 166 PetscFunctionReturn(0); 167 } 168 169 /*@ 170 DMNetworkMonitorView - Monitor function for TSMonitorSet. 171 172 Collectiveon DMNetworkMonitor 173 174 Input Parameters: 175 + monitor - DMNetworkMonitor object 176 - x - TS solution vector 177 178 Level: intermediate 179 180 .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()`, `DMNetworkMonitorAdd()` 181 @*/ 182 183 PetscErrorCode DMNetworkMonitorView(DMNetworkMonitor monitor, Vec x) 184 { 185 PetscInt varoffset, i, start; 186 const PetscScalar *xx; 187 PetscScalar *vv; 188 DMNetworkMonitorList node; 189 190 PetscFunctionBegin; 191 PetscCall(VecGetArrayRead(x, &xx)); 192 for (node = monitor->firstnode; node; node = node->next) { 193 PetscCall(DMNetworkGetGlobalVecOffset(monitor->network, node->element, ALL_COMPONENTS, &varoffset)); 194 PetscCall(VecGetArray(node->v, &vv)); 195 start = varoffset + node->start; 196 for (i = 0; i < node->nodes; i++) vv[i] = xx[start + i * node->blocksize]; 197 PetscCall(VecRestoreArray(node->v, &vv)); 198 PetscCall(VecView(node->v, node->viewer)); 199 } 200 PetscCall(VecRestoreArrayRead(x, &xx)); 201 PetscFunctionReturn(0); 202 } 203