Actual source code: ex3.c
1: static char help[] = "This example tests subnetwork coupling. \n\
2: \n\n";
4: /* T
5: Concepts: DMNetwork
6: */
7: #include <petscdmnetwork.h>
9: typedef struct{
10: PetscInt id;
11: } Comp0;
13: typedef struct{
14: PetscScalar val;
15: } Comp1;
17: int main(int argc,char ** argv)
18: {
19: PetscMPIInt size,rank;
20: DM dmnetwork;
21: PetscInt i,j,net,Nsubnet,nsubnet,ne,nv,nvar,v,ncomp,compkey0,compkey1,compkey,goffset,row;
22: PetscInt numEdges[10],*edgelist[10],asvtx,bsvtx;
23: const PetscInt *vtx,*edges;
24: PetscBool sharedv,ghost,distribute=PETSC_TRUE,test=PETSC_FALSE;
25: Vec X;
26: Comp0 comp0;
27: Comp1 comp1;
28: PetscScalar val;
30: PetscInitialize(&argc,&argv,(char*)0,help);
31: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
32: MPI_Comm_size(PETSC_COMM_WORLD,&size);
34: /* Create a network of subnetworks */
35: nsubnet = 1;
36: if (size == 1) nsubnet = 2;
38: /* Create a dmnetwork and register components */
39: DMNetworkCreate(PETSC_COMM_WORLD,&dmnetwork);
40: DMNetworkRegisterComponent(dmnetwork,"comp0",sizeof(Comp0),&compkey0);
41: DMNetworkRegisterComponent(dmnetwork,"comp1",sizeof(Comp1),&compkey1);
43: /* Set componnet values - intentionally take rank-dependent value for test */
44: comp0.id = rank;
45: comp1.val = 10.0*rank;
47: /* Set number of subnetworks, numbers of vertices and edges over each subnetwork */
48: DMNetworkSetNumSubNetworks(dmnetwork,nsubnet,PETSC_DECIDE);
49: DMNetworkGetNumSubNetworks(dmnetwork,NULL,&Nsubnet);
51: /* Input subnetworks; when size>1, process[i] creates subnetwork[i] */
52: for (i=0; i<Nsubnet; i++) numEdges[i] = 0;
53: for (i=0; i<Nsubnet; i++) {
54: if (i == 0 && (size == 1 || (rank == i && size >1))) {
55: numEdges[i] = 3;
56: PetscMalloc1(2*numEdges[i],&edgelist[i]);
57: edgelist[i][0] = 0; edgelist[i][1] = 2;
58: edgelist[i][2] = 2; edgelist[i][3] = 1;
59: edgelist[i][4] = 1; edgelist[i][5] = 3;
61: } else if (i == 1 && (size == 1 || (rank == i && size >1))) {
62: numEdges[i] = 3;
63: PetscMalloc1(2*numEdges[i],&edgelist[i]);
64: edgelist[i][0] = 0; edgelist[i][1] = 3;
65: edgelist[i][2] = 3; edgelist[i][3] = 2;
66: edgelist[i][4] = 2; edgelist[i][5] = 1;
68: } else if (i>1 && (size == 1 || (rank == i && size >1))) {
69: numEdges[i] = 3;
70: PetscMalloc1(2*numEdges[i],&edgelist[i]);
71: for (j=0; j< numEdges[i]; j++) {
72: edgelist[i][2*j] = j; edgelist[i][2*j+1] = j+1;
73: }
74: }
75: }
77: /* Add subnetworks */
78: for (i=0; i<Nsubnet; i++) {
79: PetscInt netNum = -1;
80: DMNetworkAddSubnetwork(dmnetwork,NULL,numEdges[i],edgelist[i],&netNum);
81: }
83: /* Add shared vertices -- all processes hold this info at current implementation */
84: asvtx = bsvtx = 0;
85: for (j=1; j<Nsubnet; j++) {
86: /* vertex subnet[0].0 shares with vertex subnet[j].0 */
87: DMNetworkAddSharedVertices(dmnetwork,0,j,1,&asvtx,&bsvtx);
88: }
90: /* Setup the network layout */
91: DMNetworkLayoutSetUp(dmnetwork);
93: /* Get Subnetwork(); Add nvar=1 to subnet[0] and nvar=2 to other subnets */
94: for (net=0; net<Nsubnet; net++) {
95: DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
96: for (v=0; v<nv; v++) {
97: DMNetworkIsSharedVertex(dmnetwork,vtx[v],&sharedv);
98: if (sharedv) continue;
100: if (!net) {
101: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
102: } else {
103: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
104: }
105: }
106: }
108: /* Add components and nvar to shared vertex -- owning and all ghost ranks must call DMNetworkAddComponent() */
109: DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
110: for (v=0; v<nv; v++) {
111: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
112: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
113: }
115: /* Enable runtime option of graph partition type -- must be called before DMSetUp() */
116: if (size > 1) {
117: DM plexdm;
118: PetscPartitioner part;
119: DMNetworkGetPlex(dmnetwork,&plexdm);
120: DMPlexGetPartitioner(plexdm, &part);
121: PetscPartitionerSetType(part,PETSCPARTITIONERSIMPLE);
122: PetscOptionsSetValue(NULL,"-dm_plex_csr_alg","mat"); /* for parmetis */
123: }
125: /* Setup dmnetwork */
126: DMSetUp(dmnetwork);
128: /* Redistribute the network layout; use '-distribute false' to skip */
129: PetscOptionsGetBool(NULL,NULL,"-distribute",&distribute,NULL);
130: if (distribute) {
131: DMNetworkDistribute(&dmnetwork,0);
132: DMView(dmnetwork,PETSC_VIEWER_STDOUT_WORLD);
133: }
135: /* Create a global vector */
136: DMCreateGlobalVector(dmnetwork,&X);
137: VecSet(X,0.0);
139: /* Set X values at shared vertex */
140: DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
141: for (v=0; v<nv; v++) {
142: DMNetworkIsGhostVertex(dmnetwork,vtx[v],&ghost);
143: if (ghost) continue;
145: /* only one process holds a non-ghost vertex */
146: DMNetworkGetComponent(dmnetwork,vtx[v],ALL_COMPONENTS,NULL,NULL,&nvar);
147: DMNetworkGetNumComponents(dmnetwork,vtx[v],&ncomp);
148: /* PetscPrintf(PETSC_COMM_SELF,"[%d] shared v %D: nvar %D, ncomp %D\n",rank,vtx[v],nvar,ncomp); */
149: for (j=0; j<ncomp; j++) {
150: DMNetworkGetComponent(dmnetwork,vtx[v],j,&compkey,NULL,&nvar);
151: DMNetworkGetGlobalVecOffset(dmnetwork,vtx[v],j,&goffset);
152: for (i=0; i<nvar; i++) {
153: row = goffset + i;
154: val = compkey + 1.0;
155: VecSetValues(X,1,&row,&val,INSERT_VALUES);
156: }
157: }
158: }
159: VecAssemblyBegin(X);
160: VecAssemblyEnd(X);
161: VecView(X,PETSC_VIEWER_STDOUT_WORLD);
163: /* Test DMNetworkGetSubnetwork() */
164: PetscOptionsGetBool(NULL,NULL,"-test_getsubnet",&test,NULL);
165: if (test) {
166: net = 0;
167: PetscOptionsGetInt(NULL,NULL,"-subnet",&net,NULL);
168: DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
169: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] subnet %D: nv %D, ne %D\n",rank,net,nv,ne);
170: PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
171: MPI_Barrier(PETSC_COMM_WORLD);
173: for (i=0; i<nv; i++) {
174: DMNetworkIsGhostVertex(dmnetwork,vtx[i],&ghost);
175: DMNetworkIsSharedVertex(dmnetwork,vtx[i],&sharedv);
177: DMNetworkGetNumComponents(dmnetwork,vtx[i],&ncomp);
178: if (sharedv || ghost) {
179: PetscPrintf(PETSC_COMM_SELF," [%d] v %D is shared %d, is ghost %d, ncomp %D\n",rank,vtx[i],sharedv,ghost,ncomp);
180: }
182: for (j=0; j<ncomp; j++) {
183: void* component;
184: DMNetworkGetComponent(dmnetwork,vtx[i],j,&compkey,(void**)&component,NULL);
185: if (compkey == 0) {
186: Comp0 *mycomp0;
187: mycomp0 = (Comp0*)component;
188: PetscPrintf(PETSC_COMM_SELF," [%d] v %D compkey %D, mycomp0->id %D\n",rank,vtx[i],compkey,mycomp0->id);
189: } else if (compkey == 1) {
190: Comp1 *mycomp1;
191: mycomp1 = (Comp1*)component;
192: PetscPrintf(PETSC_COMM_SELF," [%d] v %D compkey %D, mycomp1->val %g\n",rank,vtx[i],compkey,mycomp1->val);
193: }
194: }
195: }
196: }
198: /* Free work space */
199: VecDestroy(&X);
200: for (i=0; i<Nsubnet; i++) {
201: if (size == 1 || rank == i) PetscFree(edgelist[i]);
202: }
204: DMDestroy(&dmnetwork);
205: PetscFinalize();
206: return 0;
207: }
209: /*TEST
211: build:
212: requires: !single double defined(PETSC_HAVE_ATTRIBUTEALIGNED)
214: test:
215: args:
217: test:
218: suffix: 2
219: nsize: 2
220: args: -options_left no
222: test:
223: suffix: 3
224: nsize: 4
225: args: -options_left no
227: TEST*/