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) {
 53:     DMNetworkMonitorPop(*monitor);
 54:   }

 56:   PetscFree(*monitor);
 57:   return 0;
 58: }

 60: /*@
 61:   DMNetworkMonitorPop - Removes the most recently added viewer

 63:   Collective on monitor

 65:   Input Parameters:
 66: . monitor - the monitor

 68:   Level: intermediate

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

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

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

 89: /*@C
 90:   DMNetworkMonitorAdd - Adds a new viewer to monitor

 92:   Collective on monitor

 94:   Input Parameters:
 95: + monitor - the monitor
 96: . name - name of viewer
 97: . element - vertex / edge number
 98: . nodes - number of nodes
 99: . start - variable starting offset
100: . blocksize - variable blocksize
101: . xmin - xmin (or PETSC_DECIDE) for viewer
102: . xmax - xmax (or PETSC_DECIDE) for viewer
103: . ymin - ymin for viewer
104: . ymax - ymax for viewer
105: - hold - determines if plot limits should be held

107:   Level: intermediate

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

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

116: .seealso: DMNetworkMonitorCreate(), DMNetworkMonitorDestroy()
117: @*/
118: 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)
119: {
120:   PetscDrawLG          drawlg;
121:   PetscDrawAxis        axis;
122:   PetscMPIInt          rank, size;
123:   DMNetworkMonitorList node;
124:   char                 titleBuffer[64];
125:   PetscInt             vStart,vEnd,eStart,eEnd;

127:   MPI_Comm_rank(monitor->comm, &rank);
128:   MPI_Comm_size(monitor->comm, &size);

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

133:   /* Make window title */
134:   if (vStart <= element && element < vEnd) {
135:     PetscSNPrintf(titleBuffer, 64, "%s @ vertex %d [%d / %d]", name, element - vStart, rank, size-1);
136:   } else if (eStart <= element && element < eEnd) {
137:     PetscSNPrintf(titleBuffer, 64, "%s @ edge %d [%d / %d]", name, element - eStart, rank, size-1);
138:   } else {
139:     /* vertex / edge is not on local machine, so skip! */
140:     return 0;
141:   }

143:   PetscMalloc1(1, &node);

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

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

160:   node->element   = element;
161:   node->nodes     = nodes;
162:   node->start     = start;
163:   node->blocksize = blocksize;

165:   node->next         = monitor->firstnode;
166:   monitor->firstnode = node;
167:   return 0;
168: }

170: /*@
171:   DMNetworkMonitorView - Monitor function for TSMonitorSet.

173:   Collectiveon DMNetworkMonitor

175:   Input Parameters:
176: + monitor - DMNetworkMonitor object
177: - x - TS solution vector

179:   Level: intermediate

181: .seealso: DMNetworkMonitorCreate(), DMNetworkMonitorDestroy(), DMNetworkMonitorAdd()
182: @*/

184: PetscErrorCode DMNetworkMonitorView(DMNetworkMonitor monitor,Vec x)
185: {
186:   PetscInt            varoffset,i,start;
187:   const PetscScalar   *xx;
188:   PetscScalar         *vv;
189:   DMNetworkMonitorList node;

191:   VecGetArrayRead(x, &xx);
192:   for (node = monitor->firstnode; node; node = node->next) {
193:     DMNetworkGetGlobalVecOffset(monitor->network, node->element, ALL_COMPONENTS, &varoffset);
194:     VecGetArray(node->v, &vv);
195:     start = varoffset + node->start;
196:     for (i = 0; i < node->nodes; i++) {
197:       vv[i] = xx[start+i*node->blocksize];
198:     }
199:     VecRestoreArray(node->v, &vv);
200:     VecView(node->v, node->viewer);
201:   }
202:   VecRestoreArrayRead(x, &xx);
203:   return 0;
204: }