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:   PetscFunctionBeginUser;
 23:   PetscCall(PetscInitialize(&argc, &argv, NULL, NULL));
 24:   PetscCall(PetscOptionsSetValue(NULL, "-petscpartitioner_use_vertex_weights", "No"));
 25:   comm = PETSC_COMM_WORLD;
 26:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
 27:   PetscCallMPI(MPI_Comm_size(comm, &size));

 29:   PetscCall(DMNetworkCreate(PETSC_COMM_WORLD, &network));

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

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

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

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

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

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

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

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

 77:   PetscCall(PetscFree(edgelist));

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

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

 92: /*TEST

 94:    build:
 95:       requires: !complex double

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

101: TEST*/