1 static char help[] = "Test DMCreateCoordinateDM_Network, and related functions \n\n"; 2 3 #include <petscdmnetwork.h> 4 5 /* 6 CreateStarGraphEdgeList - Create a k-Star Graph Edgelist on current processor 7 Not Collective 8 9 Input Parameters: 10 . k - order of the star graph (number of edges) 11 . directin - if true direction of edges is towards the center vertex, otherwise they are directed out of the center vertex. 12 13 Output Parameters: 14 . ne - number of edges of this star graph 15 . edgelist - list of edges for this star graph, this is a one dimensional array with pairs of entries being the two vertices (in global numbering of the vertices) of each edge, 16 [first vertex of first edge, second vertex of first edge, first vertex of second edge, second vertex of second edge, etc]. 17 18 User is responsible for deallocating this memory. 19 */ 20 static PetscErrorCode StarGraphCreateEdgeList(PetscInt k, PetscBool directin, PetscInt *ne, PetscInt *edgelist[]) 21 { 22 PetscInt i; 23 24 PetscFunctionBegin; 25 *ne = k; 26 PetscCall(PetscCalloc1(2 * k, edgelist)); 27 28 if (directin) { 29 for (i = 0; i < k; i++) { 30 (*edgelist)[2 * i] = i + 1; 31 (*edgelist)[2 * i + 1] = 0; 32 } 33 } else { 34 for (i = 0; i < k; i++) { 35 (*edgelist)[2 * i] = 0; 36 (*edgelist)[2 * i + 1] = i + 1; 37 } 38 } 39 PetscFunctionReturn(PETSC_SUCCESS); 40 } 41 42 /* Simple Circular embedding of the star graph */ 43 static PetscErrorCode StarGraphSetCoordinates(DM dm, PetscReal *vcolor) 44 { 45 DM cdm; 46 Vec Coord; 47 PetscScalar *coord; 48 PetscInt vStart, vEnd, v, vglobal, compkey = 0, off, NVert; 49 PetscReal theta; 50 51 PetscFunctionBegin; 52 PetscCall(DMSetCoordinateDim(dm, 2)); 53 PetscCall(DMGetCoordinateDM(dm, &cdm)); 54 55 PetscCall(DMNetworkGetVertexRange(cdm, &vStart, &vEnd)); 56 PetscCall(DMNetworkRegisterComponent(cdm, "coordinate", sizeof(PetscReal), &compkey)); 57 for (v = vStart; v < vEnd; v++) { 58 PetscCall(DMNetworkGetGlobalVertexIndex(cdm, v, &vglobal)); 59 vcolor[v - vStart] = vglobal; 60 PetscCall(DMNetworkAddComponent(cdm, v, compkey, &vcolor[v - vStart], 2)); 61 } 62 PetscCall(DMNetworkFinalizeComponents(cdm)); 63 64 PetscCall(DMCreateLocalVector(cdm, &Coord)); 65 PetscCall(VecGetArray(Coord, &coord)); 66 PetscCall(DMNetworkGetNumVertices(cdm, NULL, &NVert)); 67 theta = 2 * PETSC_PI / (NVert - 1); 68 for (v = vStart; v < vEnd; v++) { 69 PetscCall(DMNetworkGetGlobalVertexIndex(cdm, v, &vglobal)); 70 PetscCall(DMNetworkGetLocalVecOffset(cdm, v, 0, &off)); 71 if (vglobal == 0) { 72 coord[off] = 0.0; 73 coord[off + 1] = 0.0; 74 } else { 75 /* embed on the unit circle */ 76 coord[off] = PetscCosReal(theta * (vglobal - 1)); 77 coord[off + 1] = PetscSinReal(theta * (vglobal - 1)); 78 } 79 } 80 PetscCall(VecRestoreArray(Coord, &coord)); 81 PetscCall(DMSetCoordinatesLocal(dm, Coord)); 82 PetscCall(VecDestroy(&Coord)); 83 PetscFunctionReturn(PETSC_SUCCESS); 84 } 85 86 /* 87 CreateSimpleStarGraph - Create a Distributed k-Star Graph DMNetwork with a single PetscInt component on 88 all edges and vertices, a selectable number of dofs on vertices and edges. Intended mostly to be used for testing purposes. 89 90 Input Parameters: 91 . comm - the communicator of the dm 92 . numdofvert - number of degrees of freedom (dofs) on vertices 93 . numdofedge - number of degrees of freedom (dofs) on edges 94 . k - order of the star graph (number of edges) 95 . directin - if true direction of edges is towards the center vertex, otherwise they are directed out of the center vertex 96 97 Output Parameters: 98 . newdm - The created and distributed simple Star Graph 99 */ 100 static PetscErrorCode StarGraphCreate(MPI_Comm comm, PetscInt numdofvert, PetscInt numdofedge, PetscInt k, PetscBool directin, DM *newdm) 101 { 102 DM dm; 103 PetscMPIInt rank; 104 PetscInt ne = 0, compkey, eStart, eEnd, vStart, vEnd, e, v; 105 PetscInt *edgelist = NULL, *compedge, *compvert; 106 PetscReal *vcolor; 107 108 PetscFunctionBegin; 109 PetscCall(DMNetworkCreate(comm, &dm)); 110 PetscCall(DMNetworkSetNumSubNetworks(dm, PETSC_DECIDE, 1)); 111 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 112 if (rank == 0) PetscCall(StarGraphCreateEdgeList(k, directin, &ne, &edgelist)); 113 PetscCall(DMNetworkAddSubnetwork(dm, "Main", ne, edgelist, NULL)); 114 PetscCall(DMNetworkRegisterComponent(dm, "dummy", sizeof(PetscInt), &compkey)); 115 PetscCall(DMNetworkLayoutSetUp(dm)); 116 PetscCall(PetscFree(edgelist)); 117 PetscCall(DMNetworkGetEdgeRange(dm, &eStart, &eEnd)); 118 PetscCall(DMNetworkGetVertexRange(dm, &vStart, &vEnd)); 119 PetscCall(PetscCalloc3(eEnd - eStart, &compedge, vEnd - vStart, &compvert, vEnd - vStart, &vcolor)); 120 for (e = eStart; e < eEnd; e++) { 121 compedge[e - eStart] = e; 122 PetscCall(DMNetworkAddComponent(dm, e, compkey, &compedge[e - eStart], numdofedge)); 123 } 124 for (v = vStart; v < vEnd; v++) { 125 compvert[v - vStart] = v; 126 PetscCall(DMNetworkAddComponent(dm, v, compkey, &compvert[v - vStart], numdofvert)); 127 } 128 129 PetscCall(StarGraphSetCoordinates(dm, vcolor)); 130 131 PetscCall(DMSetFromOptions(dm)); 132 PetscCall(DMSetUp(dm)); 133 PetscCall(PetscFree3(compedge, compvert, vcolor)); 134 135 *newdm = dm; 136 PetscFunctionReturn(PETSC_SUCCESS); 137 } 138 139 /* This subroutine is used in petsc/src/snes/tutorials/network/ex1.c */ 140 static PetscErrorCode CoordinatePrint(DM dm) 141 { 142 DM dmclone; 143 PetscInt cdim, v, off, vglobal, vStart, vEnd; 144 const PetscScalar *carray; 145 Vec coords; 146 MPI_Comm comm; 147 PetscMPIInt rank; 148 149 PetscFunctionBegin; 150 PetscCall(PetscObjectGetComm((PetscObject)dm, &comm)); 151 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 152 153 PetscCall(DMGetCoordinateDM(dm, &dmclone)); 154 PetscCall(DMNetworkGetVertexRange(dm, &vStart, &vEnd)); 155 PetscCall(DMGetCoordinatesLocal(dm, &coords)); 156 157 PetscCall(DMGetCoordinateDim(dm, &cdim)); 158 PetscCall(VecGetArrayRead(coords, &carray)); 159 160 PetscCall(PetscPrintf(MPI_COMM_WORLD, "\nCoordinatePrint, cdim %" PetscInt_FMT ":\n", cdim)); 161 PetscCall(PetscSynchronizedPrintf(MPI_COMM_WORLD, "[%d]\n", rank)); 162 for (v = vStart; v < vEnd; v++) { 163 PetscCall(DMNetworkGetLocalVecOffset(dmclone, v, 0, &off)); 164 PetscCall(DMNetworkGetGlobalVertexIndex(dmclone, v, &vglobal)); 165 switch (cdim) { 166 case 2: 167 PetscCall(PetscSynchronizedPrintf(MPI_COMM_WORLD, "Vertex: %" PetscInt_FMT ", x = %f y = %f \n", vglobal, (double)PetscRealPart(carray[off]), (double)PetscRealPart(carray[off + 1]))); 168 break; 169 default: 170 PetscCheck(cdim == 2, MPI_COMM_WORLD, PETSC_ERR_SUP, "Only supports Network embedding dimension of 2, not supplied %" PetscInt_FMT, cdim); 171 break; 172 } 173 } 174 PetscCall(PetscSynchronizedFlush(MPI_COMM_WORLD, NULL)); 175 PetscCall(VecRestoreArrayRead(coords, &carray)); 176 PetscFunctionReturn(PETSC_SUCCESS); 177 } 178 179 int main(int argc, char **argv) 180 { 181 DM dm; 182 PetscInt dofv = 1, dofe = 1, ne = 1; 183 PetscMPIInt rank; 184 PetscBool viewCSV = PETSC_FALSE; 185 186 PetscFunctionBeginUser; 187 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 188 PetscCallMPI(MPI_Comm_rank(MPI_COMM_WORLD, &rank)); 189 190 /* create a distributed k-Star graph DMNetwork */ 191 PetscCall(PetscOptionsGetInt(NULL, NULL, "-dofv", &dofv, NULL)); 192 PetscCall(PetscOptionsGetInt(NULL, NULL, "-dofe", &dofe, NULL)); 193 PetscCall(PetscOptionsGetInt(NULL, NULL, "-ne", &ne, NULL)); 194 PetscCall(PetscOptionsGetBool(NULL, NULL, "-viewCSV", &viewCSV, NULL)); 195 196 /* create and setup a quick R^2 embedding of the star graph */ 197 PetscCall(StarGraphCreate(PETSC_COMM_WORLD, dofv, dofe, ne, PETSC_TRUE, &dm)); 198 199 PetscCall(DMNetworkDistribute(&dm, 0)); 200 if (viewCSV) { /* CSV View of network with coordinates */ 201 PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_CSV)); 202 PetscCall(DMView(dm, PETSC_VIEWER_STDOUT_WORLD)); 203 PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD)); 204 } else { 205 PetscCall(DMView(dm, PETSC_VIEWER_STDOUT_WORLD)); 206 } 207 208 /* print or view the coordinates of each vertex */ 209 PetscCall(CoordinatePrint(dm)); 210 211 PetscCall(DMDestroy(&dm)); 212 PetscCall(PetscFinalize()); 213 } 214 215 /*TEST 216 217 test: 218 suffix: 0 219 args: -ne 4 220 221 test: 222 suffix: 1 223 nsize: 2 224 args: -ne 4 -petscpartitioner_type simple 225 226 TEST*/ 227