Actual source code: networkview.c

  1: #include <petsc/private/dmnetworkimpl.h>

  3: static PetscErrorCode DMView_Network_CSV(DM dm, PetscViewer viewer)
  4: {
  5:   DM              dmcoords;
  6:   PetscInt        nsubnets, i, subnet, nvertices, nedges, vertex, edge;
  7:   PetscInt        vertexOffsets[2], globalEdgeVertices[2];
  8:   PetscScalar     vertexCoords[2];
  9:   const PetscInt *vertices, *edges, *edgeVertices;
 10:   Vec             allVertexCoords;
 11:   PetscMPIInt     rank;
 12:   MPI_Comm        comm;

 14:   PetscFunctionBegin;
 15:   // Get the network containing coordinate information
 16:   PetscCall(DMGetCoordinateDM(dm, &dmcoords));
 17:   // Get the coordinate vector for the network
 18:   PetscCall(DMGetCoordinatesLocal(dm, &allVertexCoords));
 19:   // Get the MPI communicator and this process' rank
 20:   PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
 21:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
 22:   // Start synchronized printing
 23:   PetscCall(PetscViewerASCIIPushSynchronized(viewer));

 25:   // Write the header
 26:   PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "Type,Rank,ID,X,Y,Z,Name,Color\n"));

 28:   // Iterate each subnetwork (Note: We need to get the global number of subnets apparently)
 29:   PetscCall(DMNetworkGetNumSubNetworks(dm, NULL, &nsubnets));
 30:   for (subnet = 0; subnet < nsubnets; subnet++) {
 31:     // Get the subnetwork's vertices and edges
 32:     PetscCall(DMNetworkGetSubnetwork(dm, subnet, &nvertices, &nedges, &vertices, &edges));

 34:     // Write out each vertex
 35:     for (i = 0; i < nvertices; i++) {
 36:       vertex = vertices[i];
 37:       // Get the offset into the coordinate vector for the vertex
 38:       PetscCall(DMNetworkGetLocalVecOffset(dmcoords, vertex, ALL_COMPONENTS, vertexOffsets));
 39:       vertexOffsets[1] = vertexOffsets[0] + 1;
 40:       // Remap vertex to the global value
 41:       PetscCall(DMNetworkGetGlobalVertexIndex(dm, vertex, &vertex));
 42:       // Get the vertex position from the coordinate vector
 43:       PetscCall(VecGetValues(allVertexCoords, 2, vertexOffsets, vertexCoords));

 45:       // TODO: Determine vertex color/name
 46:       PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "Node,%" PetscInt_FMT ",%" PetscInt_FMT ",%lf,%lf,0,%" PetscInt_FMT "\n", (PetscInt)rank, vertex, (double)PetscRealPart(vertexCoords[0]), (double)PetscRealPart(vertexCoords[1]), vertex));
 47:     }

 49:     // Write out each edge
 50:     for (i = 0; i < nedges; i++) {
 51:       edge = edges[i];
 52:       PetscCall(DMNetworkGetConnectedVertices(dm, edge, &edgeVertices));
 53:       PetscCall(DMNetworkGetGlobalVertexIndex(dm, edgeVertices[0], &globalEdgeVertices[0]));
 54:       PetscCall(DMNetworkGetGlobalVertexIndex(dm, edgeVertices[1], &globalEdgeVertices[1]));
 55:       PetscCall(DMNetworkGetGlobalEdgeIndex(dm, edge, &edge));

 57:       // TODO: Determine edge color/name
 58:       PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "Edge,%" PetscInt_FMT ",%" PetscInt_FMT ",%" PetscInt_FMT ",%" PetscInt_FMT ",0,%" PetscInt_FMT "\n", (PetscInt)rank, edge, globalEdgeVertices[0], globalEdgeVertices[1], edge));
 59:     }
 60:   }
 61:   // End synchronized printing
 62:   PetscCall(PetscViewerFlush(viewer));
 63:   PetscCall(PetscViewerASCIIPopSynchronized(viewer));
 64:   PetscFunctionReturn(PETSC_SUCCESS);
 65: }

 67: #include <petscdraw.h>
 68: static PetscErrorCode DMView_Network_Matplotlib(DM dm, PetscViewer viewer)
 69: {
 70:   PetscMPIInt rank, size, rank2;
 71:   MPI_Comm    comm;
 72:   char        filename[PETSC_MAX_PATH_LEN + 1], proccall[PETSC_MAX_PATH_LEN + 500], scriptFile[PETSC_MAX_PATH_LEN + 1], streamBuffer[256];
 73:   PetscViewer csvViewer;
 74:   FILE       *processFile = NULL;
 75:   PetscBool   isnull;
 76:   PetscDraw   draw;

 78:   PetscFunctionBegin;
 79:   // Deal with the PetscDraw we are given
 80:   PetscCall(PetscViewerDrawGetDraw(viewer, 1, &draw));
 81:   PetscCall(PetscDrawIsNull(draw, &isnull));
 82:   PetscCall(PetscDrawSetVisible(draw, PETSC_FALSE));

 84:   // Clear the file name buffer so all communicated bytes are well-defined
 85:   PetscCall(PetscMemzero(filename, sizeof(filename)));

 87:   // Get the MPI communicator and this process' rank
 88:   PetscCall(PetscObjectGetComm((PetscObject)dm, &comm));
 89:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
 90:   PetscCallMPI(MPI_Comm_size(comm, &size));

 92:   // Generate and broadcast the temporary file name from rank 0
 93:   if (rank == 0) {
 94: #if defined(PETSC_HAVE_TMPNAM_S)
 95:     // Acquire a temporary file to write to and open an ASCII/CSV viewer
 96:     PetscCheck(tmpnam_s(filename, sizeof(filename)) == 0, comm, PETSC_ERR_SYS, "Could not acquire temporary file");
 97: #elif defined(PETSC_HAVE_MKSTEMP) && __STDC_VERSION__ > 199901L
 98:     size_t numChars;
 99:     // Same thing, but for POSIX systems on which tmpnam is deprecated
100:     // Note: Configure may detect mkstemp but it will not be defined if compiling for C99, so check additional defines to see if we can use it
101:     PetscCall(PetscStrcpy(filename, "/tmp/"));
102:     // Mkstemp requires us to explicitly specify part of the path, but some systems may not like putting files in /tmp/ so have an option for it
103:     PetscCall(PetscOptionsGetString(NULL, NULL, "-dmnetwork_view_tmpdir", filename, sizeof(filename), NULL));
104:     // Make sure the filename ends with a '/'
105:     PetscCall(PetscStrlen(filename, &numChars));
106:     if (filename[numChars - 1] != '/') {
107:       filename[numChars]     = '/';
108:       filename[numChars + 1] = 0;
109:     }
110:     // Perform the actual temporary file creation
111:     PetscCall(PetscStrcat(filename, "XXXXXX"));
112:     PetscCheck(mkstemp(filename) != -1, comm, PETSC_ERR_SYS, "Could not acquire temporary file");
113: #else
114:     // Same thing, but for older C versions which don't have the safe form
115:     PetscCheck(tmpnam(filename) != NULL, comm, PETSC_ERR_SYS, "Could not acquire temporary file");
116: #endif
117:     // Broadcast the filename to all other MPI ranks
118:     for (rank2 = 1; rank2 < size; rank2++) PetscCallMPI(MPI_Send(filename, PETSC_MAX_PATH_LEN, MPI_BYTE, rank2, 0, comm));
119:   } else {
120:     // Receive the file name
121:     PetscCallMPI(MPI_Recv(filename, PETSC_MAX_PATH_LEN, MPI_BYTE, 0, 0, comm, MPI_STATUS_IGNORE));
122:   }

124:   PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, filename, &csvViewer));
125:   PetscCall(PetscViewerPushFormat(csvViewer, PETSC_VIEWER_ASCII_CSV));

127:   // Use the CSV viewer to write out the local network
128:   PetscCall(DMView_Network_CSV(dm, csvViewer));

130:   // Close the viewer
131:   PetscCall(PetscViewerDestroy(&csvViewer));

133:   // Get the value of $PETSC_DIR
134:   PetscCall(PetscStrreplace(PETSC_COMM_WORLD, "${PETSC_DIR}/share/petsc/bin/dmnetwork_view.py", scriptFile, sizeof(scriptFile)));
135:   PetscCall(PetscFixFilename(scriptFile, scriptFile));
136:   // Generate the system call for 'python3 $PETSC_DIR/share/petsc/dmnetwork_view.py '
137:   PetscCall(PetscArrayzero(proccall, sizeof(proccall)));
138:   PetscCall(PetscSNPrintf(proccall, sizeof(proccall), "%s %s %s %s", PETSC_PYTHON_EXE, scriptFile, (isnull ? "-tx" : ""), filename));

140: #if defined(PETSC_HAVE_POPEN)
141:   // Perform the call to run the python script (Note: while this is called on all ranks POpen will only run on rank 0)
142:   PetscCall(PetscPOpen(PETSC_COMM_WORLD, NULL, proccall, "r", &processFile));
143:   if (processFile != NULL) {
144:     while (fgets(streamBuffer, sizeof(streamBuffer), processFile) != NULL) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "%s", streamBuffer));
145:   }
146:   PetscCall(PetscPClose(PETSC_COMM_WORLD, processFile));
147: #else
148:   // Same thing, but using the standard library for systems that don't have POpen/PClose (only run on rank 0)
149:   if (rank == 0) {
150:     PetscCheck(system(proccall) == 0, comm, PETSC_ERR_SYS, "Failed to call viewer script");
151:     // Barrier so that all ranks wait until the call completes
152:     PetscCallMPI(MPI_Barrier(PETSC_COMM_WORLD));
153:   }
154: #endif
155:   // Clean up the temporary file we used using rank 0
156:   if (rank == 0) PetscCheck(remove(filename) == 0, comm, PETSC_ERR_SYS, "Failed to delete temporary file");
157:   PetscFunctionReturn(PETSC_SUCCESS);
158: }

160: PetscErrorCode DMView_Network(DM dm, PetscViewer viewer)
161: {
162:   PetscBool         iascii, isdraw;
163:   PetscViewerFormat format;

165:   PetscFunctionBegin;
168:   PetscCall(PetscViewerGetFormat(viewer, &format));

170:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERDRAW, &isdraw));
171:   if (isdraw) {
172:     PetscCall(DMView_Network_Matplotlib(dm, viewer));
173:     PetscFunctionReturn(PETSC_SUCCESS);
174:   }

176:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
177:   if (iascii) {
178:     const PetscInt *cone, *vtx, *edges;
179:     PetscInt        vfrom, vto, i, j, nv, ne, nsv, p, nsubnet;
180:     DM_Network     *network = (DM_Network *)dm->data;
181:     PetscMPIInt     rank;

183:     PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)dm), &rank));
184:     if (format == PETSC_VIEWER_ASCII_CSV) {
185:       PetscCall(DMView_Network_CSV(dm, viewer));
186:       PetscFunctionReturn(PETSC_SUCCESS);
187:     }

189:     nsubnet = network->cloneshared->Nsubnet; /* num of subnetworks */
190:     if (!rank) {
191:       PetscCall(PetscPrintf(PETSC_COMM_SELF, "  NSubnets: %" PetscInt_FMT "; NEdges: %" PetscInt_FMT "; NVertices: %" PetscInt_FMT "; NSharedVertices: %" PetscInt_FMT ".\n", nsubnet, network->cloneshared->NEdges, network->cloneshared->NVertices,
192:                             network->cloneshared->Nsvtx));
193:     }

195:     PetscCall(DMNetworkGetSharedVertices(dm, &nsv, NULL));
196:     PetscCall(PetscViewerASCIIPushSynchronized(viewer));
197:     PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "  [%d] nEdges: %" PetscInt_FMT "; nVertices: %" PetscInt_FMT "; nSharedVertices: %" PetscInt_FMT "\n", rank, network->cloneshared->nEdges, network->cloneshared->nVertices, nsv));

199:     for (i = 0; i < nsubnet; i++) {
200:       PetscCall(DMNetworkGetSubnetwork(dm, i, &nv, &ne, &vtx, &edges));
201:       if (ne) {
202:         PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "     Subnet %" PetscInt_FMT ": nEdges %" PetscInt_FMT ", nVertices(include shared vertices) %" PetscInt_FMT "\n", i, ne, nv));
203:         for (j = 0; j < ne; j++) {
204:           p = edges[j];
205:           PetscCall(DMNetworkGetConnectedVertices(dm, p, &cone));
206:           PetscCall(DMNetworkGetGlobalVertexIndex(dm, cone[0], &vfrom));
207:           PetscCall(DMNetworkGetGlobalVertexIndex(dm, cone[1], &vto));
208:           PetscCall(DMNetworkGetGlobalEdgeIndex(dm, edges[j], &p));
209:           PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "       edge %" PetscInt_FMT ": %" PetscInt_FMT " ----> %" PetscInt_FMT "\n", p, vfrom, vto));
210:         }
211:       }
212:     }

214:     /* Shared vertices */
215:     PetscCall(DMNetworkGetSharedVertices(dm, NULL, &vtx));
216:     if (nsv) {
217:       PetscInt        gidx;
218:       PetscBool       ghost;
219:       const PetscInt *sv = NULL;

221:       PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "     SharedVertices:\n"));
222:       for (i = 0; i < nsv; i++) {
223:         PetscCall(DMNetworkIsGhostVertex(dm, vtx[i], &ghost));
224:         if (ghost) continue;

226:         PetscCall(DMNetworkSharedVertexGetInfo(dm, vtx[i], &gidx, &nv, &sv));
227:         PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "       svtx %" PetscInt_FMT ": global index %" PetscInt_FMT ", subnet[%" PetscInt_FMT "].%" PetscInt_FMT " ---->\n", i, gidx, sv[0], sv[1]));
228:         for (j = 1; j < nv; j++) PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "                                           ----> subnet[%" PetscInt_FMT "].%" PetscInt_FMT "\n", sv[2 * j], sv[2 * j + 1]));
229:       }
230:     }
231:     PetscCall(PetscViewerFlush(viewer));
232:     PetscCall(PetscViewerASCIIPopSynchronized(viewer));
233:   } else PetscCheck(iascii, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Viewer type %s not yet supported for DMNetwork writing", ((PetscObject)viewer)->type_name);
234:   PetscFunctionReturn(PETSC_SUCCESS);
235: }