Actual source code: networkmonitor.c

  1: #include <petscdmnetwork.h>
  2: #include <petscdraw.h>

  4: /*@
  5:   DMNetworkMonitorCreate - Creates a network monitor context

  7:   Collective

  9:   Input Parameters:
 10: . network - network to monitor

 12:   Output Parameters:
 13: . monitorptr - Location to put network monitor context

 15:   Level: intermediate

 17: .seealso: `DMNetworkMonitorDestroy()`, `DMNetworkMonitorAdd()`
 18: @*/
 19: PetscErrorCode DMNetworkMonitorCreate(DM network, DMNetworkMonitor *monitorptr)
 20: {
 21:   DMNetworkMonitor monitor;
 22:   MPI_Comm         comm;
 23:   PetscMPIInt      size;

 25:   PetscObjectGetComm((PetscObject)network, &comm);
 26:   MPI_Comm_size(comm, &size);

 29:   PetscMalloc1(1, &monitor);
 30:   monitor->comm      = comm;
 31:   monitor->network   = network;
 32:   monitor->firstnode = NULL;

 34:   *monitorptr = monitor;
 35:   return 0;
 36: }

 38: /*@
 39:   DMNetworkMonitorDestroy - Destroys a network monitor and all associated viewers

 41:   Collective on monitor

 43:   Input Parameters:
 44: . monitor - monitor to destroy

 46:   Level: intermediate

 48: .seealso: `DMNetworkMonitorCreate`, `DMNetworkMonitorAdd`
 49: @*/
 50: PetscErrorCode DMNetworkMonitorDestroy(DMNetworkMonitor *monitor)
 51: {
 52:   while ((*monitor)->firstnode) DMNetworkMonitorPop(*monitor);

 54:   PetscFree(*monitor);
 55:   return 0;
 56: }

 58: /*@
 59:   DMNetworkMonitorPop - Removes the most recently added viewer

 61:   Collective on monitor

 63:   Input Parameters:
 64: . monitor - the monitor

 66:   Level: intermediate

 68: .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()`
 69: @*/
 70: PetscErrorCode DMNetworkMonitorPop(DMNetworkMonitor monitor)
 71: {
 72:   DMNetworkMonitorList node;

 74:   if (monitor->firstnode) {
 75:     /* Update links */
 76:     node               = monitor->firstnode;
 77:     monitor->firstnode = node->next;

 79:     /* Free list node */
 80:     PetscViewerDestroy(&(node->viewer));
 81:     VecDestroy(&(node->v));
 82:     PetscFree(node);
 83:   }
 84:   return 0;
 85: }

 87: /*@C
 88:   DMNetworkMonitorAdd - Adds a new viewer to monitor

 90:   Collective on monitor

 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

105:   Level: intermediate

107:   Notes:
108:   This is written to be independent of the semantics associated to the variables
109:   at a given network vertex / edge.

111:   Precisely, the parameters nodes, start and blocksize allow you to select a general
112:   strided subarray of the variables to monitor.

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: {
118:   PetscDrawLG          drawlg;
119:   PetscDrawAxis        axis;
120:   PetscMPIInt          rank, size;
121:   DMNetworkMonitorList node;
122:   char                 titleBuffer[64];
123:   PetscInt             vStart, vEnd, eStart, eEnd;

125:   MPI_Comm_rank(monitor->comm, &rank);
126:   MPI_Comm_size(monitor->comm, &size);

128:   DMNetworkGetVertexRange(monitor->network, &vStart, &vEnd);
129:   DMNetworkGetEdgeRange(monitor->network, &eStart, &eEnd);

131:   /* Make window title */
132:   if (vStart <= element && element < vEnd) {
133:     PetscSNPrintf(titleBuffer, sizeof(titleBuffer), "%s @ vertex %" PetscInt_FMT " [%d / %d]", name, element - vStart, rank, size - 1);
134:   } else if (eStart <= element && element < eEnd) {
135:     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:     return 0;
139:   }

141:   PetscMalloc1(1, &node);

143:   /* Setup viewer. */
144:   PetscViewerDrawOpen(monitor->comm, NULL, titleBuffer, PETSC_DECIDE, PETSC_DECIDE, PETSC_DRAW_QUARTER_SIZE, PETSC_DRAW_QUARTER_SIZE, &(node->viewer));
145:   PetscViewerPushFormat(node->viewer, PETSC_VIEWER_DRAW_LG_XRANGE);
146:   PetscViewerDrawGetDrawLG(node->viewer, 0, &drawlg);
147:   PetscDrawLGGetAxis(drawlg, &axis);
148:   if (xmin != PETSC_DECIDE && xmax != PETSC_DECIDE) PetscDrawAxisSetLimits(axis, xmin, xmax, ymin, ymax);
149:   else PetscDrawAxisSetLimits(axis, 0, nodes - 1, ymin, ymax);
150:   PetscDrawAxisSetHoldLimits(axis, hold);

152:   /* Setup vector storage for drawing. */
153:   VecCreateSeq(PETSC_COMM_SELF, nodes, &(node->v));

155:   node->element   = element;
156:   node->nodes     = nodes;
157:   node->start     = start;
158:   node->blocksize = blocksize;

160:   node->next         = monitor->firstnode;
161:   monitor->firstnode = node;
162:   return 0;
163: }

165: /*@
166:   DMNetworkMonitorView - Monitor function for TSMonitorSet.

168:   Collectiveon DMNetworkMonitor

170:   Input Parameters:
171: + monitor - DMNetworkMonitor object
172: - x - TS solution vector

174:   Level: intermediate

176: .seealso: `DMNetworkMonitorCreate()`, `DMNetworkMonitorDestroy()`, `DMNetworkMonitorAdd()`
177: @*/

179: PetscErrorCode DMNetworkMonitorView(DMNetworkMonitor monitor, Vec x)
180: {
181:   PetscInt             varoffset, i, start;
182:   const PetscScalar   *xx;
183:   PetscScalar         *vv;
184:   DMNetworkMonitorList node;

186:   VecGetArrayRead(x, &xx);
187:   for (node = monitor->firstnode; node; node = node->next) {
188:     DMNetworkGetGlobalVecOffset(monitor->network, node->element, ALL_COMPONENTS, &varoffset);
189:     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:     VecRestoreArray(node->v, &vv);
193:     VecView(node->v, node->viewer);
194:   }
195:   VecRestoreArrayRead(x, &xx);
196:   return 0;
197: }