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