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: PetscFunctionBeginUser;
27: PetscCall(PetscInitialize(&argc, &argv, "wateroptions", help));
28: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &crank));
30: /* Create an empty network object */
31: PetscCall(DMNetworkCreate(PETSC_COMM_WORLD, &networkdm));
33: /* Register the components in the network */
34: PetscCall(DMNetworkRegisterComponent(networkdm, "edgestruct", sizeof(struct _p_EDGE_Water), &appctx.compkey_edge));
35: PetscCall(DMNetworkRegisterComponent(networkdm, "busstruct", sizeof(struct _p_VERTEX_Water), &appctx.compkey_vtx));
37: PetscCall(PetscLogStageRegister("Read Data", &stage1));
38: PetscCall(PetscLogStagePush(stage1));
39: PetscCall(PetscNew(&waterdata));
41: /* READ THE DATA */
42: if (!crank) {
43: /* READ DATA. Only rank 0 reads the data */
44: PetscCall(PetscOptionsGetString(NULL, NULL, "-waterdata", waterdata_file, sizeof(waterdata_file), NULL));
45: PetscCall(WaterReadData(waterdata, waterdata_file));
47: PetscCall(PetscCalloc1(2 * waterdata->nedge, &edgelist));
48: PetscCall(GetListofEdges_Water(waterdata, edgelist));
49: }
50: PetscCall(PetscLogStagePop());
52: PetscCall(PetscLogStageRegister("Create network", &stage2));
53: PetscCall(PetscLogStagePush(stage2));
55: /* Set numbers of nodes and edges */
56: PetscCall(DMNetworkSetNumSubNetworks(networkdm, PETSC_DECIDE, 1));
57: PetscCall(DMNetworkAddSubnetwork(networkdm, "", waterdata->nedge, edgelist, NULL));
58: if (!crank) PetscCall(PetscPrintf(PETSC_COMM_SELF, "water nvertices %" PetscInt_FMT ", nedges %" PetscInt_FMT "\n", waterdata->nvertex, waterdata->nedge));
60: /* Set up the network layout */
61: PetscCall(DMNetworkLayoutSetUp(networkdm));
63: if (!crank) PetscCall(PetscFree(edgelist));
65: /* ADD VARIABLES AND COMPONENTS FOR THE NETWORK */
66: PetscCall(DMNetworkGetSubnetwork(networkdm, 0, &nv, &ne, &vtx, &edges));
68: for (i = 0; i < ne; i++) PetscCall(DMNetworkAddComponent(networkdm, edges[i], appctx.compkey_edge, &waterdata->edge[i], 0));
70: for (i = 0; i < nv; i++) PetscCall(DMNetworkAddComponent(networkdm, vtx[i], appctx.compkey_vtx, &waterdata->vertex[i], 1));
72: /* Set up DM for use */
73: PetscCall(DMSetUp(networkdm));
75: if (!crank) {
76: PetscCall(PetscFree(waterdata->vertex));
77: PetscCall(PetscFree(waterdata->edge));
78: }
79: PetscCall(PetscFree(waterdata));
81: /* Distribute networkdm to multiple processes */
82: PetscCall(DMNetworkDistribute(&networkdm, 0));
84: PetscCall(PetscLogStagePop());
86: PetscCall(DMCreateGlobalVector(networkdm, &X));
87: PetscCall(VecDuplicate(X, &F));
89: /* HOOK UP SOLVER */
90: PetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
91: PetscCall(SNESSetDM(snes, networkdm));
92: PetscCall(SNESSetOptionsPrefix(snes, "water_"));
93: PetscCall(SNESSetFunction(snes, F, WaterFormFunction, NULL));
94: PetscCall(SNESSetFromOptions(snes));
96: PetscCall(WaterSetInitialGuess(networkdm, X));
97: /* PetscCall(VecView(X,PETSC_VIEWER_STDOUT_WORLD)); */
99: PetscCall(SNESSolve(snes, NULL, X));
100: PetscCall(SNESGetConvergedReason(snes, &reason));
102: PetscCheck(reason >= 0, PETSC_COMM_SELF, PETSC_ERR_CONV_FAILED, "No solution found for the water network");
103: /* PetscCall(VecView(X,PETSC_VIEWER_STDOUT_WORLD)); */
105: PetscCall(VecDestroy(&X));
106: PetscCall(VecDestroy(&F));
107: PetscCall(SNESDestroy(&snes));
108: PetscCall(DMDestroy(&networkdm));
109: PetscCall(PetscFinalize());
110: return 0;
111: }
113: /*TEST
115: build:
116: depends: waterreaddata.c waterfunctions.c
117: requires: !complex double defined(PETSC_HAVE_ATTRIBUTEALIGNED)
119: test:
120: args: -water_snes_converged_reason -options_left no -fp_trap 0
121: localrunfiles: wateroptions sample1.inp
122: output_file: output/water.out
123: requires: double !complex defined(PETSC_HAVE_ATTRIBUTEALIGNED)
125: test:
126: suffix: 2
127: nsize: 3
128: args: -water_snes_converged_reason -options_left no -fp_trap 0
129: localrunfiles: wateroptions sample1.inp
130: output_file: output/water.out
131: requires: double !complex defined(PETSC_HAVE_ATTRIBUTEALIGNED)
133: TEST*/