Actual source code: water.c

  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: #include "water.h"
  7: #include <petscdmnetwork.h>

  9: int main(int argc,char ** argv)
 10: {
 11:   char             waterdata_file[PETSC_MAX_PATH_LEN] = "sample1.inp";
 12:   WATERDATA        *waterdata;
 13:   AppCtx_Water     appctx;
 14: #if defined(PETSC_USE_LOG)
 15:   PetscLogStage    stage1,stage2;
 16: #endif
 17:   PetscMPIInt      crank;
 18:   DM               networkdm;
 19:   PetscInt         *edgelist = NULL;
 20:   PetscInt         nv,ne,i;
 21:   const PetscInt   *vtx,*edges;
 22:   Vec              X,F;
 23:   SNES             snes;
 24:   SNESConvergedReason reason;

 26:   PetscInitialize(&argc,&argv,"wateroptions",help);
 27:   MPI_Comm_rank(PETSC_COMM_WORLD,&crank);

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

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

 36:   PetscLogStageRegister("Read Data",&stage1);
 37:   PetscLogStagePush(stage1);
 38:   PetscNew(&waterdata);

 40:   /* READ THE DATA */
 41:   if (!crank) {
 42:     /* READ DATA. Only rank 0 reads the data */
 43:     PetscOptionsGetString(NULL,NULL,"-waterdata",waterdata_file,sizeof(waterdata_file),NULL);
 44:     WaterReadData(waterdata,waterdata_file);

 46:     PetscCalloc1(2*waterdata->nedge,&edgelist);
 47:     GetListofEdges_Water(waterdata,edgelist);
 48:   }
 49:   PetscLogStagePop();

 51:   PetscLogStageRegister("Create network",&stage2);
 52:   PetscLogStagePush(stage2);

 54:   /* Set numbers of nodes and edges */
 55:   DMNetworkSetNumSubNetworks(networkdm,PETSC_DECIDE,1);
 56:   DMNetworkAddSubnetwork(networkdm,"",waterdata->nedge,edgelist,NULL);
 57:   if (!crank) {
 58:     PetscPrintf(PETSC_COMM_SELF,"water nvertices %D, nedges %D\n",waterdata->nvertex,waterdata->nedge);
 59:   }

 61:   /* Set up the network layout */
 62:   DMNetworkLayoutSetUp(networkdm);

 64:   if (!crank) {
 65:     PetscFree(edgelist);
 66:   }

 68:   /* ADD VARIABLES AND COMPONENTS FOR THE NETWORK */
 69:   DMNetworkGetSubnetwork(networkdm,0,&nv,&ne,&vtx,&edges);

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

 75:   for (i = 0; i < nv; i++) {
 76:     DMNetworkAddComponent(networkdm,vtx[i],appctx.compkey_vtx,&waterdata->vertex[i],1);
 77:   }

 79:   /* Set up DM for use */
 80:   DMSetUp(networkdm);

 82:   if (!crank) {
 83:     PetscFree(waterdata->vertex);
 84:     PetscFree(waterdata->edge);
 85:   }
 86:   PetscFree(waterdata);

 88:   /* Distribute networkdm to multiple processes */
 89:   DMNetworkDistribute(&networkdm,0);

 91:   PetscLogStagePop();

 93:   DMCreateGlobalVector(networkdm,&X);
 94:   VecDuplicate(X,&F);

 96:   /* HOOK UP SOLVER */
 97:   SNESCreate(PETSC_COMM_WORLD,&snes);
 98:   SNESSetDM(snes,networkdm);
 99:   SNESSetOptionsPrefix(snes,"water_");
100:   SNESSetFunction(snes,F,WaterFormFunction,NULL);
101:   SNESSetFromOptions(snes);

103:   WaterSetInitialGuess(networkdm,X);
104:   /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */

106:   SNESSolve(snes,NULL,X);
107:   SNESGetConvergedReason(snes,&reason);

110:   /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */

112:   VecDestroy(&X);
113:   VecDestroy(&F);
114:   SNESDestroy(&snes);
115:   DMDestroy(&networkdm);
116:   PetscFinalize();
117:   return 0;
118: }

120: /*TEST

122:    build:
123:       depends: waterreaddata.c waterfunctions.c
124:       requires: !complex double defined(PETSC_HAVE_ATTRIBUTEALIGNED)

126:    test:
127:       args: -water_snes_converged_reason -options_left no
128:       localrunfiles: wateroptions sample1.inp
129:       output_file: output/water.out
130:       requires: double !complex defined(PETSC_HAVE_ATTRIBUTEALIGNED)

132:    test:
133:       suffix: 2
134:       nsize: 3
135:       args: -water_snes_converged_reason -options_left no
136:       localrunfiles: wateroptions sample1.inp
137:       output_file: output/water.out
138:       requires: double !complex defined(PETSC_HAVE_ATTRIBUTEALIGNED)

140: TEST*/