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: {
20: PetscMPIInt size,rank;
21: DM dmnetwork;
22: PetscInt i,j,net,Nsubnet,nsubnet,ne,nv,nvar,v,ncomp,compkey0,compkey1,compkey,goffset,row;
23: PetscInt numEdges[10],*edgelist[10],asvtx,bsvtx;
24: const PetscInt *vtx,*edges;
25: PetscBool sharedv,ghost,distribute=PETSC_TRUE,test=PETSC_FALSE;
26: Vec X;
27: Comp0 comp0;
28: Comp1 comp1;
29: PetscScalar val;
31: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
32: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
33: MPI_Comm_size(PETSC_COMM_WORLD,&size);
35: /* Create a network of subnetworks */
36: nsubnet = 1;
37: if (size == 1) nsubnet = 2;
39: /* Create a dmnetwork and register components */
40: DMNetworkCreate(PETSC_COMM_WORLD,&dmnetwork);
41: DMNetworkRegisterComponent(dmnetwork,"comp0",sizeof(Comp0),&compkey0);
42: DMNetworkRegisterComponent(dmnetwork,"comp1",sizeof(Comp1),&compkey1);
44: /* Set componnet values - intentionally take rank-dependent value for test */
45: comp0.id = rank;
46: comp1.val = 10.0*rank;
48: /* Set number of subnetworks, numbers of vertices and edges over each subnetwork */
49: DMNetworkSetNumSubNetworks(dmnetwork,nsubnet,PETSC_DECIDE);
50: DMNetworkGetNumSubNetworks(dmnetwork,NULL,&Nsubnet);
52: /* Input subnetworks; when size>1, process[i] creates subnetwork[i] */
53: for (i=0; i<Nsubnet; i++) numEdges[i] = 0;
54: for (i=0; i<Nsubnet; i++) {
55: if (i == 0 && (size == 1 || (rank == i && size >1))) {
56: numEdges[i] = 3;
57: PetscMalloc1(2*numEdges[i],&edgelist[i]);
58: edgelist[i][0] = 0; edgelist[i][1] = 2;
59: edgelist[i][2] = 2; edgelist[i][3] = 1;
60: edgelist[i][4] = 1; edgelist[i][5] = 3;
62: } else if (i == 1 && (size == 1 || (rank == i && size >1))) {
63: numEdges[i] = 3;
64: PetscMalloc1(2*numEdges[i],&edgelist[i]);
65: edgelist[i][0] = 0; edgelist[i][1] = 3;
66: edgelist[i][2] = 3; edgelist[i][3] = 2;
67: edgelist[i][4] = 2; edgelist[i][5] = 1;
69: } else if (i>1 && (size == 1 || (rank == i && size >1))) {
70: numEdges[i] = 3;
71: PetscMalloc1(2*numEdges[i],&edgelist[i]);
72: for (j=0; j< numEdges[i]; j++) {
73: edgelist[i][2*j] = j; edgelist[i][2*j+1] = j+1;
74: }
75: }
76: }
78: /* Add subnetworks */
79: for (i=0; i<Nsubnet; i++) {
80: PetscInt netNum = -1;
81: DMNetworkAddSubnetwork(dmnetwork,NULL,numEdges[i],edgelist[i],&netNum);
82: }
84: /* Add shared vertices -- all processes hold this info at current implementation */
85: asvtx = bsvtx = 0;
86: for (j=1; j<Nsubnet; j++) {
87: /* vertex subnet[0].0 shares with vertex subnet[j].0 */
88: DMNetworkAddSharedVertices(dmnetwork,0,j,1,&asvtx,&bsvtx);
89: }
91: /* Setup the network layout */
92: DMNetworkLayoutSetUp(dmnetwork);
94: /* Get Subnetwork(); Add nvar=1 to subnet[0] and nvar=2 to other subnets */
95: for (net=0; net<Nsubnet; net++) {
96: DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
97: for (v=0; v<nv; v++) {
98: DMNetworkIsSharedVertex(dmnetwork,vtx[v],&sharedv);
99: if (sharedv) {
100: #if 0
101: /* current version requirs all processess add componenets and nvar at the shared vertices! */
102: if (net == 0) {
103: //printf("[%d] net %d, v[%d]=%d is a shared vertex, add comp0[0]\n",rank,net,v,vtx[v]);
104: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
105: } else if (net == 1) {
106: //printf("[%d] net %d, v[%d]=%d is a shared vertex, add comp1\n",rank,net,v,vtx[v]);
107: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
108: }
109: #endif
110: continue;
111: }
113: if (!net) {
114: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
115: } else {
116: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
117: }
118: }
119: }
121: /* Add components and nvar to shared vertex -- only owner of the vertex does this! */
122: DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
123: for (v=0; v<nv; v++) {
124: DMNetworkIsGhostVertex(dmnetwork,vtx[v],&ghost);
125: if (ghost) continue;
126: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
127: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
128: }
130: /* Enable runtime option of graph partition type -- must be called before DMSetUp() */
131: if (size > 1) {
132: DM plexdm;
133: PetscPartitioner part;
134: DMNetworkGetPlex(dmnetwork,&plexdm);
135: DMPlexGetPartitioner(plexdm, &part);
136: PetscPartitionerSetType(part,PETSCPARTITIONERSIMPLE);
137: PetscOptionsSetValue(NULL,"-dm_plex_csr_via_mat","true"); /* for parmetis */
138: }
140: /* Setup dmnetwork */
141: DMSetUp(dmnetwork);
143: /* Redistribute the network layout; use '-distribute false' to skip */
144: PetscOptionsGetBool(NULL,NULL,"-distribute",&distribute,NULL);
145: if (distribute) {
146: DMNetworkDistribute(&dmnetwork,0);
147: DMView(dmnetwork,PETSC_VIEWER_STDOUT_WORLD);
148: }
150: /* Create a global vector */
151: DMCreateGlobalVector(dmnetwork,&X);
152: VecSet(X,0.0);
154: /* Set X values at shared vertex */
155: DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
156: for (v=0; v<nv; v++) {
157: DMNetworkIsGhostVertex(dmnetwork,vtx[v],&ghost);
158: if (ghost) continue;
160: /* only one process holds a non-ghost vertex */
161: DMNetworkGetComponent(dmnetwork,vtx[v],ALL_COMPONENTS,NULL,NULL,&nvar);
162: DMNetworkGetNumComponents(dmnetwork,vtx[v],&ncomp);
163: /* PetscPrintf(PETSC_COMM_SELF,"[%d] shared v %D: nvar %D, ncomp %D\n",rank,vtx[v],nvar,ncomp); */
164: for (j=0; j<ncomp; j++) {
165: DMNetworkGetComponent(dmnetwork,vtx[v],j,&compkey,NULL,&nvar);
166: DMNetworkGetGlobalVecOffset(dmnetwork,vtx[v],j,&goffset);
167: for (i=0; i<nvar; i++) {
168: row = goffset + i;
169: val = compkey + 1.0;
170: VecSetValues(X,1,&row,&val,INSERT_VALUES);
171: }
172: }
173: }
174: VecAssemblyBegin(X);
175: VecAssemblyEnd(X);
176: VecView(X,PETSC_VIEWER_STDOUT_WORLD);
178: /* Test DMNetworkGetSubnetwork() */
179: PetscOptionsGetBool(NULL,NULL,"-test_getsubnet",&test,NULL);
180: if (test) {
181: net = 0;
182: PetscOptionsGetInt(NULL,NULL,"-subnet",&net,NULL);
183: DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
184: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] subnet %D: nv %D, ne %D\n",rank,net,nv,ne);
185: PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
186: MPI_Barrier(PETSC_COMM_WORLD);
188: for (i=0; i<nv; i++) {
189: DMNetworkIsGhostVertex(dmnetwork,vtx[i],&ghost);
190: DMNetworkIsSharedVertex(dmnetwork,vtx[i],&sharedv);
192: DMNetworkGetNumComponents(dmnetwork,vtx[i],&ncomp);
193: if (sharedv || ghost) {
194: PetscPrintf(PETSC_COMM_SELF," [%d] v %D is shared %d, is ghost %d, ncomp %D\n",rank,vtx[i],sharedv,ghost,ncomp);
195: }
197: for (j=0; j<ncomp; j++) {
198: void* component;
199: DMNetworkGetComponent(dmnetwork,vtx[i],j,&compkey,(void**)&component,NULL);
200: if (compkey == 0) {
201: Comp0 *mycomp0;
202: mycomp0 = (Comp0*)component;
203: PetscPrintf(PETSC_COMM_SELF," [%d] v %D compkey %D, mycomp0->id %D\n",rank,vtx[i],compkey,mycomp0->id);
204: } else if (compkey == 1) {
205: Comp1 *mycomp1;
206: mycomp1 = (Comp1*)component;
207: PetscPrintf(PETSC_COMM_SELF," [%d] v %D compkey %D, mycomp1->val %g\n",rank,vtx[i],compkey,mycomp1->val);
208: }
209: }
210: }
211: }
213: /* Free work space */
214: VecDestroy(&X);
215: for (i=0; i<Nsubnet; i++) {
216: if (size == 1 || rank == i) {PetscFree(edgelist[i]);}
217: }
219: DMDestroy(&dmnetwork);
220: PetscFinalize();
221: return ierr;
222: }
224: /*TEST
226: build:
227: requires: !single double defined(PETSC_HAVE_ATTRIBUTEALIGNED)
229: test:
230: args:
232: test:
233: suffix: 2
234: nsize: 2
235: args: -options_left no
237: test:
238: suffix: 3
239: nsize: 4
240: args: -options_left no
242: TEST*/