Actual source code: ex4.c
1: static char help[] = "This example tests subnetwork coupling with zero size components. \n\n";
3: #include <petscdmnetwork.h>
5: int main(int argc,char ** argv)
6: {
7: PetscMPIInt size,rank;
8: DM dmnetwork;
9: PetscInt i,j,net,Nsubnet,ne,nv,nvar,v,goffset,row,compkey0,compkey1,compkey;
10: PetscInt *numEdges,**edgelist,asvtx[2],bsvtx[2];
11: const PetscInt *vtx,*edges;
12: PetscBool ghost,distribute=PETSC_TRUE,sharedv;
13: Vec X;
14: PetscScalar val;
16: PetscInitialize(&argc,&argv,(char*)0,help);
17: MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
18: MPI_Comm_size(PETSC_COMM_WORLD,&size);
20: /* Create a network of subnetworks */
21: if (size == 1) Nsubnet = 2;
22: else Nsubnet = (PetscInt)size;
23: PetscCalloc2(Nsubnet,&numEdges,Nsubnet,&edgelist);
25: /* when size>1, process[i] creates subnetwork[i] */
26: for (i=0; i<Nsubnet; i++) {
27: if (i == 0 && (size == 1 || (rank == i && size >1))) {
28: numEdges[i] = 3;
29: PetscMalloc1(2*numEdges[i],&edgelist[i]);
30: edgelist[i][0] = 0; edgelist[i][1] = 1;
31: edgelist[i][2] = 1; edgelist[i][3] = 2;
32: edgelist[i][4] = 2; edgelist[i][5] = 3;
34: } else if (i == 1 && (size == 1 || (rank == i && size >1))) {
35: numEdges[i] = 3;
36: PetscMalloc1(2*numEdges[i],&edgelist[i]);
37: edgelist[i][0] = 0; edgelist[i][1] = 1;
38: edgelist[i][2] = 1; edgelist[i][3] = 2;
39: edgelist[i][4] = 2; edgelist[i][5] = 3;
41: } else if (i>1 && (size == 1 || (rank == i && size >1))) {
42: numEdges[i] = 3;
43: PetscMalloc1(2*numEdges[i],&edgelist[i]);
44: for (j=0; j< numEdges[i]; j++) {
45: edgelist[i][2*j] = j; edgelist[i][2*j+1] = j+1;
46: }
47: }
48: }
50: /* Create a dmnetwork */
51: DMNetworkCreate(PETSC_COMM_WORLD,&dmnetwork);
53: /* Register zero size componets to get compkeys to be used by DMNetworkAddComponent() */
54: DMNetworkRegisterComponent(dmnetwork,"comp0",0,&compkey0);
55: DMNetworkRegisterComponent(dmnetwork,"comp1",0,&compkey1);
57: /* Set number of subnetworks, numbers of vertices and edges over each subnetwork */
58: DMNetworkSetNumSubNetworks(dmnetwork,PETSC_DECIDE,Nsubnet);
60: for (i=0; i<Nsubnet; i++) {
61: PetscInt netNum = -1;
62: DMNetworkAddSubnetwork(dmnetwork,NULL,numEdges[i],edgelist[i],&netNum);
63: }
65: /* Add shared vertices -- all processes hold this info at current implementation
66: net[0].0 -> net[j].0, j=0,...,Nsubnet-1
67: net[0].1 -> net[j].1, j=0,...,Nsubnet-1 */
68: asvtx[0] = bsvtx[0] = 0;
69: asvtx[1] = bsvtx[1] = 1;
70: for (j=Nsubnet-1; j>=1; j--) {
71: DMNetworkAddSharedVertices(dmnetwork,0,j,2,asvtx,bsvtx);
72: }
74: /* Setup the network layout */
75: DMNetworkLayoutSetUp(dmnetwork);
77: /* Get Subnetwork(); Add nvar=1 to subnet[0] and nvar=2 to other subnets */
78: for (net=0; net<Nsubnet; net++) {
79: DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
80: for (v=0; v<nv; v++) {
81: DMNetworkIsSharedVertex(dmnetwork,vtx[v],&sharedv);
82: if (sharedv) continue;
84: if (!net) {
85: /* Set nvar = 2 for subnet0 */
86: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,NULL,2);
87: } else {
88: /* Set nvar = 1 for other subnets */
89: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,NULL,1);
90: }
91: }
92: }
94: /* Add nvar to shared vertex -- owning and all ghost ranks must call DMNetworkAddComponent() */
95: DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
96: for (v=0; v<nv; v++) {
97: DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,NULL,2);
98: DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,NULL,1);
99: }
101: /* Enable runtime option of graph partition type -- must be called before DMSetUp() */
102: if (size > 1) {
103: DM plexdm;
104: PetscPartitioner part;
105: DMNetworkGetPlex(dmnetwork,&plexdm);
106: DMPlexGetPartitioner(plexdm, &part);
107: PetscPartitionerSetType(part,PETSCPARTITIONERSIMPLE);
108: PetscOptionsSetValue(NULL,"-dm_plex_csr_alg","mat"); /* for parmetis */
109: }
111: /* Setup dmnetwork */
112: DMSetUp(dmnetwork);
114: /* Redistribute the network layout; use '-distribute false' to skip */
115: PetscOptionsGetBool(NULL,NULL,"-distribute",&distribute,NULL);
116: if (distribute) {
117: DMNetworkDistribute(&dmnetwork,0);
118: DMView(dmnetwork,PETSC_VIEWER_STDOUT_WORLD);
119: }
121: /* Create a global vector */
122: DMCreateGlobalVector(dmnetwork,&X);
123: VecSet(X,0.0);
125: /* Set X values at shared vertex */
126: DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
127: for (v=0; v<nv; v++) {
128: DMNetworkIsGhostVertex(dmnetwork,vtx[v],&ghost);
129: if (ghost) continue;
131: /* only one process holds a non-ghost vertex */
132: DMNetworkGetComponent(dmnetwork,vtx[v],ALL_COMPONENTS,NULL,NULL,&nvar);
133: DMNetworkGetGlobalVecOffset(dmnetwork,vtx[v],ALL_COMPONENTS,&goffset);
134: for (i=0; i<nvar; i++) {
135: row = goffset + i;
136: val = (PetscScalar)rank + 1.0;
137: VecSetValues(X,1,&row,&val,ADD_VALUES);
138: }
140: DMNetworkGetComponent(dmnetwork,vtx[v],1,&compkey,NULL,&nvar);
141: DMNetworkGetGlobalVecOffset(dmnetwork,vtx[v],compkey,&goffset);
142: for (i=0; i<nvar; i++) {
143: row = goffset + i;
144: val = 1.0;
145: VecSetValues(X,1,&row,&val,ADD_VALUES);
146: }
147: }
148: VecAssemblyBegin(X);
149: VecAssemblyEnd(X);
150: VecView(X,PETSC_VIEWER_STDOUT_WORLD);
152: /* Free work space */
153: VecDestroy(&X);
154: for (i=0; i<Nsubnet; i++) {
155: if (size == 1 || rank == i) PetscFree(edgelist[i]);
156: }
157: PetscFree2(numEdges,edgelist);
158: DMDestroy(&dmnetwork);
159: PetscFinalize();
160: return 0;
161: }
163: /*TEST
165: build:
166: requires: !single double defined(PETSC_HAVE_ATTRIBUTEALIGNED)
168: test:
169: args:
171: test:
172: suffix: 2
173: nsize: 2
174: args: -options_left no -petscpartitioner_type simple
176: test:
177: suffix: 3
178: nsize: 4
179: args: -options_left no -petscpartitioner_type simple
181: TEST*/