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