Actual source code: ex3.c

  1: static char help[] = "This example tests subnetwork coupling. \n\
  2:               \n\n";

  4: #include <petscdmnetwork.h>

  6: typedef struct{
  7:   PetscInt id;
  8: } Comp0;

 10: typedef struct{
 11:   PetscScalar val;
 12: } Comp1;

 14: int main(int argc,char ** argv)
 15: {
 16:   PetscMPIInt    size,rank;
 17:   DM             dmnetwork;
 18:   PetscInt       i,j,net,Nsubnet,nsubnet,ne,nv,nvar,v,ncomp,compkey0,compkey1,compkey,goffset,row;
 19:   PetscInt       numEdges[10],*edgelist[10],asvtx,bsvtx;
 20:   const PetscInt *vtx,*edges;
 21:   PetscBool      sharedv,ghost,distribute=PETSC_TRUE,test=PETSC_FALSE;
 22:   Vec            X;
 23:   Comp0          comp0;
 24:   Comp1          comp1;
 25:   PetscScalar    val;

 27:   PetscInitialize(&argc,&argv,(char*)0,help);
 28:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 29:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 31:   /* Create a network of subnetworks */
 32:   nsubnet = 1;
 33:   if (size == 1) nsubnet = 2;

 35:   /* Create a dmnetwork and register components */
 36:   DMNetworkCreate(PETSC_COMM_WORLD,&dmnetwork);
 37:   DMNetworkRegisterComponent(dmnetwork,"comp0",sizeof(Comp0),&compkey0);
 38:   DMNetworkRegisterComponent(dmnetwork,"comp1",sizeof(Comp1),&compkey1);

 40:   /* Set componnet values - intentionally take rank-dependent value for test */
 41:   comp0.id  = rank;
 42:   comp1.val = 10.0*rank;

 44:   /* Set number of subnetworks, numbers of vertices and edges over each subnetwork */
 45:   DMNetworkSetNumSubNetworks(dmnetwork,nsubnet,PETSC_DECIDE);
 46:   DMNetworkGetNumSubNetworks(dmnetwork,NULL,&Nsubnet);

 48:   /* Input subnetworks; when size>1, process[i] creates subnetwork[i] */
 49:   for (i=0; i<Nsubnet; i++) numEdges[i] = 0;
 50:   for (i=0; i<Nsubnet; i++) {
 51:     if (i == 0 && (size == 1 || (rank == i && size >1))) {
 52:       numEdges[i] = 3;
 53:       PetscMalloc1(2*numEdges[i],&edgelist[i]);
 54:       edgelist[i][0] = 0; edgelist[i][1] = 2;
 55:       edgelist[i][2] = 2; edgelist[i][3] = 1;
 56:       edgelist[i][4] = 1; edgelist[i][5] = 3;

 58:     } else if (i == 1 && (size == 1 || (rank == i && size >1))) {
 59:       numEdges[i] = 3;
 60:       PetscMalloc1(2*numEdges[i],&edgelist[i]);
 61:       edgelist[i][0] = 0; edgelist[i][1] = 3;
 62:       edgelist[i][2] = 3; edgelist[i][3] = 2;
 63:       edgelist[i][4] = 2; edgelist[i][5] = 1;

 65:     } else if (i>1 && (size == 1 || (rank == i && size >1))) {
 66:       numEdges[i] = 3;
 67:       PetscMalloc1(2*numEdges[i],&edgelist[i]);
 68:       for (j=0; j< numEdges[i]; j++) {
 69:         edgelist[i][2*j] = j; edgelist[i][2*j+1] = j+1;
 70:       }
 71:     }
 72:   }

 74:   /* Add subnetworks */
 75:   for (i=0; i<Nsubnet; i++) {
 76:     PetscInt netNum = -1;
 77:     DMNetworkAddSubnetwork(dmnetwork,NULL,numEdges[i],edgelist[i],&netNum);
 78:   }

 80:   /* Add shared vertices -- all processes hold this info at current implementation */
 81:   asvtx = bsvtx = 0;
 82:   for (j=1; j<Nsubnet; j++) {
 83:     /* vertex subnet[0].0 shares with vertex subnet[j].0 */
 84:     DMNetworkAddSharedVertices(dmnetwork,0,j,1,&asvtx,&bsvtx);
 85:   }

 87:   /* Setup the network layout */
 88:   DMNetworkLayoutSetUp(dmnetwork);

 90:   /* Get Subnetwork(); Add nvar=1 to subnet[0] and nvar=2 to other subnets */
 91:   for (net=0; net<Nsubnet; net++) {
 92:     DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
 93:     for (v=0; v<nv; v++) {
 94:       DMNetworkIsSharedVertex(dmnetwork,vtx[v],&sharedv);
 95:       if (sharedv) continue;

 97:       if (!net) {
 98:         DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
 99:       } else {
100:         DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
101:       }
102:     }
103:   }

105:   /* Add components and nvar to shared vertex -- owning and all ghost ranks must call DMNetworkAddComponent() */
106:   DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
107:   for (v=0; v<nv; v++) {
108:     DMNetworkAddComponent(dmnetwork,vtx[v],compkey0,&comp0,1);
109:     DMNetworkAddComponent(dmnetwork,vtx[v],compkey1,&comp1,2);
110:   }

112:   /* Enable runtime option of graph partition type -- must be called before DMSetUp() */
113:   if (size > 1) {
114:     DM               plexdm;
115:     PetscPartitioner part;
116:     DMNetworkGetPlex(dmnetwork,&plexdm);
117:     DMPlexGetPartitioner(plexdm, &part);
118:     PetscPartitionerSetType(part,PETSCPARTITIONERSIMPLE);
119:     PetscOptionsSetValue(NULL,"-dm_plex_csr_alg","mat"); /* for parmetis */
120:   }

122:   /* Setup dmnetwork */
123:   DMSetUp(dmnetwork);

125:   /* Redistribute the network layout; use '-distribute false' to skip */
126:   PetscOptionsGetBool(NULL,NULL,"-distribute",&distribute,NULL);
127:   if (distribute) {
128:     DMNetworkDistribute(&dmnetwork,0);
129:     DMView(dmnetwork,PETSC_VIEWER_STDOUT_WORLD);
130:   }

132:   /* Create a global vector */
133:   DMCreateGlobalVector(dmnetwork,&X);
134:   VecSet(X,0.0);

136:   /* Set X values at shared vertex */
137:   DMNetworkGetSharedVertices(dmnetwork,&nv,&vtx);
138:   for (v=0; v<nv; v++) {
139:     DMNetworkIsGhostVertex(dmnetwork,vtx[v],&ghost);
140:     if (ghost) continue;

142:     /* only one process holds a non-ghost vertex */
143:     DMNetworkGetComponent(dmnetwork,vtx[v],ALL_COMPONENTS,NULL,NULL,&nvar);
144:     DMNetworkGetNumComponents(dmnetwork,vtx[v],&ncomp);
145:     /* PetscPrintf(PETSC_COMM_SELF,"[%d] shared v %D: nvar %D, ncomp %D\n",rank,vtx[v],nvar,ncomp); */
146:     for (j=0; j<ncomp; j++) {
147:       DMNetworkGetComponent(dmnetwork,vtx[v],j,&compkey,NULL,&nvar);
148:       DMNetworkGetGlobalVecOffset(dmnetwork,vtx[v],j,&goffset);
149:       for (i=0; i<nvar; i++) {
150:         row = goffset + i;
151:         val = compkey + 1.0;
152:         VecSetValues(X,1,&row,&val,INSERT_VALUES);
153:       }
154:     }
155:   }
156:   VecAssemblyBegin(X);
157:   VecAssemblyEnd(X);
158:   VecView(X,PETSC_VIEWER_STDOUT_WORLD);

160:   /* Test DMNetworkGetSubnetwork() */
161:   PetscOptionsGetBool(NULL,NULL,"-test_getsubnet",&test,NULL);
162:   if (test) {
163:     net = 0;
164:     PetscOptionsGetInt(NULL,NULL,"-subnet",&net,NULL);
165:     DMNetworkGetSubnetwork(dmnetwork,net,&nv,&ne,&vtx,&edges);
166:     PetscSynchronizedPrintf(PETSC_COMM_WORLD,"[%d] subnet %D: nv %D, ne %D\n",rank,net,nv,ne);
167:     PetscSynchronizedFlush(PETSC_COMM_WORLD,PETSC_STDOUT);
168:     MPI_Barrier(PETSC_COMM_WORLD);

170:     for (i=0; i<nv; i++) {
171:       DMNetworkIsGhostVertex(dmnetwork,vtx[i],&ghost);
172:       DMNetworkIsSharedVertex(dmnetwork,vtx[i],&sharedv);

174:       DMNetworkGetNumComponents(dmnetwork,vtx[i],&ncomp);
175:       if (sharedv || ghost) {
176:         PetscPrintf(PETSC_COMM_SELF,"  [%d] v %D is shared %d, is ghost %d, ncomp %D\n",rank,vtx[i],sharedv,ghost,ncomp);
177:       }

179:       for (j=0; j<ncomp; j++) {
180:         void* component;
181:         DMNetworkGetComponent(dmnetwork,vtx[i],j,&compkey,(void**)&component,NULL);
182:         if (compkey == 0) {
183:           Comp0  *mycomp0;
184:           mycomp0 = (Comp0*)component;
185:           PetscPrintf(PETSC_COMM_SELF,"  [%d] v %D compkey %D, mycomp0->id %D\n",rank,vtx[i],compkey,mycomp0->id);
186:         } else if (compkey == 1) {
187:           Comp1  *mycomp1;
188:           mycomp1 = (Comp1*)component;
189:           PetscPrintf(PETSC_COMM_SELF,"  [%d] v %D compkey %D, mycomp1->val %g\n",rank,vtx[i],compkey,mycomp1->val);
190:         }
191:       }
192:     }
193:   }

195:   /* Free work space */
196:   VecDestroy(&X);
197:   for (i=0; i<Nsubnet; i++) {
198:     if (size == 1 || rank == i) PetscFree(edgelist[i]);
199:   }

201:   DMDestroy(&dmnetwork);
202:   PetscFinalize();
203:   return 0;
204: }

206: /*TEST

208:    build:
209:       requires: !single double defined(PETSC_HAVE_ATTRIBUTEALIGNED)

211:    test:
212:       args:

214:    test:
215:       suffix: 2
216:       nsize: 2
217:       args: -options_left no

219:    test:
220:       suffix: 3
221:       nsize: 4
222:       args: -options_left no

224: TEST*/