Actual source code: ex10.c
1: /*
2: Simple example demonstrating that creating a one sub-network DMNetwork in parallel.
3: */
5: #include <petscdmnetwork.h>
7: int main(int argc,char ** argv)
8: {
9: PetscErrorCode ierr;
10: DM network;
11: PetscMPIInt size,rank;
12: MPI_Comm comm;
13: PetscInt e,ne,nv,v;
14: PetscInt *edgelist = NULL;
15: const PetscInt *nodes,*edges;
16: DM plex;
17: PetscSection section;
18: PetscInt Ne,Ni;
19: PetscInt nodeOffset,k = 2,nedge;
21: PetscInitialize(&argc,&argv,NULL,NULL);if (ierr) return ierr;
22: PetscOptionsSetValue(NULL,"-petscpartitioner_use_vertex_weights","No");
23: comm = PETSC_COMM_WORLD;
24: MPI_Comm_rank(comm,&rank);
25: MPI_Comm_size(comm,&size);
27: DMNetworkCreate(PETSC_COMM_WORLD,&network);
29: Ne = 2;
30: Ni = 1;
31: nodeOffset = (Ne+Ni)*rank; /* The global node index of the first node defined on this process */
33: /* There are three nodes on each rank and two edges. The edges only connect nodes on the given rank */
34: nedge = k * Ni;
36: PetscCalloc1(2*nedge,&edgelist);
37: edgelist[0] = nodeOffset + 0;
38: edgelist[1] = nodeOffset + 2;
39: edgelist[2] = nodeOffset + 1;
40: edgelist[3] = nodeOffset + 2;
42: DMNetworkSetNumSubNetworks(network,PETSC_DECIDE,1);
43: DMNetworkAddSubnetwork(network,"Subnetwork 1",nedge,edgelist,NULL);
44: DMNetworkLayoutSetUp(network);
46: /* Add components and variables for the network */
47: DMNetworkGetSubnetwork(network,0,&nv,&ne,&nodes,&edges);
48: for (e = 0; e < ne; e++) {
49: /* The edges have no degrees of freedom */
50: DMNetworkAddComponent(network,edges[e],-1,NULL,1);
51: }
52: for (v = 0; v < nv; v++) {
53: DMNetworkAddComponent(network,nodes[v],-1,NULL,2);
54: }
56: DMSetUp(network);
57: PetscPrintf(PETSC_COMM_WORLD,"Network after DMSetUp:\n");
58: DMView(network,PETSC_VIEWER_STDOUT_WORLD);
59: DMNetworkGetPlex(network,&plex);
60: /* DMView(plex,PETSC_VIEWER_STDOUT_WORLD); */
61: DMGetLocalSection(plex,§ion);
62: PetscSectionView(section,PETSC_VIEWER_STDOUT_WORLD);
64: PetscFree(edgelist);
66: DMNetworkDistribute(&network,0);
67: PetscPrintf(PETSC_COMM_WORLD,"Network after DMNetworkDistribute:\n");
68: DMView(network,PETSC_VIEWER_STDOUT_WORLD);
69: DMNetworkGetPlex(network,&plex);
70: /* DMView(plex,PETSC_VIEWER_STDOUT_WORLD); */
71: DMGetLocalSection(plex,§ion);
72: PetscSectionView(section,PETSC_VIEWER_STDOUT_WORLD);
74: DMDestroy(&network);
75: PetscFinalize();
76: return ierr;
77: }
79: /*TEST
81: build:
82: requires: !complex double
84: test:
85: nsize: 2
87: TEST*/