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