Actual source code: ex10.c

  1: /*
  2:     Simple example demonstrating creating a one sub-network DMNetwork in parallel.

  4:     In this example vertices 0 and 1 are not connected to any edges.
  5: */

  7: #include <petscdmnetwork.h>

  9: int main(int argc,char ** argv)
 10: {
 11:   DM                network;
 12:   PetscMPIInt       size,rank;
 13:   MPI_Comm          comm;
 14:   PetscInt          e,ne,nv,v,ecompkey,vcompkey;
 15:   PetscInt          *edgelist = NULL;
 16:   const PetscInt    *nodes,*edges;
 17:   DM                plex;
 18:   PetscSection      section;
 19:   PetscInt          Ne,Ni;
 20:   PetscInt          nodeOffset,k = 2,nedge;

 22:   PetscInitialize(&argc,&argv,NULL,NULL);
 23:   PetscOptionsSetValue(NULL,"-petscpartitioner_use_vertex_weights","No");
 24:   comm = PETSC_COMM_WORLD;
 25:   MPI_Comm_rank(comm,&rank);
 26:   MPI_Comm_size(comm,&size);

 28:   DMNetworkCreate(PETSC_COMM_WORLD,&network);

 30:   /* Register zero size componets to get compkeys to be used by DMNetworkAddComponent() */
 31:   DMNetworkRegisterComponent(network,"ecomp",0,&ecompkey);
 32:   DMNetworkRegisterComponent(network,"vcomp",0,&vcompkey);

 34:   Ne = 2;
 35:   Ni = 1;
 36:   nodeOffset = (Ne+Ni)*rank;   /* The global node index of the first node defined on this process */

 38:   /* There are three nodes on each rank and two edges. The edges only connect nodes on the given rank */
 39:   nedge = k * Ni;

 41:   if (rank == 0) {
 42:     nedge = 1;
 43:     PetscCalloc1(2*nedge,&edgelist);
 44:     edgelist[0] = nodeOffset + 2;
 45:     edgelist[1] = nodeOffset + 3;
 46:   } else {
 47:     nedge = 2;
 48:     PetscCalloc1(2*nedge,&edgelist);
 49:     edgelist[0] = nodeOffset + 0;
 50:     edgelist[1] = nodeOffset + 2;
 51:     edgelist[2] = nodeOffset + 1;
 52:     edgelist[3] = nodeOffset + 2;
 53:   }

 55:   DMNetworkSetNumSubNetworks(network,PETSC_DECIDE,1);
 56:   DMNetworkAddSubnetwork(network,"Subnetwork 1",nedge,edgelist,NULL);
 57:   DMNetworkLayoutSetUp(network);

 59:   PetscPrintf(PETSC_COMM_WORLD,"Network after DMNetworkLayoutSetUp:\n");
 60:   DMView(network,PETSC_VIEWER_STDOUT_WORLD);

 62:   /* Add components and variables for the network */
 63:   DMNetworkGetSubnetwork(network,0,&nv,&ne,&nodes,&edges);
 64:   for (e = 0; e < ne; e++) {
 65:     /* The edges have no degrees of freedom */
 66:     DMNetworkAddComponent(network,edges[e],ecompkey,NULL,1);
 67:   }
 68:   for (v = 0; v < nv; v++) {
 69:     DMNetworkAddComponent(network,nodes[v],vcompkey,NULL,2);
 70:   }

 72:   DMSetUp(network);
 73:   DMNetworkGetPlex(network,&plex);
 74:   /* DMView(plex,PETSC_VIEWER_STDOUT_WORLD); */
 75:   DMGetLocalSection(plex,&section);
 76:   PetscSectionView(section,PETSC_VIEWER_STDOUT_WORLD);

 78:   PetscFree(edgelist);

 80:   DMNetworkDistribute(&network,0);
 81:   PetscPrintf(PETSC_COMM_WORLD,"\nNetwork after DMNetworkDistribute:\n");
 82:   DMView(network,PETSC_VIEWER_STDOUT_WORLD);
 83:   DMNetworkGetPlex(network,&plex);
 84:   /* DMView(plex,PETSC_VIEWER_STDOUT_WORLD); */
 85:   DMGetLocalSection(plex,&section);
 86:   PetscSectionView(section,PETSC_VIEWER_STDOUT_WORLD);

 88:   DMDestroy(&network);
 89:   PetscFinalize();
 90:   return 0;
 91: }

 93: /*TEST

 95:    build:
 96:       requires: !complex double

 98:    test:
 99:       nsize: 2
100:       args: -petscpartitioner_type simple

102: TEST*/