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: }