1 static char help[] = "Check if DMClone for DMNetwork Correctly Shallow Clones Topology Only \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 */
CreateStarGraphEdgeList(PetscInt k,PetscBool directin,PetscInt * ne,PetscInt * edgelist[])20 PetscErrorCode CreateStarGraphEdgeList(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, aselectable 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 */
CreateSimpleStarGraph(MPI_Comm comm,PetscInt numdofvert,PetscInt numdofedge,PetscInt k,PetscBool directin,DM * newdm)56 PetscErrorCode CreateSimpleStarGraph(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(CreateStarGraphEdgeList(k, directin, &ne, &edgelist));
68 PetscCall(DMNetworkAddSubnetwork(dm, "Main", ne, edgelist, NULL));
69 PetscCall(DMNetworkRegisterComponent(dm, "unused", 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(DMViewFromOptions(dm, NULL, "-dm_view"));
85 PetscCall(DMSetUp(dm));
86 PetscCall(PetscFree2(compedge, compvert));
87 PetscCall(DMNetworkDistribute(&dm, 0));
88 *newdm = dm;
89 PetscFunctionReturn(PETSC_SUCCESS);
90 }
91
main(int argc,char ** argv)92 int main(int argc, char **argv)
93 {
94 DM dm, dmclone, plex;
95 PetscInt e, eStart, eEnd, ndofs, ndofsprev;
96 PetscInt *compprev, *comp, compkey;
97 PetscInt dofv = 1, dofe = 1, ne = 1;
98 PetscSection sec;
99 Vec vec;
100
101 PetscFunctionBeginUser;
102 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
103 /* create a distributed k-Star graph DMNetwork */
104 PetscCall(PetscOptionsGetInt(NULL, NULL, "-dofv", &dofv, NULL));
105 PetscCall(PetscOptionsGetInt(NULL, NULL, "-dofe", &dofe, NULL));
106 PetscCall(PetscOptionsGetInt(NULL, NULL, "-ne", &ne, NULL));
107 PetscCall(CreateSimpleStarGraph(PETSC_COMM_WORLD, dofv, dofe, ne, PETSC_TRUE, &dm));
108 PetscCall(DMNetworkGetEdgeRange(dm, &eStart, &eEnd));
109
110 /* check if cloning changed any component */
111 if (eStart < eEnd) PetscCall(DMNetworkGetComponent(dm, eStart, 0, NULL, (void **)&compprev, &ndofsprev));
112 PetscCall(DMClone(dm, &dmclone));
113 if (eStart < eEnd) {
114 PetscCall(DMNetworkGetComponent(dm, eStart, 0, NULL, (void **)&comp, &ndofs));
115 PetscCheck(*comp == *compprev, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Cloning changed the Original, comp (previous) : %" PetscInt_FMT " comp (now) : %" PetscInt_FMT, *compprev, *comp);
116 PetscCheck(ndofsprev == ndofs, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Cloning changed the Original, ndofs (previous) : %" PetscInt_FMT " ndofs (now) : %" PetscInt_FMT, ndofsprev, ndofs);
117 }
118
119 /* register new components to the clone and add a unused component to every point */
120 PetscCall(DMNetworkRegisterComponent(dmclone, "unusedclone", sizeof(PetscInt), &compkey));
121 PetscCall(DMNetworkGetEdgeRange(dmclone, &eStart, &eEnd));
122 PetscCall(PetscMalloc1(eEnd - eStart, &comp));
123 for (e = eStart; e < eEnd; e++) {
124 comp[e - eStart] = e;
125 PetscCall(DMNetworkAddComponent(dmclone, e, compkey, &comp[e - eStart], 2));
126 }
127 PetscCall(DMNetworkFinalizeComponents(dmclone));
128 PetscCall(PetscFree(comp));
129
130 PetscCall(PetscPrintf(PETSC_COMM_WORLD, " dm: \n"));
131 PetscCall(DMView(dm, PETSC_VIEWER_STDOUT_WORLD));
132 PetscCall(DMNetworkGetPlex(dm, &plex));
133 PetscCall(DMGetLocalSection(plex, &sec));
134 PetscCall(PetscSectionView(sec, PETSC_VIEWER_STDOUT_WORLD));
135
136 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n dmclone: \n"));
137 PetscCall(DMView(dmclone, PETSC_VIEWER_STDOUT_WORLD));
138 PetscCall(DMNetworkGetPlex(dmclone, &plex));
139 PetscCall(DMGetLocalSection(plex, &sec));
140 PetscCall(PetscSectionView(sec, PETSC_VIEWER_STDOUT_WORLD));
141
142 /* create Vectors */
143 PetscCall(DMCreateGlobalVector(dm, &vec));
144 PetscCall(VecSet(vec, 1.0));
145 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n dm vec:\n"));
146 PetscCall(VecView(vec, PETSC_VIEWER_STDOUT_WORLD));
147 PetscCall(VecDestroy(&vec));
148
149 PetscCall(DMCreateGlobalVector(dmclone, &vec));
150 PetscCall(VecSet(vec, 2.0));
151 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "\n dmclone vec:\n"));
152 PetscCall(VecView(vec, PETSC_VIEWER_STDOUT_WORLD));
153 PetscCall(VecDestroy(&vec));
154
155 PetscCall(DMDestroy(&dm));
156 PetscCall(DMDestroy(&dmclone));
157 PetscCall(PetscFinalize());
158 }
159
160 /*TEST
161
162 test:
163 suffix: 0
164 args:
165
166 test:
167 suffix: 1
168 nsize: 2
169 args: -dofv 2 -dofe 2 -ne 2
170
171 TEST*/
172