1: static char help[] = "This example demonstrates the use of DMNetwork interface for solving a steady-state water network model.\n\
  2:                       The water network equations follow those used for the package EPANET\n\
  3:                       The data file format used is from the EPANET package (https://www.epa.gov/water-research/epanet).\n\
  4:                       Run this program: mpiexec -n <n> ./water\n\\n";

  6: /* T
  7:    Concepts: DMNetwork
  8:    Concepts: PETSc SNES solver
  9: */

 11: #include "water.h"
 12:  #include <petscdmnetwork.h>

 14: int main(int argc,char ** argv)
 15: {
 16:   PetscErrorCode   ierr;
 17:   char             waterdata_file[PETSC_MAX_PATH_LEN]="sample1.inp";
 18:   WATERDATA        *waterdata;
 19:   AppCtx_Water     appctx;
 20:   PetscLogStage    stage1,stage2;
 21:   PetscMPIInt      crank;
 22:   DM               networkdm;
 23:   PetscInt         *edgelist = NULL;
 24:   PetscInt         nv,ne,i;
 25:   const PetscInt   *vtx,*edges;
 26:   Vec              X,F;
 27:   SNES             snes;
 28:   SNESConvergedReason reason;

 30:   PetscInitialize(&argc,&argv,"wateroptions",help);if (ierr) return ierr;
 31:   MPI_Comm_rank(PETSC_COMM_WORLD,&crank);

 33:   /* Create an empty network object */
 34:   DMNetworkCreate(PETSC_COMM_WORLD,&networkdm);

 36:   /* Register the components in the network */
 37:   DMNetworkRegisterComponent(networkdm,"edgestruct",sizeof(struct _p_EDGE_Water),&appctx.compkey_edge);
 38:   DMNetworkRegisterComponent(networkdm,"busstruct",sizeof(struct _p_VERTEX_Water),&appctx.compkey_vtx);

 40:   PetscLogStageRegister("Read Data",&stage1);
 41:   PetscLogStagePush(stage1);
 42:   PetscNew(&waterdata);

 44:   /* READ THE DATA */
 45:   if (!crank) {
 46:     /* READ DATA. Only rank 0 reads the data */
 47:     PetscOptionsGetString(NULL,NULL,"-waterdata",waterdata_file,PETSC_MAX_PATH_LEN-1,NULL);
 48:     WaterReadData(waterdata,waterdata_file);

 50:     PetscCalloc1(2*waterdata->nedge,&edgelist);
 51:     GetListofEdges_Water(waterdata,edgelist);
 52:   }
 53:   PetscLogStagePop();

 55:   PetscLogStageRegister("Create network",&stage2);
 56:   PetscLogStagePush(stage2);

 58:   /* Set numbers of nodes and edges */
 59:   DMNetworkSetSizes(networkdm,1,&waterdata->nvertex,&waterdata->nedge,0,NULL);
 60:   if (!crank) {
 61:     PetscPrintf(PETSC_COMM_SELF,"water nvertices %D, nedges %D\n",waterdata->nvertex,waterdata->nedge);
 62:   }

 64:   /* Add edge connectivity */
 65:   DMNetworkSetEdgeList(networkdm,&edgelist,NULL);

 67:   /* Set up the network layout */
 68:   DMNetworkLayoutSetUp(networkdm);

 70:   if (!crank) {
 71:     PetscFree(edgelist);
 72:   }

 74:   /* ADD VARIABLES AND COMPONENTS FOR THE NETWORK */
 75:   DMNetworkGetSubnetworkInfo(networkdm,0,&nv,&ne,&vtx,&edges);

 77:   for (i = 0; i < ne; i++) {
 78:     DMNetworkAddComponent(networkdm,edges[i],appctx.compkey_edge,&waterdata->edge[i]);
 79:   }

 81:   for (i = 0; i < nv; i++) {
 82:     DMNetworkAddComponent(networkdm,vtx[i],appctx.compkey_vtx,&waterdata->vertex[i]);
 83:     /* Add number of variables */
 84:     DMNetworkAddNumVariables(networkdm,vtx[i],1);
 85:   }

 87:   /* Set up DM for use */
 88:   DMSetUp(networkdm);

 90:   if (!crank) {
 91:     PetscFree(waterdata->vertex);
 92:     PetscFree(waterdata->edge);
 93:   }
 94:   PetscFree(waterdata);

 96:   /* Distribute networkdm to multiple processes */
 97:   DMNetworkDistribute(&networkdm,0);

 99:   PetscLogStagePop();

101:   DMCreateGlobalVector(networkdm,&X);
102:   VecDuplicate(X,&F);

104:   /* HOOK UP SOLVER */
105:   SNESCreate(PETSC_COMM_WORLD,&snes);
106:   SNESSetDM(snes,networkdm);
107:   SNESSetOptionsPrefix(snes,"water_");
108:   SNESSetFunction(snes,F,WaterFormFunction,NULL);
109:   SNESSetFromOptions(snes);

111:   WaterSetInitialGuess(networkdm,X);
112:   /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */

114:   SNESSolve(snes,NULL,X);
115:   SNESGetConvergedReason(snes,&reason);
116:   if (reason < 0) {
117:     SETERRQ(PETSC_COMM_SELF,0,"No solution found for the water network");
118:   }
119:   /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */

121:   VecDestroy(&X);
122:   VecDestroy(&F);
123:   SNESDestroy(&snes);
124:   DMDestroy(&networkdm);
125:   PetscFinalize();
126:   return ierr;
127: }

129: /*TEST

131:    build:
132:       depends: waterreaddata.c waterfunctions.c
133:       requires: !complex double define(PETSC_HAVE_ATTRIBUTEALIGNED)

135:    test:
136:       args: -water_snes_converged_reason -options_left no
137:       localrunfiles: wateroptions sample1.inp
138:       output_file: output/water.out
139:       requires: double !complex define(PETSC_HAVE_ATTRIBUTEALIGNED)

141:    test:
142:       suffix: 2
143:       nsize: 3
144:       args: -water_snes_converged_reason -options_left no
145:       localrunfiles: wateroptions sample1.inp
146:       output_file: output/water.out
147:       requires: double !complex define(PETSC_HAVE_ATTRIBUTEALIGNED)

149: TEST*/