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: /* 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: char waterdata_file[PETSC_MAX_PATH_LEN] = "sample1.inp";
17: WATERDATA *waterdata;
18: AppCtx_Water appctx;
19: #if defined(PETSC_USE_LOG)
20: PetscLogStage stage1,stage2;
21: #endif
22: PetscMPIInt crank;
23: DM networkdm;
24: PetscInt *edgelist = NULL;
25: PetscInt nv,ne,i;
26: const PetscInt *vtx,*edges;
27: Vec X,F;
28: SNES snes;
29: SNESConvergedReason reason;
31: PetscInitialize(&argc,&argv,"wateroptions",help);
32: MPI_Comm_rank(PETSC_COMM_WORLD,&crank);
34: /* Create an empty network object */
35: DMNetworkCreate(PETSC_COMM_WORLD,&networkdm);
37: /* Register the components in the network */
38: DMNetworkRegisterComponent(networkdm,"edgestruct",sizeof(struct _p_EDGE_Water),&appctx.compkey_edge);
39: DMNetworkRegisterComponent(networkdm,"busstruct",sizeof(struct _p_VERTEX_Water),&appctx.compkey_vtx);
41: PetscLogStageRegister("Read Data",&stage1);
42: PetscLogStagePush(stage1);
43: PetscNew(&waterdata);
45: /* READ THE DATA */
46: if (!crank) {
47: /* READ DATA. Only rank 0 reads the data */
48: PetscOptionsGetString(NULL,NULL,"-waterdata",waterdata_file,sizeof(waterdata_file),NULL);
49: WaterReadData(waterdata,waterdata_file);
51: PetscCalloc1(2*waterdata->nedge,&edgelist);
52: GetListofEdges_Water(waterdata,edgelist);
53: }
54: PetscLogStagePop();
56: PetscLogStageRegister("Create network",&stage2);
57: PetscLogStagePush(stage2);
59: /* Set numbers of nodes and edges */
60: DMNetworkSetNumSubNetworks(networkdm,PETSC_DECIDE,1);
61: DMNetworkAddSubnetwork(networkdm,"",waterdata->nedge,edgelist,NULL);
62: if (!crank) {
63: PetscPrintf(PETSC_COMM_SELF,"water nvertices %D, nedges %D\n",waterdata->nvertex,waterdata->nedge);
64: }
66: /* Set up the network layout */
67: DMNetworkLayoutSetUp(networkdm);
69: if (!crank) {
70: PetscFree(edgelist);
71: }
73: /* ADD VARIABLES AND COMPONENTS FOR THE NETWORK */
74: DMNetworkGetSubnetwork(networkdm,0,&nv,&ne,&vtx,&edges);
76: for (i = 0; i < ne; i++) {
77: DMNetworkAddComponent(networkdm,edges[i],appctx.compkey_edge,&waterdata->edge[i],0);
78: }
80: for (i = 0; i < nv; i++) {
81: DMNetworkAddComponent(networkdm,vtx[i],appctx.compkey_vtx,&waterdata->vertex[i],1);
82: }
84: /* Set up DM for use */
85: DMSetUp(networkdm);
87: if (!crank) {
88: PetscFree(waterdata->vertex);
89: PetscFree(waterdata->edge);
90: }
91: PetscFree(waterdata);
93: /* Distribute networkdm to multiple processes */
94: DMNetworkDistribute(&networkdm,0);
96: PetscLogStagePop();
98: DMCreateGlobalVector(networkdm,&X);
99: VecDuplicate(X,&F);
101: /* HOOK UP SOLVER */
102: SNESCreate(PETSC_COMM_WORLD,&snes);
103: SNESSetDM(snes,networkdm);
104: SNESSetOptionsPrefix(snes,"water_");
105: SNESSetFunction(snes,F,WaterFormFunction,NULL);
106: SNESSetFromOptions(snes);
108: WaterSetInitialGuess(networkdm,X);
109: /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */
111: SNESSolve(snes,NULL,X);
112: SNESGetConvergedReason(snes,&reason);
115: /* VecView(X,PETSC_VIEWER_STDOUT_WORLD); */
117: VecDestroy(&X);
118: VecDestroy(&F);
119: SNESDestroy(&snes);
120: DMDestroy(&networkdm);
121: PetscFinalize();
122: return 0;
123: }
125: /*TEST
127: build:
128: depends: waterreaddata.c waterfunctions.c
129: requires: !complex double defined(PETSC_HAVE_ATTRIBUTEALIGNED)
131: test:
132: args: -water_snes_converged_reason -options_left no
133: localrunfiles: wateroptions sample1.inp
134: output_file: output/water.out
135: requires: double !complex defined(PETSC_HAVE_ATTRIBUTEALIGNED)
137: test:
138: suffix: 2
139: nsize: 3
140: args: -water_snes_converged_reason -options_left no
141: localrunfiles: wateroptions sample1.inp
142: output_file: output/water.out
143: requires: double !complex defined(PETSC_HAVE_ATTRIBUTEALIGNED)
145: TEST*/