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