Seismic-Reliability-Based Optimal Layout of a Water Distribution Network

We proposed an economic, cost-constrained optimal design of a water distribution system (WDS) that maximizes seismic reliability while satisfying pressure constraints. The model quantifies the seismic reliability of a WDS through a series of procedures: stochastic earthquake generation, seismic intensity attenuation, determination of the pipe failure status (normal, leakage, and breakage), pipe failure modeling in hydraulic simulation, and negative pressure treatment. The network’s seismic reliability is defined as the ratio of the available quantity of water to the required water demand under stochastic earthquakes. The proposed model allows no pipe option in decisions, making it possible to identify seismic-reliability-based optimal layout for a WDS. The model takes into account the physical impact of earthquake events on the WDS, which ultimately affects the network’s boundary conditions (e.g., failure level of pipes). A well-known benchmark network, the Anytown network, is used to demonstrate the proposed model. The network’s optimal topology and pipe layouts are determined from a series of optimizations. The results show that installing large redundant pipes degrades the system’s seismic reliability because the pipes will cause a large rupture opening under failure. Our model is a useful tool to find the optimal pipe layout that maximizes system reliability under earthquakes.


Introduction
An earthquake occurs as a result of a sudden release of energy in the earth's crust.The released massive energy creates seismic waves, which cause deformation of the water distribution system (WDS) of which components are mainly connected beneath the ground.The damage of an earthquake is the multiple failure of system components.For example, for an earthquake occurred in Kobe, Japan in 1995, many pipes were ruptured, causing water to flow out of the system, while some of the pump stations stopped working completely due to power outage.However, the likelihood of concurrence of multiple failures within a system is low under normal conditions.Therefore, a different strategy should be adopted for WDS design of the regions with the risk of earthquake occurrence, such as Japan, Korean peninsula, and West Coast cities of the United States.
During the last two decades, many optimal design approaches have been proposed for WDSs.The early works have mainly taken into account the economic cost as a single objective [1][2][3][4][5].Thereafter, some studies have included the system's performance index within the optimization framework.Lansey et al. [6] developed a chance-constrained least-cost model in which the capacity reliability is defined as the probability that stochastic pressures are equal to or higher than pressure requirement and are constrained at a certain level while minimizing cost.The nodal demands, pressure requirements, and pipe roughness coefficients were assumed to be uncertain while optimizing the pipe size.Several later studies also used chance-constrained optimization models for WDS design [7,8].Kapelan et al. [9] was one of the earliest studies that proposed the multi-objective optimal design of WDS.Since then, most design approaches have adopted two objective functions to minimize economic cost and to maximize system reliability [10,11].The various system reliability indices were suggested to reflect the uncertainties of pipe roughness and system demands when sizing the system.Component failures have been rarely considered in this stream of optimal WDS design.
Su et al. [12] is one of the few studies that considered component failure in WDS design.In their study, reliability was defined as the probability of water demand provision under pipe break (failure) conditions and is used as a constraint for the least-cost design of a WDS.The minimum cut-set (i.e., the most critical set of pipes) was identified by closing the pipe individually for calculating the system reliability.Note that the pipe failure conditions were not considered in hydraulic simulation within the optimization framework, since the computational intensity is overwhelming.
The seismic hazard assessment models have been developed recently.HAZUS [13] was the first model to assess the economic losses of an infrastructure by earthquake.HAZUS is mainly intended to assess the seismic damage but not detailed simulation of the system behavior.Later, the Mid-America Earthquake Center developed a seismic impact assessment model and investigated interdependencies between water and power systems [14].
Early earthquake studies in the water domain focused on investigating individual components' physical behavior under earthquakes rather than quantifying the system-wide performance by modeling the WDS and earthquakes [15][16][17].In a recent study, Fragiadakis et al. [18] proposed a seismic reliability assessment model of a WDS using survival curves of pipes based on general seismic assessment standards and American Lifelines Alliance (ALA) [19] guidelines.However, detailed hydraulic simulations were not conducted in this paper.Later studies began proposing methodologies to evaluate seismic reliability with hydraulic simulations using well-known hydraulic solvers and seismic simulations [20][21][22][23][24][25][26][27][28][29].
The latest and most popular assessment model is the graphical iterative response analysis for flow following earthquakes (GIRAFFE [30]) developed by a research team at Cornell University.GIRAFFE can simulate various pipe leakage and breakage conditions by using EPANET [31].The model's graphical user interface helps to visualize the model results, which is compatible with other geographic information system tools.After the initial development in 2008, the model has been improved and validated through many case studies [27,[32][33][34][35][36].However, GIRAFFE uses a controversial approach to treat negative pressure that removes the nodes' negative pressure and connected pipes.This process is repeated until the negative-pressure nodes are no longer produced.This approach can be time-consuming because the system file must be revised iteratively.
Yoo et al. [37] developed seismic reliability assessment model under stochastic earthquake events.The model quantifies the seismic reliability of a WDS through a series of procedures: stochastic earthquake generation, seismic intensity attenuation, determination of the pipe failure status (normal, leakage, and breakage), pipe failure modeling in hydraulic simulation, and negative pressure treatment.
To the authors' best knowledge, this study is the first attempt to develop an economic, cost-constrained optimal WDS design approach that takes into account seismic reliability based on detailed hydraulic simulations.The proposed model maximizes the system's seismic reliability while satisfying the constraints on economic cost and node pressure requirements.The physical impacts of a seismic wave to WDS components are simulated to determine the failure conditions.The seismic reliability is defined as the ratio of the supplied water to the required demand under stochastic earthquake events.A well-known benchmark network, the Anytown network, is used for the applications of the optimal design and layout maximizing the seismic reliability.

Methodology
In this study, a seismic reliability-based WDS optimization approach is proposed.Figure 1 shows the flowchart of the proposed WDS optimal design model which consists of two sub-models: seismic reliability estimation and optimization model.Seismic reliability estimation model (SREM) first generates stochastic earthquake events by using earthquake generation module.Here, the locations and magnitude of stochastic earthquakes are determined and peak ground acceleration (PGA) reached at each pipe is calculated (described in Section 2.1).Note that no hydraulic calculation is performed from this module.This module also determines the failure mode of pipe where each pipe is classified as normal, leakage, or breakage (Section 2.2).

Methodology
In this study, a seismic reliability-based WDS optimization approach is proposed.Figure 1 shows the flowchart of the proposed WDS optimal design model which consists of two sub-models: seismic reliability estimation and optimization model.Seismic reliability estimation model (SREM) first generates stochastic earthquake events by using earthquake generation module.Here, the locations and magnitude of stochastic earthquakes are determined and peak ground acceleration (PGA) reached at each pipe is calculated (described in Section 2.1).Note that no hydraulic calculation is performed from this module.This module also determines the failure mode of pipe where each pipe is classified as normal, leakage, or breakage (Section 2.2).The failure modes are transferred to hydraulic calculation module where nodal pressures under the failure conditions (i.e., the boundary conditions to the WDS system equation) are calculated (Figure 1 and Section 2.3).Pipe leakage and breakage are modeled by the emitter in EPANET.In order to better simulate the conditions, a semi pressure-driven analysis (semi-PDA) is proposed.After an initial run of EPANET (a hydraulic simulator for DDA), the negative-pressure nodes' demand is set to zero, and the second run is made (Section 2.4).
The resulting hydraulics are provided to seismic reliability calculation module to quantify seismic reliability which is defined as the ratio of available supply to the required demand under stochastic earthquake events.The proposed WDS optimal design model maximizes the system's seismic reliability with constraints on economic cost and pressure requirements (Figure 1 and Section 2.5).Harmony search algorithm (HSA) is used for optimization (Figure 1) [38,39].In the last section of the Methodology, Todini [40]'s resilience indicator is defined.This indicator is not used in the proposed model but used for a postoptimization analysis.
The following subsections describe the details of the methodology, such as earthquake simulation, objective function, optimization approach, and semi-PDA method.Note that the following subsections are in the order that Figure 1 represents.

Earthquake Intensity Attenuation
The proposed model considers the physical characteristics of the WDS's components to determine the failure modes given earthquake intensity.Therefore, supplemental information such as the types of pipes and surrounding soil characteristics should be provided.Network information such as the network layout, nodal demands and coordinates, pipe diameters and lengths, and sizes of the pump and tanks are entered from the EPANET input file.The nodal coordinates are used to calculate the location of a pipe that is used for calculating the distance from the epicenter.
Earthquake intensity is attenuated as it travels from the epicenter.In the model, three attenuation equations were adopted to quantify the damped energy of seismic waves.Kawashima et al. [41] relates the distance from the epicenter and earthquake magnitude to PGA (cm/s 2 ) based on the 90 earthquakes that have occurred in Japan as  The failure modes are transferred to hydraulic calculation module where nodal pressures under the failure conditions (i.e., the boundary conditions to the WDS system equation) are calculated (Figure 1 and Section 2.3).Pipe leakage and breakage are modeled by the emitter in EPANET.In order to better simulate the conditions, a semi pressure-driven analysis (semi-PDA) is proposed.After an initial run of EPANET (a hydraulic simulator for DDA), the negative-pressure nodes' demand is set to zero, and the second run is made (Section 2.4).
The resulting hydraulics are provided to seismic reliability calculation module to quantify seismic reliability which is defined as the ratio of available supply to the required demand under stochastic earthquake events.The proposed WDS optimal design model maximizes the system's seismic reliability with constraints on economic cost and pressure requirements (Figure 1 and Section 2.5).Harmony search algorithm (HSA) is used for optimization (Figure 1) [38,39].In the last section of the Methodology, Todini [40]'s resilience indicator is defined.This indicator is not used in the proposed model but used for a postoptimization analysis.
The following subsections describe the details of the methodology, such as earthquake simulation, objective function, optimization approach, and semi-PDA method.Note that the following subsections are in the order that Figure 1 represents.

Earthquake Intensity Attenuation
The proposed model considers the physical characteristics of the WDS's components to determine the failure modes given earthquake intensity.Therefore, supplemental information such as the types of pipes and surrounding soil characteristics should be provided.Network information such as the network layout, nodal demands and coordinates, pipe diameters and lengths, and sizes of the pump and tanks are entered from the EPANET input file.The nodal coordinates are used to calculate the location of a pipe that is used for calculating the distance from the epicenter.Earthquake intensity is attenuated as it travels from the epicenter.In the model, three attenuation equations were adopted to quantify the damped energy of seismic waves.Kawashima et al. [41] relates the distance from the epicenter and earthquake magnitude to PGA (cm/s 2 ) based on the 90 earthquakes that have occurred in Japan as PGA " 403.8 ˆ10 0.265M ˆpR `30q ´1.218  (1) where M = earthquake magnitude (e.g., M7 stands for earthquake magnitude 7), and R = the distance from the epicenter (km).Lee and Cho [42] presented a similar PGA function that was developed using the historical data of South Korea as log PGA " ´1.83 `0.386M ´logR ´0.0015R Instead of using the horizontal distance from the epicenter, Baag et al. [43] relate the Euclidean distance from focus (∆) and earthquake magnitude to PGA as lnPGA " 0.40 `1.2M ´0.76ln∆ ´0.0094∆ The average of the three calculated PGAs (from Equation (1) to Equation ( 3)) was used in our model.

Determination of Pipe Failure Mode
The fragility of the pipes was determined based on the pipe repair rate (RR, number of pipe repairs/km) and the PGA.In this study, the RR equation suggested by Isoyama et al. [44] and adopted by ALA [19] was used.Isoyama et al. related the RR to the PGA by multiplying some correction factors denoting the pipe types, diameters and soil characteristics (Table 1) as where C1, C2, C3, and C4 represent the correction factors according to the pipe diameter, pipe material, topography, and soil liquefaction, respectively.Because of the lack of information, this study only considered the correction factor of the pipe diameter, C1, in Equation ( 4), while assuming other correction factors are equal to one.As seen in Table 1, for the smaller pipe, C1 increases and results in the higher probability of damage given the same PGA.The post-earthquake pipe status was classified into normal, leakage, or breakage conditions.Pipe leakage is defined as a small rupture on the pipe wall or at the joint through which continuous minor water loss occurs.Pipe breakage is defined as a complete separation of the pipe into two pieces that causes complete loss of transportation ability.ALA [19] suggested that the probability of pipe breakage (P f breakage,i ) is a function of the RR and the length of the pipe (L i , km) expressed as As a rule of thumb, the probability of pipe leakage (P f leakage ,i ) is assumed to be five times higher than the probability of pipe breakage [19].If the cumulative probability of the two conditions is less than 1.0 (e.g., P fbreakage,i = 0.1 and P fleakage,i = 0.5), the complementary probability is assigned for the normal condition (P normal = 0.4).However, if the cumulative probability of the leakage and breakage is greater than 1 (e.g., P fbreakage,i = 0.18 and P fleakage,i = 0.90), the model assumes that the normal condition (without any damage) is not available, and the two failure probabilities are normalized to make the sum of 1.0 (P fbreakage,i = 0.18/(0.18+ 0.90) = 0.17 and P fleakage,i = 0.90/(0.18+ 0.90) = 0.83).

Pipe Failure Modeling
Once the status of each pipe is determined, the information is entered into a network solver, EPANET, for hydraulic simulation (i.e., solving for WDS system equations consisting of conservations of mass and energy).To simulate the pipe leakage, pressure-dependent flows (Q l ) are assigned to the closest node to the damaged pipe using the following equation: where C D = the discharge coefficient (= ˆ2g γ w ˙α A) (the emitter coefficient in EPANET), where g = gravitational acceleration and α = the exponent of the power function; γ w = the specific weight of water; A = the opening area of the damaged pipe; and p = pressure at the closest node.Lambert [45] conducted an experimental study to investigate the accuracy of the power function model in Equation ( 6) and provided a guideline on using exponent α based on the pipe material and level of leakage.It was suggested that the exponent of 0.5 would be used to simulate large detectable leakages in metal pipes; the exponent of 1.0 (linear relationship between discharge and pressure) would be used if no information on the pipe material is available.Puchovsky [46] theoretically derived a discharge coefficient and validated it using sprinkler data.
The shape and area of the opening on ruptured pipe varies depending on pipe material, origin of the pipe damage, and direction of external force.Large opening area is more likely to occur in the failure of large pipes than smaller pipes.To that end, it was assumed that the total opening area is equivalent to 10% of the entire cross-sectional pipe area during leaks.In case of breakage, the entire cross-sectional area was used as the opening area.
If a pipe was tagged as leaking, it was modeled in a hydraulic simulation, as illustrated in Figure 2a, and the discharge coefficient was assigned to the downstream node.The pipe breakage was modeled as shown in Figure 2b.The discharge coefficient was assigned to the upper node in flow direction.Then, the broken pipe was set to "closed" to disconnect the water flow.The demand of the node connected to a broken pipe was modified to consider the degraded delivery capacity due to disconnection.In the model, the nodal demands of both end junctions of the broken pipe were reduced by degrees of node (DoN).Here, DoN is defined as the number of connections/edges that a node has to other nodes in a network.As shown in Figure 2b, the base demand of the upstream node is reduced by 25% (= 1/DoN = 1/4), while the downstream node is reduced by 33.3% (= 1/DoN = 1/3).

Negative Pressure Treatment
Hydraulic analysis approaches for water distribution network can be classified as a DDA and a pressure-driven analysis (PDA).DDA assumes that the demand at individual node is satisfied regardless of the associated pressure, while PDA assumes that nodal discharge is dependent on the pressure head expressed by the head-outflow relationship (HOR).DDA often generates negative pressure when simulating abnormal conditions, such as multiple pipe breaks.PDA, on the other hand, provides more realistic hydraulics avoiding negative pressures.WaterGEMS [47], WDNetXL [48], and WaterNetGen [49] are well-known programs equipped with PDA option based on HOR.However, the system-specific HOR should be provided, which is spatially and temporally inconsistent and also operational dependent.Therefore, for PDA techniques used in pipe network analyses, the analyzer should assume a HOR for the whole network.Given that hydraulic states of pipe networks can vary greatly with changes in HORs, analysis results cannot guarantee high accuracy.
Quasi-PDA methods suppress the occurrence of negative pressure through repetitive DDA analyses.In general, if negative pressure occurs from the first DDA run, better hydraulic calculation results can be derived through quasi-PDA methods by resetting the nodal demands and the components' status.
Ballantyne et al. [21] used the KYPIPE model [50] for hydraulic analyses.KYPIPE, a DDA program similar to EPANET, can also produce negative pressures during simulation of pipeline destruction.Ballantyne et al. calculated system reliability by assuming that water could not be supplied to points with negative pressures, but the hydraulic results around the negative-pressure node still did not have realistic values.The method of Shinozuka et al. [22] and Shi [27] completely removed negative-pressure nodes and the connected pipelines from the original network, but the process required time-consuming repetitive hydraulic analyses and regeneration of input files.
Here, a quasi-PDA approach is proposed for realistic hydraulic simulation under multiple failure conditions to avoid the occurrence of negative pressure.That is, if negative pressure occurred after DDA simulation using EPANET, the updated base demand (after modification

Negative Pressure Treatment
Hydraulic analysis approaches for water distribution network can be classified as a DDA and a pressure-driven analysis (PDA).DDA assumes that the demand at individual node is satisfied regardless of the associated pressure, while PDA assumes that nodal discharge is dependent on the pressure head expressed by the head-outflow relationship (HOR).DDA often generates negative pressure when simulating abnormal conditions, such as multiple pipe breaks.PDA, on the other hand, provides more realistic hydraulics avoiding negative pressures.WaterGEMS [47], WDNetXL [48], and WaterNetGen [49] are well-known programs equipped with PDA option based on HOR.However, the system-specific HOR should be provided, which is spatially and temporally inconsistent and also operational dependent.Therefore, for PDA techniques used in pipe network analyses, the analyzer should assume a HOR for the whole network.Given that hydraulic states of pipe networks can vary greatly with changes in HORs, analysis results cannot guarantee high accuracy.
Quasi-PDA methods suppress the occurrence of negative pressure through repetitive DDA analyses.In general, if negative pressure occurs from the first DDA run, better hydraulic calculation results can be derived through quasi-PDA methods by resetting the nodal demands and the components' status.
Ballantyne et al. [21] used the KYPIPE model [50] for hydraulic analyses.KYPIPE, a DDA program similar to EPANET, can also produce negative pressures during simulation of pipeline destruction.Ballantyne et al. calculated system reliability by assuming that water could not be supplied to points with negative pressures, but the hydraulic results around the negative-pressure node still did not have realistic values.The method of Shinozuka et al. [22] and Shi [27] completely removed negative-pressure nodes and the connected pipelines from the original network, but the process required time-consuming repetitive hydraulic analyses and regeneration of input files.Here, a quasi-PDA approach is proposed for realistic hydraulic simulation under multiple failure conditions to avoid the occurrence of negative pressure.That is, if negative pressure occurred after DDA simulation using EPANET, the updated base demand (after modification described in Chapter 2.3) of the negative-pressure node is set to zero and DDA is repeated.If negative pressure reappears, the pressure of the relevant node is assumed to zero.Compared to Shinozuka et al. [22] and Shi [27] approach, the proposed approach saves overhead processing time by avoiding regeneration of the EPANET input file for the new configuration.In addition to the negative pressure treatment, pressure-dependent water supply is considered in calculating the quantity of available water at each node.The available nodal demand at node j (Q avl,j ) is estimated as Q new,j ˆd P j P min if 0 ă P j ă P min Q new,j if P j ě P min (7) where Q new,j = the updated base demand with pipe breakage consideration and negative pressure treatment at node j; P j = the pressure head at node j; and P min " the minimum pressure requirement.

Seismic Reliability Indicator
Various surrogate measures of WDS reliability have been proposed including capacity reliability [6,7], resilience [40], robustness [10,11], availability [51].Bao and Mays [52] proposed three formulations of WDS system reliability which can be calculated from nodal reliability values: minimum nodal reliability, arithmetic mean reliability, and flow-weighted mean reliability.The first one concerns on the worst nodal reliability while the latter two values indicate system-wide reliability.In order to represent post-earthquake system performance in the water supply, a seismic reliability measure should be able to reflect system-wide water availability rather than local level of system performance.
To quantify system-wide seismic reliability, a new reliability indicator is proposed herein.The system seismic reliability (S S ) is defined as the ratio of the total available system demand to the total required system demand: where m = total number of demand nodes; and Q req,j = the required demand at node j.
A water system planner would intend to design the most reliable network for earthquakes under a given budget condition.To that end, the proposed model here is a single-objective optimal design model that maximizes the system's seismic reliability with a constraint on economic cost.The optimization model is formulated as follows: Maximize F " S s s.t.CC ď CC given (9) where CC = the pipe construction cost; CC given = the available budget for pipe system.Equation ( 9) is valid for each set of available commercial pipe diameters and in conditions where the minimum pressure head is guaranteed.The pipe construction cost accounting for the pipe material and installation is expressed as: where uc pD i q = the unit cost of the pipe with diameter D i determined for the ith pipe (USD/m); and L i = the length of the ith pipe.

Resilience Indicator
In this study, Todini [40]'s resilience is quantified from a range of designs of the Anytown network and compared with seismic reliability.Todini [40] introduced a resilience index to quantify the resilience of looped network.Resilience is defined as surplus power available within the network as a percentage of net input power as: Power i γ ´řm j Q req,j H req (11) where H j = total head of the jth node; Q k = flow provided by reservoir k (m 3 /s); H k = head at reservoir k; Power i = power of the ith pump (Nm/s); γ = specific weight of water (N/m 3 ); nr and np = number of reservoir and pumps, respectively.Note that this indicator was used only for a postoptimization analysis.Todini's resilience is one of the most popular and widely used surrogate measures of WDS reliability.Farmani et al.
[53] investigated the trade-off between economic cost and the resilience for a rehabilitation problem of the Anytown network.Prasad and Park [54] proposed a multiobjective optimization approach to minimize cost and maximize modified version of Todini's resilience indicator.Recently, Gheisi and Naser [55] compared entropy-based reliability, Todini's resilience, and three modified versions of Todini's resilience in the twenty-two potential pipe layouts of a hypothetical network.

Study Network
The proposed optimization approach is applied for optimal design of a well-known benchmark WDS, the Anytown network.The network, which was firstly published by Walski et al. [56], was modified by Jung et al. [11] for pipe only and pipe/pump designs that minimize total cost and maximize system robustness.The original benchmark network was modified by removing the two tanks and the connected riser pipes to be solely supplied by a single reservoir with a fixed source head.The fixed source head is elevated from 3 m (10 ft) to 73.2 m (240 ft).A peaking factor of 1.8 is applied to given average-based demand (2005 average daily use) to create the daily peak condition.
Other modifications of the study network were made for application purposes.First, this study suggests new pipe layouts and sizes assuming there are no existing pipes in the system.Pump and tank design are not considered.Second, in addition to 10 commercial pipe sizes (152, 203, 254, 305, 356, 406, 457, 508, 610, and 762 mm), zero diameter can be selected suggesting no pipe installation for the potential path.The unit costs of the commercial pipes are adopted from Walski et al. [56].
To generate random epicenters, a 2 ˆ2 grid was created and laid on the study network as seen in Figure 3.The rectangular boundary of the grid is defined by the four end nodes: north, south, east, and west.Total 900 earthquakes are generated and consistent number of earthquakes are assigned at each corner of the grid (marked as "x" in Figure 3).The earthquake intensity reached at the pipe is a function of the earthquake magnitude and the Euclidean distance from the epicenter.Note the earthquake magnitude of M4 was generated from each corner of the grid.The minimum pressure requirement of 28.12 m (40 psi) should be satisfied under based demand condition and is also considered in Equation (7).and west.Total 900 earthquakes are generated and consistent number of earthquakes are assigned at each corner of the grid (marked as "x" in Figure 3).The earthquake intensity reached at the pipe is a function of the earthquake magnitude and the Euclidean distance from the epicenter.Note the earthquake magnitude of M4 was generated from each corner of the grid.The minimum pressure requirement of 28.12 m (40 psi) should be satisfied under based demand condition and is also considered in Equation (7).HSA was used to find an optimal solution for the pipe sizing problem.HSA was inspired by musical performance process and widely used for WDS optimizations.While the applications encompass from pipe network design to pump scheduling, HSA was proven to be generally outperform other algorithms such as genetic algorithm [57][58][59].
Several assumptions and simplifications were made in this study: (1) the seismic damages of the pump, tank, reservoir, and valve are not considered; (2) the earthquake's focal depth (tectonics) is assumed to be 10 km; (3) the coordinate of the center of the pipe is used as a reference point for calculating the distance from the epicenter; (4) Pipes are cast iron pipes; and (5) the Anytown network is laid on the alluvial plain with no liquefaction.Based on the assumption (4) and ( 5), all correction factors except C1 in the Equation (4) (C2-4) are equal to one.

Application Results
First, we applied two different commercial pipe sets for the pipe sizing of the study network.This analysis investigates the WDS's layout changes with respect to the seismic reliability increase.The same total cost constraint was applied for the two designs and the resultant pipe layouts were compared.Then, seven designs with a fixed layout and different redundancy levels were compared for the systems' seismic reliability.Finally, discussions and suggestions on improving WDS seismic reliability were provided at the end of this section.

Different Available Pipe Sizing Options
The first application is intended to investigate the impact of seismic reliability on the network's layout.Two sets of commercial pipe sizes are assumed to be available for two design cases.In Case 1, all commercial pipes (152, 203, 254, 305, 356, 406, 457, 508, 610, and 762 mm) except the zero size option are available.In Case 2, zero diameter is available in addition to the commercial sizes in Case 1. Considering the same cost constraint (CC given = 18 million (M) USD) in Equation ( 9) provided a platform for consistent comparison of the resulting designs.
Figure 4 shows the optimal seismic reliability values of the two case designs.The corresponding optimal pipe layouts are presented in Figure 5. Contrary to expectations, the seismic reliability decreased with the availability of more pipe sizes.In Case 1, at least a 152 mm pipe should be installed because no pipe option is not available.By comparing Figures 5a,b and 6, we can observe that 152 mm pipes were installed in Case 1 for the link at which no pipe was constructed in Case 2. The 152 mm pipe, which is the most vulnerable pipe to earthquakes, almost always causes failure.Being able to have no pipe instead of a 152 mm pipe increased the system's seismic reliability by 0.2 (a 20% increase in the amount of available water during an earthquake).The network layout becomes smaller from Case 2 to Case 1 (Figure 5a,b), while the overall pipe diameters also increase (Figure 6).
optimal pipe layouts are presented in Figure 5. Contrary to expectations, the seismic reliability decreased with the availability of more pipe sizes.In Case 1, at least a 152 mm pipe should be installed because no pipe option is not available.By comparing Figure 5a,b and Figure 6, we can observe that 152 mm pipes were installed in Case 1 for the link at which no pipe was constructed in Case 2. The 152 mm pipe, which is the most vulnerable pipe to earthquakes, almost always causes failure.Being able to have no pipe instead of a 152 mm pipe increased the system's seismic reliability by 0.2 (a 20% increase in the amount of available water during an earthquake).The network layout becomes smaller from Case 2 to Case 1 (Figure 5a,b), while the overall pipe diameters also increase (Figure 6).Note that this is very different from what we observed in a traditional capacity reliability-based design.In the context of the capacity reliability field, it was believed that having additional paths Note that this is very different from what we observed in a traditional capacity reliability-based design.In the context of the capacity reliability field, it was believed that having additional paths and more loops would result in the increase of system reliability and redundancy [11,60].However, the result of this study indicates that a different strategy should be available under the conditions where WDS component failures are affected by strong external forces (i.e., earthquakes) and the components' physical characteristics (pipe's probability of failure as a function of RR). and more loops would result in the increase of system reliability and redundancy [11,60].However, the result of this study indicates that a different strategy should be available under the conditions where WDS component failures are affected by strong external forces (i.e., earthquakes) and the components' physical characteristics (pipe's probability of failure as a function of RR).This plateau in seismic reliability was also observed in the optimal pipe designs of the Anytown network by applying different total cost constraints and using the available pipe sizes in Case 1. Figure 7 shows the Pareto optimal solutions' total cost and seismic reliability.The marginal cost and more loops would result in the increase of system reliability and redundancy [11,60].However, the result of this study indicates that a different strategy should be available under the conditions where WDS component failures are affected by strong external forces (i.e., earthquakes) and the components' physical characteristics (pipe's probability of failure as a function of RR).This plateau in seismic reliability was also observed in the optimal pipe designs of the Anytown network by applying different total cost constraints and using the available pipe sizes in Case 1. Figure 7 shows the Pareto optimal solutions' total cost and seismic reliability.The marginal cost becomes infinite for the designs whose cost is greater than 16.5 M USD.A reliability increase can no longer be achieved once a sufficient investment is made.Although the designs greater than 16.5 M This plateau in seismic reliability was also observed in the optimal pipe designs of the Anytown network by applying different total cost constraints and using the available pipe sizes in Case 1. Figure 7 shows the Pareto optimal solutions' total cost and seismic reliability.The marginal cost becomes infinite for the designs whose cost is greater than 16.5 M USD.A reliability increase can no longer be achieved once a sufficient investment is made.Although the designs greater than 16.5 M USD have more large pipes compared to the solutions less than 16.5 M USD (the number of pipes equal to or larger than 508 mm is between 7 and 10, while the designs less than 16.5 M USD have three to five pipes), no benefit of having larger pipes was obtained with respect to seismic reliability.For effective and economical improvement of WDS seismic reliability, the threshold investment for a network should first be identified.

Constant Layout with a Single Pipe Sizing Option
The impacts of having large pipes are also investigated through the seismic reliability evaluation of seven uniform designs.The Design 1 has 305 mm for all pipes in the study network.Design 2, 3, 4, 5, 6, and 7 have 356, 406, 457, 508, 610, and 762 mm, respectively, for all pipes.The seismic reliabilities of the seven designs are shown in Figure 8.For comparison, Todini's resilience is also calculated from the seven designs and plotted in Figure 8.While there is a large increase in seismic reliability from the 457 mm design (Design 4) to the 508 mm design (Design 5), reliability decreases from the 508 mm design to the 762 mm design (Design 7).This explains why we observed a plateau in seismic reliability in Figure 7. On the other hand, resilience (a traditional reliability measure) consistently increases with increasing pipe sizes.The marginal cost of improving resilience increases substantially for a resilience value of 0.3-0.8 and stabilized for a value higher than 0.8.As the pipe sizes decrease, the correction factor for RR (Equation ( 4)) increases, finally resulting in high pipe breakage and leakage probability.However, although the pipes failed, the resulting impact is not significant to the system's seismic reliability because the calculated discharge

Constant Layout with a Single Pipe Sizing Option
The impacts of having large pipes are also investigated through the seismic reliability evaluation of seven uniform designs.The Design 1 has 305 mm for all pipes in the study network.Design 2, 3, 4, 5, 6, and 7 have 356, 406, 457, 508, 610, and 762 mm, respectively, for all pipes.The seismic reliabilities of the seven designs are shown in Figure 8.For comparison, Todini's resilience is also calculated from the seven designs and plotted in Figure 8.While there is a large increase in seismic reliability from the 457 mm design (Design 4) to the 508 mm design (Design 5), reliability decreases from the 508 mm design to the 762 mm design (Design 7).This explains why we observed a plateau in seismic reliability in Figure 7. On the other hand, resilience (a traditional reliability measure) consistently increases with increasing pipe sizes.The marginal cost of improving resilience increases substantially for a resilience value of 0.3-0.8 and stabilized for a value higher than 0.8.

Constant Layout with a Single Pipe Sizing Option
The impacts of having large pipes are also investigated through the seismic reliability evaluation of seven uniform designs.The Design 1 has 305 mm for all pipes in the study network.Design 2, 3, 4, 5, 6, and 7 have 356, 406, 457, 508, 610, and 762 mm, respectively, for all pipes.The seismic reliabilities of the seven designs are shown in Figure 8.For comparison, Todini's resilience is also calculated from the seven designs and plotted in Figure 8.While there is a large increase in seismic reliability from the 457 mm design (Design 4) to the 508 mm design (Design 5), reliability decreases from the 508 mm design to the 762 mm design (Design 7).This explains why we observed a plateau in seismic reliability in Figure 7. On the other hand, resilience (a traditional reliability measure) consistently increases with increasing pipe sizes.The marginal cost of improving resilience increases substantially for a resilience value of 0.3-0.8 and stabilized for a value higher than 0.8.As the pipe sizes decrease, the correction factor for RR (Equation ( 4)) increases, finally resulting in high pipe breakage and leakage probability.However, although the pipes failed, the resulting impact is not significant to the system's seismic reliability because the calculated discharge coefficient, which is a function of the pipe's cross-sectional area (Table 2), is not large.However, as seen in Table 2, large pipes such as 508, 610, and 762 mm have a smaller failure probability As the pipe sizes decrease, the correction factor for RR (Equation ( 4)) increases, finally resulting in high pipe breakage and leakage probability.However, although the pipes failed, the resulting impact is not significant to the system's seismic reliability because the calculated discharge coefficient, which is a function of the pipe's cross-sectional area (Table 2), is not large.However, as seen in Table 2, large pipes such as 508, 610, and 762 mm have a smaller failure probability compared to small pipes, but the resulting discharge coefficient is more than 10 times as big as that of small pipes.The failure effects are more significant for system seismic reliability compared to small pipes.

Random Designs
Finally, random designs that satisfy pressure requirements are generated from the proposed model to confirm the aforementioned conclusion.Figure 9 shows the profiles of seismic reliability and resilience of many random designs.We can clearly see that there is an inflection point around the total cost of 23 M USD and seismic reliability around 0.4, after which the overall seismic reliability decreases (Figure 9a).Because the total cost is a direct function of the pipe diameter, this plot indicates that installing large pipes does not always guarantee an increase in seismic reliability.The printed solutions are all suboptimal solutions that are dominated by the Pareto solutions found in Figure 7. On the other hand, installing large pipes resulted in an upward trend in resilience (Figure 8b).

Random Designs
Finally, random designs that satisfy pressure requirements are generated from the proposed model to confirm the aforementioned conclusion.Figure 9 shows the profiles of seismic reliability and resilience of many random designs.We can clearly see that there is an inflection point around the total cost of 23 M USD and seismic reliability around 0.4, after which the overall seismic reliability decreases (Figure 9a).Because the total cost is a direct function of the pipe diameter, this plot indicates that installing large pipes does not always guarantee an increase in seismic reliability.The printed solutions are all suboptimal solutions that are dominated by the Pareto solutions found in Figure 7. On the other hand, installing large pipes resulted in an upward trend in resilience (Figure 8b).Under normal failure condition, having more additional paths and installing large pipes throughout the system are beneficial with respect to system reliability (i.e., ability to supply required quantity of water) and redundancy (i.e., level of pressure redundancy).However, we observed this does not apply under earthquake failures.Earthquake deforms WDS components and their original function.Once failed under earthquakes, large pipes, which can deliver large volume of water under normal condition, help accelerate water loss out of the system.This study was the first attempt to take into account such irregular hydraulic behavior in the WDS modeling and design and highlight the need to find the optimal pipe layout considering the system's performances under the two different states.Under normal failure condition, having more additional paths and installing large pipes throughout the system are beneficial with respect to system reliability (i.e., ability to supply required quantity of water) and redundancy (i.e., level of pressure redundancy).However, we observed this does not apply under earthquake failures.Earthquake deforms WDS components and their original function.Once failed under earthquakes, large pipes, which can deliver large volume of water under normal condition, help accelerate water loss out of the system.This study was the first attempt to take into account such irregular hydraulic behavior in the WDS modeling and design and highlight the need to find the optimal pipe layout considering the system's performances under the two different states.

Summary and Conclusions
Earthquakes can cause many simultaneous failures of components throughout WDSs, which result in very different and severe conditions compared to normal failure conditions that water communities usually deal with (e.g., single pipe failure).In this study, an economic, cost-constrained optimal design approach of a WDS was proposed to maximize seismic reliability.The network's seismic reliability is defined as the ratio of the available supply to the required water demand during stochastic earthquakes.The physical relationship between the earthquake intensity and the WDS components' vulnerability characteristics was defined and utilized for quantifying seismic reliability.To overcome the limitation of DDA in modeling earthquake failures, negative pressures were assumed to have realistic hydraulic calculations.Then, we investigated the seismic reliability improvement with respect to the economic investment in the system and the associated topological and pipe diameter changes.A traditional benchmark network, the Anytown network, was used to demonstrate the approach.Random earthquakes around the study network were generated for reliability quantification.
The results were quite different from what we generally observed from traditional reliability-based design (e.g., resilience-based design).First, allowing redundant small pipes in the system did not help improve seismic reliability.Those pipes cause frequent failures because of their low durability to earthquake forces.Second, having too-large pipes also degrades the system reliability during earthquakes.Compared to small pipes that are 152 mm to 305 mm, a big pipe (762 mm) has a lower failure probability; however, the failure's influence is very significant to the system once it fails.The large pipe's cross-sectional area will help release more water out of the system.Therefore, the first step to efficiently improve the WDS seismic reliability should be to identify the most appropriate pipe sizes for a system.For example, from the Anytown network, having a uniform 508 mm for all pipes provides the highest system reliability of 0.38 among seven uniform designs.
This study has several limitations that future research must address.First, this model is sensitive to correction factors for the pipe and the assumption of the opening area of the pipe under failure.Intensive sensitivity analysis should be conducted with experimental studies to simulate the most realistic failure behavior of WDS components.Second, this study only considers the pipe failures while pump, valve, tank, and reservoir failure can occur during an earthquake.Therefore, various failure types could be included in future studies.The proposed semi-PDA approach and EPANet can be replaced with a PDA-based network solver in order to simulate more realistic behavior of WDS under earthquake events.Third, while this study focuses on the system's reliability right after an earthquake occurs, post-earthquake recovery strategies should also be investigated to enhance the overall WDS's reliability/resilience.This could be found mostly in the context of operation and management.Fourth, the proposed semi-PDA approach and EPANet can be replaced with a PDA-based network solver in order to simulate more realistic behavior of WDS under earthquake events.Finally, interdependencies among multiple lifeline infrastructures (e.g., the water, power, transportation, and communication systems) can help improve and recover an individual infrastructure's reliability during an earthquake.Therefore, more efforts should be made to identify these interdependencies.

Figure 1 .
Figure 1.Flowchart of the proposed water distribution system (WDS) optimal design model.

Figure 1 .
Figure 1.Flowchart of the proposed water distribution system (WDS) optimal design model.

Figure 2 .
Figure 2. A schematic diagram to describe pipe damage modeling: (a) leakage condition and (b) breakage condition.

Figure 2 .
Figure 2. A schematic diagram to describe pipe damage modeling: (a) leakage condition and (b) breakage condition.

Figure 4 .
Figure 4. Maximum seismic reliability values of two case designs with different pipe sizing options.

Figure 4 .
Figure 4. Maximum seismic reliability values of two case designs with different pipe sizing options.

Figure 5 .
Figure 5. Pipe layout comparison for the solutions obtained from Cases 1 and 2 (Figure 4); pipe diameters are in mm; the thicker and darker pipe is larger.(a) Case 1; (b) Case 2.

Figure 6 .
Figure 6.Histogram of pipes by the pipe diameter.

Figure 5 .
Figure 5. Pipe layout comparison for the solutions obtained from Cases 1 and 2 (Figure 4); pipe diameters are in mm; the thicker and darker pipe is larger.(a) Case 1; (b) Case 2.

Figure 5 .
Figure 5. Pipe layout comparison for the solutions obtained from Cases 1 and 2 (Figure 4); pipe diameters are in mm; the thicker and darker pipe is larger.(a) Case 1; (b) Case 2.

Figure 6 .
Figure 6.Histogram of pipes by the pipe diameter.

Figure 6 .
Figure 6.Histogram of pipes by the pipe diameter.

Figure 8 .
Figure 8. Seismic reliability of the seven uniform designs; all pipes are 508 mm in Case 5.

Figure 8 .
Figure 8. Seismic reliability of the seven uniform designs; all pipes are 508 mm in Case 5.

Figure 8 .
Figure 8. Seismic reliability of the seven uniform designs; all pipes are 508 mm in Case 5.

Figure 9 .
Figure 9. (a) Seismic reliability of randomly generated solutions (black dot) and Pareto solutions shown in Figure 7 (blue diamond) and (b) resilience of the same random designs.

Figure 9 .
Figure 9. (a) Seismic reliability of randomly generated solutions (black dot) and Pareto solutions shown in Figure 7 (blue diamond) and (b) resilience of the same random designs.

Table 2 .
Correction factors and discharge coefficients (Equation (6)) for the pipe sizes considered.

Table 2 .
Correction factors and discharge coefficients (Equation (6)) for the pipe sizes considered.