Probabilistic modeling of sustainable urban drainage systems

Sustainable Urban Drainage Systems (SUDS) include strategies and solutions for distributed stormwater management and control. They are strongly encouraged, especially in highly urbanized areas that suffer the combined effect of impervious surfaces and the increase in extreme rainfall events due to urbanization growth and climate change. Their integration into traditional urban drainage systems can mitigate flood risk and pollution of receiving water bodies. The main goal of SUDS is to restore the natural water balance by increasing infiltration and evapotranspiration processes and promoting rainwater harvesting and reuse. This paper proposes an analytical-probabilistic approach for SUDS modeling applicable to different systems. Developed equations allow estimating the runoff and residual storage probability for evaluating the efficiency of the storage volume both in terms of flood control and, depending on SUDS type, in terms of emptying time or water needs supply. The modeling considers the possibility of consecutive chained rainfalls; this feature is relevant for SUDS, often characterized by low outflow rates. Relating characteristic parameters to a probabilistic level (the Average Return Interval, ARI) makes the formulas interesting to be used in the design practice. An application to two case studies confirmed the goodness of the proposed method.


Introduction
The combination of more frequent extreme events due to climate change, and the increase of impervious surfaces due to urbanization, causes the alteration of the natural water balance. The alternation of high-intensity rainfalls to water scarcity periods involves heavy damage and issues, especially in the most fragile areas. In the natural water balance, precipitation is largely evapotranspirated from vegetation and soil and infiltrated into underlying layers, and only a limited percentage produces runoff. In highly urbanized contexts, the natural water balance is upset since only a small amount of precipitation can be infiltrated and evapotranspirated; moreover, wastewater discharges and poor-quality runoffs increase even for potable waters imported to satisfy urban water needs. SUDS are best practices that aim to reduce the flow of stormwater reaching combined and separate sewers, restoring the natural water balance. They are based on Nature-Based Solutions (NBS), promoting infiltration and evapotranspiration processes, and on stormwater harvesting and reuse, reducing flood risk and potable water consumption. SUDS allow meeting the hydrologic invariance, for which both discharged peak flow rates and runoff volumes can be limited to the pre-urbanization condition. If considering only traditional urban drainage systems (with detention tanks and over-size pipes for runoff control), only the peak flow rates get the pre-urbanization conditions (hydraulic invariance target). SUDS belong to three main categories: infiltration systems, which include permeable pavements, swales, rain gardens, soakaways, infiltration trenches, and infiltration basins; evapotranspiration systems, which include green roofs, retention basins, wetlands, and rainwater harvesting systems (Fig. 1).
SUDS are not only best practices for flood risk mitigation; they involve several benefits to the surrounding environment, such as runoff quality enhancement, health, wellbeing, amenity, recreation, thermal comfort, education, and biodiversity.
Depending on the kind of SUDS, they achieve different targets. Infiltration systems must limit runoff and be emptied in about 24-48 h to restore the storage capacity. On the contrary, rainwater harvesting systems and evapotranspiration systems such as green roofs benefit from residual storage during dry periods: the first to satisfy water supply requirements and the second to guarantee a residual humidity for vegetation survival without irrigation; this feature must coexist with the need to limit runoffs into the drainage system.
In literature, there are several studies on SUDS that evaluate their performance and analyze their characteristic parameters (Berndtsson 2010;Carter and Rasmussen 2007;Hakimdavar et al. 2014;Herrera et al. 2018;Lee 2019;Li et al. 2017;Marchioni et al. 2022;Newman et al. 2013;Palermo et al. 2019;Palla et al. 2012). Traditional methods are based on continuous simulation, design approaches, or experimental formulas. In the last decades, analyticalprobabilistic approaches have often been proposed as a trade-off to join the simplicity of design storm methods and the accuracy of the continuous simulation. When compared with design storm methods, analytical-probabilistic approaches relate design variables to a probabilistic level (linked to an Average Return Interval, ARI) and consider the possibility of residual storage from previous rainfall events, i.e., a partial pre-filling of the storage volume. When compared with continuous simulation, analytical-probabilistic approaches only need the average values of rainfall depth, rainfall duration, and interevent time as input variables, with no need for complete series of recorded rainfalls, which are often not available in the medium-long term.
This paper proposes an analytical-probabilistic approach suitable for the SUDS modeling with different configurations and typologies. Developed equations allow estimating both the runoff probability and the probability of residual storage from previous events. The combination of these two features allows a reliable approach for SUDS design, relating the need for runoff control with that of restoring the whole storage volume or, on the contrary, to have sufficient storage to satisfy water demand.
Literature about analytical-probabilistic approaches generally considered only a couple of chained consecutive rainfall events. Raimondi and Becciu (2014) concluded that this assumption is usually adequate only for small storage volumes characterized by great outflow rates. For SUDS, often characterized by low release flow rates, or when strict limitations on discharges into the downstream drainage system occur, considering more than two consecutive chained rainfall events leads to more accurate results. Equations developed in this paper include the possibility of residual storage from more than one previous rainfall event.
To relate runoff probability and residual storage probability to design parameters (e.g., the growing medium layer thickness of green roofs, the storage volume of rainwater harvesting systems, the substrate layer depth of permeable pavements, etc.) makes the proposed equations useful in the design practice. Two case studies in Genova and Milano (Italy), characterized by similar climate conditions about average annual rainfall and temperature but different rainfall distribution throughout the year, were used as tests. The comparison between results obtained by applying the Fig. 1 Comparison between traditional urban drainage systems and SUDS derived formulas with the results obtained through the continuous simulation of observed data has confirmed the goodness of the proposed approach.

Methodology
Stormwater control through SUDS increases mainly three components of the urban water balance: infiltration, evapotranspiration, and rainwater harvesting and reuse. The three main categories of SUDS, identified in Fig. 1, are distinguished by the predominance of one of the said three components, each characterized by a specific outflow. In Fig. 2, for all the three systems, rainfall and runoff volumes are identified respectively with the parameters h and v.
The output volume for a unit of area (Q) was assumed: • Q = F (the infiltration volume) for infiltration systems, • Q = ET (evapotranspiration volume) for evapotranspiration systems, • Q = Y (yielded volume) for rainwater harvesting systems.
The dashed line means that the different types of SUDS can work alone and in combination to increase their potential and effectiveness. For example, rainwater harvesting systems can integrate green roofs known as Green Blue Roofs (Busker et al. 2022) or infiltration systems (Raimondi and Becciu 2014).
One of the model hypotheses is to assume a constant outflow rate (q) over the considered time interval. The infiltration rate is often assumed equal to the evapotranspiration rate at saturation. The hypothesis is precautionary for emptying time and runoff estimation, which are the main parameters to evaluate since infiltration systems must be empty quickly and limit floods.
The evapotranspiration rate is generally assumed equal to the actual monthly evapotranspiration rate. Using the potential evapotranspiration rate would be precautionary when vegetation survival without irrigation is considered but would lead to underestimating the runoff probability. For rainwater harvesting systems, the hypothesis of constant water supply-demand is usually suitable for regular uses (e.g., WC toilet flushing) or for a limited time scale (e.g., monthly irrigation demand).
The input variable to the analytical probabilistic model is rainfall. For the identification of independent events from the continuous record of rainfalls, a minimum interevent time known as Inter Event Time Definition (IETD) was defined (USEPA 1986). If the interevent time between two consecutive rainfall events is minor of IETD, the two rainfalls are joined together into a single event; on the contrary, they are considered independent.
Main rainfall characteristics, rainfall depth ( h) , rainfall duration ( ) and interevent time (d) were assumed to be independent and exponentially distributed. The effects of these two simplifying hypotheses were deepened and accepted by different studies in the literature (Eagleson 1978;Adams et al. 1986;Adams and Papa 2000;Bedient and Huber 1992).
The bias due to their use is negligible compared to the complexity reduction of the analytical derivation. Regarding the hypothesis of exponential distribution of rainfall variables, Bacchi et al. (2008) tested that for the most Italian catchments, the Weibull distribution better fits the frequency distribution of observed data. However, the benefits of its use in the analytical probabilistic model are not justifiable since it would make the integration more complex. Moreover, the double-exponential distribution fits the frequency distribution of observed data better than the exponential distribution; it is quite easy to integrate, but final expressions are longer and more complex. In addition, its use only brings little improvement in the accuracy of results (Becciu and Raimondi 2012). The exponential probability density functions of rainfall depth, rainfall duration, and interevent time are expressed by: (1) where: = 1∕ h , = 1∕ and = 1∕ d − IETD . h is the average rainfall depth, is the average rainfall duration and d is the average interevent time.
The estimate of the runoff and the residual storage probability distribution considers the possibility of consecutive chained rainfalls (Fig. 3).
In the following modeling, a generic SUDS of volume w with outflow rate q was considered.
The first step is the definition of the analytical equation expressing the water content w in the system at the end of the (i-1) rainfall event: The subscripts "u" and "e" refer respectively to the end and the beginning of the rainfall event. The subscripts "i" and "i-1" identify the order number of the event in the stochastic rainfall series, in the specific the (i-1) th and the i th events. For i = 1 , Eq. (4) results: w e,0 = 0 , the storage capacity was assumed empty at the beginning of the analysis.
At the beginning of the i th rainfall event, the water content in the storage volume can be assumed equal to: ( The runoff at the end of the i th event ( v i ), valid for i > 0, results: For Condition 1 : no runoff at the end of the ( i − 1) th event, residual storage at the beginning of the i th rainfall, and runoff at its end occur. For Condition 2 : no runoff at the end of the (i − 1) th event, no residual storage at the beginning of the i th event, and runoff at its end occur. For Condition 3 : runoff at the end of the ( i − 1) th event, no residual storage at the beginning of the i th event, and runoff at its end occur. For Condition 4 : runoff at the end of the (i − 1) th event, residual storage at the beginning of the i th event, and runoff at its end occur. The runoff volume for i = 0 ( v 0 ) is: Applying the probabilistic approach to the analytical Eqs. (7) and (8) allows the estimation of the runoff probability distribution P v . Two different conditions were distinguished: -w∕q ≤ IETD , independent rainfall events without residual storage from previous rainfalls (storage volume fully available at the beginning of the considered event); it is the case of small storage volumes with great outflow rates, allowing emptying time less than IETD. -w∕q > IETD , consecutive chained events with residual storage from previous rainfalls at the beginning of the  When considering IETD of a few hours, as often happens in small urban catchments where SUDS are suitable, there is a high probability of residual storage from previous rainfalls at the beginning of the considered one.
For the first condition ( w∕q ≤ IETD ), the runoff probability is: where = +q⋅ and v is the runoff threshold.
For the second condition ( w∕q > IETD ), the runoff probability is: N is the number of consecutive chained rainfall events. The probability of residual storage in the system at the end of the interevent time ( P r ), was estimated by applying the probabilistic approach to the analytical Eq. (6), which defines the water content at the beginning of the i th rainfall event. Also, in this case, two different conditions were considered: -N = 2 , a couple of consecutive chained rainfall events, -N > 2 , more than two consecutive chained events.
For the first condition ( N = 2 ), the residual storage probability from the first rainfall at the beginning of the second one ( P r ) results: When the second condition occurs (N > 2), the residual storage probability from previous rainfalls ( P r ) is expressed by: . Both Eqs. (9) and (10) for estimating the runoff probability and Eqs. (11) and (12) for estimating residual storage probability allow associating the design parameters to a probabilistic level (related to the ARI). Set the level equations are solved making explicit the selected design parameter (i.e., the design volume).
The flowchart in Fig. 4 summarizes the procedure for SUDS design through the analytical-probabilistic method.   . Despite these analogies, rainfall distribution throughout the year is different (Fig. 5). Dry summers and rainy autumns characterize the rainfall regime of Genova, where the average rainfall in October is five times greater than in July. In Milano, the dry season is winter, and the difference between the driest month (January) and the rainiest one (October) is less pronounced. The modeling considered the rainfall series recorded at the Milano-Monviso gauge station in 34 years  and those recorded at the Genova-Villa Cambiaso gauge station in 24 years . Genova is rainier than Milano, with a mean of 56 events per year compared to 32 events recorded at Milano.
In the analysis, the Initial Abstraction (IA) was assumed equal to 2 [mm], neglecting all rainfall depths lower than the threshold value in the series. A minimum dry period (IETD) equal to 6 [hours] was assumed, for both the case studies, to identify independent rainfalls from the continuous record. Table 1 reports the average values, per event, of the three rainfall variables used in the modeling (rainfall depth h, rainfall duration θ, and interevent time d).
Rainfall depth is higher in Genova than in Milano, whereas rainfall duration is very similar; this means higher average rainfall intensity in Genova than in Milano. Moreover, in Milano, the average interevent time is longer than in Genova, following the lower number of rainfalls per year. Table 2 reports the coefficient of variation (V) of the three rainfall variables for both cities to test the assumption of their exponential distribution.
The hypothesis of exponential distribution resulted suitable only for rainfall duration (for which the coefficient of variation tends to be one). It is an expected result; as discussed in "Methodology" section, the exponential distribution is preferable, despite other probability distributions better fitting the frequency distribution of observed data. It is easy to integrate, and its use limits the final expressions complexity, with no significant effects on the quality of results. Table 3 reports the correlation indexes among the rainfall variables to test the hypothesis of independence.  [-] 0.69 0.61 ρ d,θ [-] 0.10 0.01

Fig. 6
Runoff probability P v (w) and frequency F v (w), at varying the storage volume ( w ) for Genova (GE) and Milano (MI) Fig. 7 Residual storage probability P r (w) and frequency F r (w), at varying storage volumes (w) for Genova (GE) and Milano (MI) 1 3 The correlation between rainfall depth and interevent time and between interevent time and rainfall duration is negligible, whereas the correlation between rainfall depth and duration is significant. Copula functions were introduced in the last decades in the hydrologic research to overcome the correlation among rainfall variables (Abdollahi et al. 2019); however, they are used in this study since they involve an increase in complexity that is not justifiable by a significant improvement in results. The outflow rate was set equal to q = 0.36 mm∕hr = 1[l∕ s ⋅ ha imp ] . The runoff threshold v in Eqs. (9) and (10) and the minimum water content w in Eqs. (11) and (12) were assumed equal to zero. The storage capacity w was varied between 0 and 250 [mm]. Figure 6 compares the runoff probability distribution P v (w), obtained by applying the proposed probabilistic Eqs. (9) and (10), with the runoff frequency distribution F v (w ), obtained through the continuous simulation of the observed data Eq. (7). The term P v (w) corresponds to: -P v1 (w) , for w∕q ≤ IETD (runoff considering a single event without residual storage from previous rainfalls). -P vN (w) for w∕q > IETD (consecutive chained events with residual storage).
The analysis was carried out for the considered cities, Milano (MI) and Genova (GE).
For both case studies, the results of the proposed equations show a good agreement with those of the continuous simulation of observed data. The runoff probability decreases as the storage volume increases and is higher for Genova than for Milano. It is explained by the higher average rainfall depth and the lower interevent time characterizing the rainfall regime of Genova, which also involves a higher number of consecutive chained events ( N = 4 for Genova and N = 2 for Milano). Figure 7 compares the residual storage probability distribution P r (w ), obtained by applying the proposed probabilistic Eqs. (11) and (12) to the residual storage frequency distribution F r (w ), obtained through the continuous simulation of the observed data Eq. (6). The term P r (w) corresponds to: -P r2 (w) , for N = 2 , if only a couple of rainfall events are considered. -P rN (w) , for N > 2 , if more than two consecutive chained rainfall events are considered.
The analysis considered both Milano (MI) and Genova (GE). The number of consecutive chained rainfall events was still N = 4 for Genova and N = 2 for Milano.
The residual storage probability is higher for Genova than Milano because of the higher average rainfall depth and the lower average interevent time characterizing the rainfall regime of Genova. The proposed equations for estimating residual storage probability well fit the frequency distribution. In both cases, the storage volume slightly affects the residual storage probability.
Table 4 (referred to Genova) and Table 5 (referred to Milano) report the values of the storage volume associated with four different runoff ARIs (T v ), usually used in practice: T v = 5-10-20-50 [years].
The storage volume (w) associated with the runoff ARI (T v ), related to runoff probability by the equation T v = 1∕P v (w ), was estimated by Eqs. (9) and (10). It increases as the ARI (T v ) grows. To achieve the same runoff ARI (T v ), the rainfall regime of Genova requires larger storage volumes than Milano. The ratio (R) between storage volumes required in Genova and Milano, with the same ARI, is reported in Table 6.
The ratio (R) increases as ARI (T v ) grows.

Conclusions
This paper proposes an analytical-probabilistic approach for estimating runoff and residual storage probability applicable to different types of SUDS. The system type and principle (infiltration, evapotranspiration, or rainwater harvesting) isn't significant in the modeling scheme aiming to simulate the storage process. In SUDS design is fundamental to consider both runoff and residual storage probability. These two factors are related to the system failure in terms of flooding occurrence, vegetation survival, and water demand fulfillment. The derived equations allow taking into account the possibility of consecutive chained rainfall events, which is a meaningful feature for systems with low flow release, as SUDS often are. The application to case studies, characterized by two different rainfall regimes, confirmed the goodness of the suggested method. Comparison with results from continuous simulation confirmed the accuracy and reliability of developed equations. A clear advantage of using these equations is the possibility of achieving the accuracy of the continuous simulation with the simplicity of typical design storm methods. Finally, the type of hydrological information needed for applying these equations is, in most cases, easier and cheaper to acquire than long series of rainfall records, not always available.

Symbols
The following symbols are used in this paper (*volume per unit area): Q: Output volume*; F: Infiltration volume*; ET: Evapotranspiration volume*; Y: Volume yielded from the rainwater harvesting system*; h: Rainfall depth; v: Runoff volume*; q: Outflow rate; θ: Rainfall duration; d: Interevent time; μ h : Average rainfall depth; μ θ : Average rainfall duration; μ d : Average interevent time; w: Storage volume*; w e,i : Water content at the beginning of the i th rainfall event; w u,i : Water content at the end of the i th rainfall event; v : Runoff volume threshold*; N : Number of consecutive chained rainfall events; w : Water content threshold; V h : Coefficients of variation of rainfall depth; V θ : Coefficients of variation of rainfall duration; V d : Coefficients of variation of interevent time; ρ h,d : Correlation index among rainfall depth and interevent time; ρ θ,h : Correlation index among rainfall duration and rainfall depth; ρ d,θ : Correlation index among interevent time and rainfall duration; P v : Runoff probability; P v2 : Runoff probability considering a couple of rainfall events; P vN : Runoff probability considering a chain of N consecutive rainfall events; P r : Residual storage probability; P r2 : Residual storage probability considering a couple of rainfall events; P rN : Residual storage probability considering a chain of N consecutive rainfall events;