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