Friction decoupling and loss of rotational invariance in flooding models
Luca Cozzolino, Giada Varra, Luigi Cimorelli, Domenico Pianese, Renata Della Morte
This work is licensed under the Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License. To view a copy of this license, visit http://creativecommons.org/licenses/by-nc-nd/4.0/ or send a letter to Creative Commons, PO Box 1866, Mountain View, CA 94042, USA.
Friction decoupling and loss of rotational invariance in flooding models
Luca Cozzolino , Giada Varra , Luigi Cimorelli , Domenico Pianese , Renata Della Morte Abstract
Friction decoupling, i.e. the computation of friction vector components making separate use of the corresponding velocity components, is common in staggered grid models of the SWE simplifications (Zero-Inertia and Local Inertia Approximation), due to the programming simplicity and to the consequent calculations speed-up. In the present paper, the effect of friction decoupling has been studied from the theoretical and numerical point of view. First, it has been found that friction vector decoupling causes the reduction of the computed friction force and the rotation of the friction force vector. Second, it has been demonstrated that decoupled-friction models lack of rotational invariance, i.e. model results depend on the alignment of the reference framework. These theoretical results have been confirmed by means of numerical experiments. On this basis, it is evident that the decoupling of the friction vector causes a major loss of credibility of the corresponding mathematical and numerical models. Despite the modest speed-up of decoupled-friction computations, classic coupled-friction models should be preferred in every case. Ass. Prof., Ph. D., Dept. of Engrg., Parthenope Univ., Centro Direzionale di Napoli – Is. C4, 80143 Napoli, Italy. E-mail: [email protected] Ph. D. Stud., Dept. of Engrg., Parthenope Univ., Centro Direzionale di Napoli – Is. C4, 80143 Napoli, Italy. E-mail: [email protected] Sen. Res., Ph. D., Dept. of Civil, Architectural and Environmental Engrg., Federico II Univ., via Claudio 21, 80125 Napoli, Italy. E-mail: [email protected] Full Prof., Dept. of Civil, Architectural and Environmental Engrg., Federico II Univ., via Claudio 21, 80125 Napoli, Italy. E-mail: [email protected] Full Prof., Ph. D., Dept. of Engrg., Parthenope Univ., Centro Direzionale di Napoli – Is. C4, 80143 Napoli, Italy. E-mail: [email protected]
Subject headings : Friction; Flooding; Shallow-water Equations; Zero-Inertia model; Local-Inertia model; Rotational invariance.
Corresponding author: Luca Cozzolino E-mail address: [email protected] Address: Dipartimento di Ingegneria, Università degli Studi di Napoli Parthenope, Isola C4, 80143 Napoli (Italy) Phone: +390815476723. . Introduction
The two-dimensional Shallow water Equations (SWE) are a very popular mathematical model for flooding simulations in urbanized and non-urbanized areas, and a variety of numerical methods are now available to tackle their approximate solution. Among the others, Finite Volume schemes with interface Riemann solvers (Toro 2001) have been successful for their ability to manage flow discontinuities such as shocks and wet-dry frontiers. In addition, smart numerical approaches in Finite Volume schemes allow satisfying properties such as depth positivity, well-balancing, non-linear stability (Bouchut 2004). Despite their success, these methods are computationally expensive, and researchers are currently making efforts along different research directions to speed up computations. Among the others, it is possible to cite research fields such as parallel computing (Delis and Mathioudakis 2009), Porous Shallow water Equations in urbanized areas (Varra et al. 2020), Finite Difference schemes avoiding the approximation of a Riemann problem (Stelling and Duinmeijer 2003), and the adoption of simplified forms of the SWE model such as the Local Inertia Approximation (LInA, Bates et al. 2010) and the Zero-Inertia approximation (Xanthopoulos and Koutitas 1976). Simplified versions of the SWE model are obtained by canceling the equation terms that are negligible with respect to the others in appropriate conditions. For example, the LInA model arises from cancelation of the convective acceleration terms (de Almeida and Bates 2013), while the Zero-Inertia model is obtained by neglecting both local and convective accelerations (Xanthopoulos and Koutitas 1976). However, the friction forces acting at the flow bed are important in most applications for physical credibility of the models, and for this reason they are always retained. In two-dimensional flooding models, the friction forces per unit of weight are represented by the dimensionless friction vector , ,
Tf f x f y
S S S , where S f,x and S f,y are the component of S f along the horizontal coordinate axes x and y , respectively, while T is the matrix transpose symbol. If bed oughness is isotropic, the friction vector is parallel to the flow velocity, implying that S f,x and S f,y are coupled because they both depend on the components of the velocity along x and y (Viero and Valipour 2017). Due to the availability of raster digital elevation models as topographic input, the LInA and the Zero-Inertia models are often solved on square staggered grids (Bates and De Roo 2000). The two horizontal components of velocity are evaluated in different points of the grid, implying that the evaluation of friction is accomplished by means of some form of averaging on a large computational stencil (Hromadka and Chintu Lai 1985, Galbiati and Savi 1995, Wasantha Lal 1998, Bradbrook et al. 2004, Panday and Huyakorn 2004, Yu and Lane 2006, de Almeida and Bates 2013, Leandro et al. 2014). Alternative to the averaging, acceleration of staggered grid algorithms can be obtained by decoupling the components of S f , i.e. computing S f,x and S f,y by separately using the x and y components of the velocity. Until now, the decoupled friction approach has been favored by many Authors (Bates and De Roo 2000, Horritt and Bates 2001, 2002, Hunter et al. 2005, 2007, Fewtrell et al. 2008, Dottori and Todini 2011, Neal et al. 2011, 2012, 2018, de Almeida et al. 2012, Falter et al. 2013, Yamazaki et al. 2013, Coulthard et al. 2013, Mendicino et al. 2015, Savage et al. 2016, Nguyen et al. 2016, Mateo et al. 2017, Sridharan et al. 2020) because it is easily programmed, more computationally efficient (Hunter et al. 2007), and numerical results have been found satisfactory when the topographic control is dominant (Bates and De Roo 2000, Horritt and Bates 2001). Nonetheless, it has been recognized that the decoupled-friction approach does not faithfully reflect the structure of the models derived from first principles (Hunter et al. 2007), and issues arise in urbanized areas (Hunter et al. 2007). It is believed that the calibration of friction coefficients may reduce the discrepancies between coupled- and decoupled-friction models (Horritt and Bates 2001), but the question if the errors introduced in decoupled-friction models can be compensated by calculations speed-up remains still open. ue to the lack of consensus about the viability of decoupled-friction approach, in the present paper the question is tackled from the theoretical and numerical point of view. In particular, the coupled and decoupled definitions of the friction vector are compared, and the influence of these definitions on the mathematical structure of flooding models is explored, demonstrating that the decoupled friction definition causes the loss of rotational invariance of the models. This theoretical finding allows the interpretation of the numerical results, which confirm the lack of rotational invariance of decoupled-friction models. The paper is organized as follows: in Section 2, the friction vector is considered from a theoretical point of view; in Section 3, the LInA and the Zero-Inertia models with coupled and decoupled friction are presented, and four numerical schemes from the literature are used to solve a case study with radial symmetry; the theoretical and numerical findings are discussed in Section 4 while, finally, conclusions are drawn in Section 5.
2. Analysis of the friction term
If the friction force is parallel to the flow velocity, it is possible to write (1) , , , yxf x f f y f qqS S S Sq q . In Eq. (1), the meaning of the symbols is as follows: x q hu and y q hv are the components along x and y , respectively, of the unit-width discharge vector Tx y q q q ; h is the flow depth; u and v are the component along x and y , respectively, of the depth-averaged flow velocity; f f x f x S S S is the modulus of S f and x y q q q is the modulus of q . Despite its empirical nature, the Manning’s formula has been justified theoretically (Gioia and Bombardelli 2002) and it is valid in a large variety of circumstances, including vegetated beds with high submergence (Gualtieri et al. 2018). For this reason, the expression (2)
22 10 3 2 2 2 10 3 f M M x y n h n q q h
S q is usually adopted for the evaluation of the friction vector modulus. Eq. (2), where n M is the Manning’s friction coefficient, implies that the components of the friction vector are (3) , M x x y M y x yf x f y n q q q n q q qS Sh h . The inspection of Eq. (3) clearly shows that the components of the friction vector are coupled because they simultaneously depend on q x and q y . In the literature of simplified flooding models, it has become popular the definition of a decoupled friction vector , , ' ' ' Tf f x f y
S S S , where the components are calculated separately in a strictly one-dimensional fashion. In the case of the Manning’s formula, the corresponding components are (4)
22, ,10 3 10 3 ' , '
M y yM x xf x f y n q qn q qS Sh h . Notably, the decoupled friction vector ' f S is not parallel to the flow vector q , because , ' f x S is not proportional to q x but to x x q q and , ' f y S is not proportional to q y but to y y q q . In the present section, the coupled and decoupled definitions of the friction vector will be compared, and the influence of these definitions on the mathematical structure of flooding models will be explored. The modulus of the decoupled friction vector is (5) ' ' '
M x yM yM xf f x f x n q qn qn qS S h h h S , which differs from that of coupled friction vector in Eq. (2). If is the angle between the friction vector S f and the x axis, the corresponding tangent is , , tan f y f x y x S S q q . This allows calculating the ratio ' f f S S between the moduli of the vectors S ’ f and S f as (6) ' tan 1 tantan 1 tan f x y x xx y x xf q q q qq q q q SS . The ratio ' f f S S ranges in the interval [ , 1] for [0, 2 ] (see Figure 1a), and this means that the friction force calculated by means of the decoupled formulas is systematically smaller than the coupled friction vector forces. The maximum error is attained for angles = (2 n + 1) /4 with n = 0, 1, 2, 3, i.e. when the vector S f is aligned with the bisectors of the reference framework quadrants. The average value of ' f f S S is 7) ' 1 tan 1 0.862 tan 1 ff E d SS , showing that the average error between the vector moduli is of about 14%. If is the angle between the friction vector S ’ f and the x axis, the corresponding tangent is (8) ,, 'tan tan tan' y yf yf x x x q qSS q q , which differs from tan . This means that the decoupled friction vector is rotated with respect to the coupled friction vector by the angle - . The absolute deviation ranges in the interval [0, 0.071 ] (see Figure 1b), and the corresponding average value is (9) E d , which is small and corresponding to 8.28° degrees. From this discussion, the apparent conclusions are that: - the modulus of S ’ f is smaller than the modulus of S f , and this may require some calibration of the friction coefficient; - the vectors S f and S ’ f , which coincide only for n with n = 0, 1, 2, 3, are quite similar in the average; - the substitution of S f with S ’ f is safe. These provisional conclusions seem to justify the use of the decoupled friction S ’ f instead of the friction S f in flooding mathematical and numerical models. In the next sub-section, it will be hown that the substitution of S f with S ’ f induces unexpected structural changes in flooding models that are strongly in contrast with the preceding observations. [Insert Figure 1 about here] In Newton’s mechanics, the galilean transformations (rotation of the reference framework by an angle , translation with uniform velocity, shifting of the origin) preserve intervals of time and distance between simultaneous events (Arnold 1989). This implies that flooding models are congruent with Newton’s mechanics only if their structure remains invariant under rotation of the reference framework. While SWE simplifications such as LInA and Zero-Inertia models make often use of the decoupled friction vector, this never happens in the literature of SWE model itself. Nonetheless, the preservation of the rotational invariance for the SWE model only is considered here, because the corresponding LInA and Zero-Inertia simplifications trivially inherit this property. Under the assumption that Eqs. (1) and (2) are used for the friction vector definition , the SWE can be written as (10)
222 7 32 22 7 3
02 2 yx x y M xx xy x y y M y qqht x y q q n qq qgh zgh gt x h y h x hq q q q n qgh zgh gt x h y h y h qq , here t is the time, h is the flow depth, z is the bed elevation, and g = 9.81 m/s is the gravity acceleration, while all the other symbols have been already defined. The first of Eq. (10) represents the mass conservation, while the second and the third represent the momentum conservation along x and y , respectively. In the following, the invariance of Eq. (10) to the frame reference rotation will be easily verified. Let OXY be a reference framework rotated with respect to the reference framework
Oxy by the angle . The relation between the coordinates in the two reference frameworks is: (11) cos sin , sin cos X x y Y x y . If
TX Y
Q Q Q is the unit-width discharge vector in the rotated framework, the relation between Q and q is (12) cos sinsin cos x X Yy X Y q Q Qq Q Q . From Eq. (12) it is evident that x y X Y q q Q Q q Q , confirming that the modulus of the velocity vector remains unchanged under a rotation of the reference framework. Recalling that the relation between the derivatives in the two frameworks is (13) cos sinsin cos x X Yy X Y , the Eq. (10) can be rewritten in the rotated framework (see Appendix A) as: (14)
222 7 3222 7 3
02 2
X Y M XX X X Y M YY X Y Y
Q Qht X Y n QQ Q Q Qgh zgh gt X h Y h X hn QQ Q Q Qgh zgh gt X h Y h Y h QQ . The comparison between Eq. (10) and Eq. (14) shows that the structure of SWE model with coupled friction remains unchanged under the rotation framework, because the rotated systems of Eq (14) can be obtained from Eq. (10) by simply substituting X to x , Y to y , Q X to q x , and Q Y to q y . Interestingly, the reference framework rotation preserves the structure of the coupled friction terms because the modulus of the friction vector and its parallelism with the flow velocity vector are preserved. These findings imply that the computation results in the rotated framework are indistinguishable from the computation results in the original reference framework, congruently with the Newton’s mechanics and with our physical intuition. In other words, all the reference frameworks are equally good for the simulation of floodings with coupled-friction models. Is the preceding result applicable to the case of decoupled friction? The answer is no. With respect to the Oxy reference framework, the SWE model with decoupled friction can be written as (15)
222 7 3222 7 3
02 2 yx x y M x xx x M y yy x y y qqht x y q q n q qq qgh zgh gt x h y h x hn q qq q q qgh zgh gt x h y h y h . The rotation of the reference framework by the angle supplies the rotated system 16)
222 1 7 3222 2 7 3
X Y M X YX X X Y M X YY X Y Y
Q Qht X Y n f Q QQ Q Q Qgh zgh gt X h Y h X hn f Q QQ Q Q Qgh zgh gt X h Y h Y h , where the functions , , X Y f Q Q and , , X Y f Q Q are defined as (17) , , cos cos sin cos sinsin sin cos sin cos, , sin cos sin cos sincos sin cos sin cos X Y X Y X YX Y X YX X Y X YX Y X Y f Q Q Q Q Q QQ Q Q Qf Q V Q Q Q QQ Q Q Q . The comparison between Eqs. (15) and (16) shows that the rotation of the reference framework causes the structural change of the decoupled-friction SWE model. In particular, the two components of the friction in the rotated framework are not proportional to
X X
Q Q and
Y Y
Q Q while Eq. (16) reduces to Eq. (15) only for n with n = 0, 1, 2, 3. The lack of rotational invariance implies that the results of the SWE model with decoupled friction depend on the choice of the reference framework . This characteristic is in contrast not only with Newton’s mechanics but also with the physical intuition, and it represents a major loss of meaning of the model. In this case, it is evident that calculations are physically meaningful only if the flow is strictly one-dimensional and parallel to one of the coordinate axes. t is possible to demonstrate that the LInA and the Zero-Inertia models inherit the rotational invariance property of the SWE model. In particular, the LInA and the Zero-Inertia models with coupled friction terms preserve the rotational invariance, while the corresponding decoupled-friction models do not satisfy this fundamental property. This is easily proven by applying the procedure of Sub-section 2.2 to the simplified models with coupled and decoupled friction, but the complete demonstration is omitted here for the sake of brevity. It follows that the results of LInA and Zero-Inertia models with coupled friction do not depend on the reference framework, while the same models with decoupled friction suffer from dependency on the choice of the reference framework . The computational aspects of this discouraging conclusion will be explored in the next Section.
3. Flooding computation with coupled/decoupled friction
In the present Section, the LInA and the Zero-Inertia mathematical models are presented and, for each friction version (coupled/decoupled) of these models, a numerical algorithm from the literature (Wasantha Lal 1998, Bates and De Roo 2000, de Almeida et al. 2012, de Almeida and Bates 2013) is briefly described. Finally, the results corresponding to a special case study are presented for each numerical model. For the sake of simplicity, details about the treatment of wetting-drying fronts are omitted because the case study refers to a fully wet flow field.
In the present Sub-section, the mathematical and numerical LInA and Zero-Inertia models are presented. Despite the controversy about viability of the LInA model in flooding applications Cozzolino et al. 2019), this model is considered here because it is often applied after the decoupling of friction vector components.
The LInA model is obtained from the SWE model by neglecting the convective acceleration terms. Assuming that discontinuities are absent in the flow field, the LInA model with coupled friction is written as (de Almeida and Bates 2013) (18) yx M xxy M y qqht x y n qq gh gt x hq n qgh gt y h qq , where h z is the free-surface elevation. For the approximate solution of Eq. (18), the numerical scheme on staggered grid by de Almeida and Bates (2013) is considered here. After that the physical domain is subdivided in square cells with sides of length s that are aligned with the x and y axes, the numerical approximation of h is collocated at the center of the cell (black dots in Figure 2). Correspondingly, the numerical approximations of q x and q y are collocated at the center of interfaces which are normal to the axis x (grey squares in Figure 2) and to the axis y (grey triangles), respectively. The component x q of the unit-width discharge at the interface ( i -1/2, j ) between the cells ( i -1, j ) and ( i , j ) is updated by means of (19) , 1 2, , 1 2, , 3 2, 1 2, , 1,1, 1 2, 7 32 1 2, 1 2, n n n n n nx i j x i j x i j i j i j i jnx i j n nM i j i j q q q gh t sq g tn h q . In Eq. (19), the meaning of the symbols is as follows: , 1 2, nx i j q is a numerical approximation of x q at the interface ( i -1/2, j ) at the time level t n ; is a relaxation factor chosen in the interval [0, 1]; n n t t t is the time step; , , , n ni j i j i j h z is a numerical approximation of the free-surface elevation in the cell ( i , j ), where , ni j h and , i j z are the cell-averaged flow depth and bed elevation, respectively; ni j h is a numerical approximation of h at the interface ( i -1/2, j ); finally, n n ni j x i j y i j q q q is an approximation of q at the interface ( i -1/2, j ). Notably, the numerical counterparts of h and q y are not defined at the interface ( i -1/2, j ), which is parallel to the y axis. For this reason, the approximations (20) max , max , n n ni j i j i j i j i j h z z and the average (21) , 1 2, , 1, 1 2 , 1, 1 2 , , 1 2 , , 1 2 n n n n ny i j y i j y i j y i j y i j q q q q q are used. A similar formulation, with obvious changes, is used for updating the unit-width discharge q y at cell interfaces. Once that the discharges are updated, it is possible to update the cell-averaged flow depth , ni j h by means of (22) n n n n n ni j i j x i j x i j y i j y i j t th h q q q qs s . [Insert Figure 2 about here] Note that the calculation of
1, 1 2, nx i j q involves only , 1 2, nx i j q when = 1, while , 3 2, nx i j q and , 1 2, nx i j q are also involved when < 1. In other words, the numerical stencil for the calculation of
1, 1 2, nx i j q is expanded in the x direction when < 1. The LInA model with decoupled friction can be written as (de Almeida et al. 2012) (23) yx M x xx M y yy qqht x y n q qq gh gt x hn q qq gh gt y h . The numerical scheme for its solution (de Almeida et al. 2012) is identical to the scheme used for the case of coupled friction, but now q x is updated without considering the contribution of q y in the friction calculation, and q y is updated neglecting the contribution of q x . For this reason, Eq. (19) is substituted by (24) , 1 2, , 1 2, , 3 2, 1 2, , 1,1, 1 2, 7 32 , 1 2, 1 2, n n n n n nx i j x i j x i j i j i j i jnx i j n nM x i j i j q q q gh t sq g tn q h , while Eq. (21) is not used. Due to the smaller number of operations, the algorithm with decoupled friction is faster than the algorithm with coupled friction. The Zero-Inertia model is obtained from the SWE model by neglecting both the local and the convective acceleration terms. The model version with coupled friction can be written as (Xanthopoulos and Koutitas 1976) (25) yxM xM y qqht x yn qx hn qy h qq . For computational purposes, it is convenient to rewrite the system in the form (Hromadka and Chintu Lai 1985) (26) yxx My M qqht x yhq xn hq yn SS where T x y S is the free-surface gradient and x y S is the free-surface gradient modulus . Following Wasantha Lal (1998), the component x q of the unit-width discharge at the interface ( i -1/2, j ) between the cells ( i -1, j ) and ( i , j ) at the time level t n is calculated by means of 27) n n ni j i j i jnx i j nM i j hq sn S , where (28) n n n n n ni j i j i j i j i j i jni j s s s S is a numerical approximation of S at the interface ( i -1/2, j ), while (29)
1, ,1 2, n ni j i jni j h hh is a numerical approximation of h at the same interface. Similar formulas are used for the calculation of the unit-width discharges on the other interfaces, allowing the update of h in the cell ( i , j ) by means of the explicit formula (29)
1, , , 1 2, , 1 2, , , 1 2 , , 1 2 n n n n n ni j i j x i j x i j y i j y i j t th h q q q qs s . The Zero-Inertia model with decoupled friction can be written as 30) yxM x xM y y qqht x yn q qx hn q qy h , which can be rewritten as (Hunter et al. 2005) (31) yxx My M qqht x yhq n x xhq n y y , where sgn ( ) represents the signum function. Following Bates and De Roo (2000), the component x q of the unit-width discharge at the interface ( i -1/2, j ) between the cells ( i -1, j ) and ( i , j ) at the time level t n is calculated by means of (31) sgn n n n n ni j i j i j i j i jnx i j M hq n s s , where Eq. (20) is used to approximate the flow depth ni j h at the interface ( i -1/2, j ). Similar formulas, with obvious changes, are used for the calculation of the unit-width discharges on the interfaces with normal y , allowing the update of h in the cell ( i , j ) by means of Eq. (29). he numerical scheme with decoupled friction is faster than the scheme with coupled friction because the approximation of the complete free-surface gradient at cell interfaces is avoided. The mathematical and numerical models of Sub-section 3.1 are evaluated considering a case study inspired by that introduced in Wasantha Lal (1998). The physical domain consists of a square area with side L = 160930 m characterised by uniform bed elevation ( z = 0 m) and uniform Manning’s roughness coefficient ( n M = 0.03 s m -1/3 ). The origin O of the reference framework Oxy is at the centre of the domain, with axes x and y parallel to the square sides. The initial conditions are characterised by sinusoidal free-surface elevation (32) x y x yx y x y , while wall boundary conditions are considered at the limits of the domain. The expected solution of the problem remains radially symmetric with time and the fronts of the propagating wave are expected to keep a circular shape for t > 0, because the initial conditions are radially symmetric while roughness and bed elevation are uniform. Due to the symmetry of initial conditions with respect to the x and y axes, it is sufficient to consider the square numerical domain with side l = L /2 = 80465 m for x y
0, with origin O of the reference framework located at the lower left vertex of the domain (see for example Figure 3). The computer programs for the numerical models compared here are written in Microsoft Visual asic 6.0, then compiled and run on a laptop computer with Vista Home Basic operative system, Intel T2080 processor (1.73 GHz clock), and 2.00 Gb RAM. The results of the Zero-Inertia model with coupled friction (numerical scheme of sub-section 3.1.3) are presented at time t = 60000 s in the left column of Figure 3 for three different densities of the numerical grid. The corresponding uniform grid spacing and constant time steps used in the computations are detailed in Table 1. The free-surface elevation for the coarse grid computations is represented in the upper panel (Figure 3a). The inspection of this panel shows that, congruently with expectation, the free-surface elevation propagates symmetrically along all the directions originating from O , exhibiting good isotropy properties. The grid refinement in Figure 3c (medium grid) and Figure 3e (refined grid) confirms the radial symmetry of the solution, congruently with the physical intuition. The results of the corresponding decoupled friction model (numerical scheme of sub-section 3.1.4) are presented at time t = 60000 s in the right column of Figure 3. A closer inspection of Figure 3b shows that the effect of friction decoupling is the distortion of the propagating wave, which progressively tends to a sort of smoothed square, because the free-surface contour lines tend to be parallel to the coordinate axes. This implies that a sort of preferential line of wave propagation arises along the bisector of the numerical domain, where the flow is faster. Also in the case of Zero-Inertia model with decoupled friction, the distorted solution does not change significantly with the refinement of the grid in Figure 3d and 3f, meaning that the structural issue of the decoupled-friction model cannot be solved by refining the grid. By convention, we set the limit of the expanding wave at the free-surface elevation = 0.31 m, which is also the contour line with lower value plotted in Figure 3. The computations show that the area with km , while the corresponding area supplied by the decoupled friction model is 2.25 km . The ncrease of the area occupied by the expanding wave is of about 6.2 %, implying that the decoupled-friction wave is faster than the coupled-friction one. [Insert Figure 3 about here] [Insert Table 1 about here] The LInA model with coupled friction (numerical scheme of sub-section 3.1.1) is run with relaxation parameter = 1.0, using the time step and grid spacing of Table 1 to allow a fair comparison with the Zero-Inertia model. The initial velocity must be specified, and its value has been set to zero everywhere for the sake of simplicity. The corresponding results at time t = 60000 s are presented in the left column of Figure 4. The inspection of Figure 4a (coarse grid), Figure 4c (medium grid), and Figure 4e (refined grid), confirms that the LInA model with coupled friction satisfies the radial symmetry of the problem, similarly to the case of the Zero-Inertia model with coupled friction. The calculations are repeated with the decoupled-friction LInA model (numerical scheme of sub-section 3.1.2) with relaxation parameter = 1.0, and the corresponding results are represented in the right column of Figure 4 for the grids of Table 1. Again, the use of a decoupled-friction approach distorts the propagating wave fronts, which assume a smoothed square shape, while a preferential propagation line is individuated along the bisector of the numerical domain. In conclusion, the LInA model with decoupled friction exhibits the same anisotropy phenomenon of the Zero-Inertia model with decoupled friction, and this structural issue cannot be solved by means of grid refinement. The area with km , while the corresponding area supplied by the decoupled friction model is 2.26 m . Also for the LInA model, the decoupled-friction wave is faster than the coupled-friction one, because the increase of the area occupied by the expanding wave is of about 5.6%. [Insert Figure 4 about here] The use of a constant time step allows the comparison of the algorithm run times, which are reported in Table 2 for the computations on the medium size grid. It is evident that these run times are only indicative, and may vary depending on the case study at hand, on the optimization of the algorithm construction, on the programming language, and on the hardware used. The inspection of Table 2 first shows that, for given constant time step, the LInA model is computationally more expensive than the Zero-Inertia model. This is expected because the quantities h , q x , and q y , are advanced in time in the LInA model, while only h evolves in the Zero-Inertia model. The comparison between coupled and decoupled models show that the Zero-Inertia model computations with decoupled friction allows a 1.25 times speedup over the corresponding coupled-friction computations, while the LInA model with decoupled friction allows a 1.17 times speedup over the LInA model with coupled friction. It is possible to compare these results with the speedups obtained in the literature by means of parallel computing. Leandro et al. (2014) obtained up to a 1.7 times speedup for the Zero-Inertia model on 12 parallel cores, while Artichowicz and Gaşiorowski (19) obtained up to 4.86 times speedup on a 8 cores architecture. Similarly, Neal et al. (2018) reported up to 60 times speed up for an optimized version of the LInA model on 16 cores. It is evident the speedup obtained by decoupling the friction vector is modest with respect to the speed up obtained by parallel computing. [Insert Table 2 about here]
4. Discussion
In the present Section, the theoretical results of Section 2 are discussed considering the numerical results of Section 3. An additional remark is made about solution distortion caused by the decoupling introduced by numerics.
It is commonly believed (Horritt and Bates 2001) that the Manning’s friction coefficient can be adjusted to improve the results of decoupled-friction models. In order to understand the limits of this procedure, a calibration exercise has been accomplished by searching the value of n M that makes the decoupled-friction LInA model supply a moving wave area equal to the one supplied by the coupled-friction model at time t = 60000 s on the refined grid. After some attempt, this value has been found to be n M = 0.0325 s m -1/3 , which differs by 8% with the original value n M = 0.03 s m -1/3 . The results are compared in Figure 5, where the contour corresponding to = 0.31 m for the coupled-friction model is plotted with a thin black line, while the same contour for the decoupled-friction model is plotted with a thick dashed line. The inspection of the figure clearly shows that, despite the two curves enclose the same area, the shape of the decoupled-friction moving wave remains different from the one of the coupled-friction wave. In other words, friction calibration cannot reduce the systematic difference between wave propagation celerities on x / y axes and quadrant bisector. [Insert Figure 5 about here] he influence of the reference framework rotation on the model results is not an academic question and it has practical relevance. Alignment of the rectangular grid (i.e. of the reference framework) with the main direction of flow propagation is a common modelling practice (Bradbrook et al. 2004), because this procedure simplifies the setting of inflow and outflow boundary conditions. The connection between the grid alignment and the loss of rotational invariance studied in Sub-section 2.2 can be elucidated by considering Figure 6. In the left panel (Figure 6a), the solution of the LInA model with decoupled friction on refined grid is represented considering the reference framework Oxy where the x axis is aligned along the direction West-East (W-E) and the y axis is aligned along the direction South-North (S-N). In this case, the stretching direction of the solution (dashed line) is parallel to the direction from South-West to North-East. Consider now the reference framework OXY (Figure 6b), which is rotated by the angle = /4 with respect to the reference framework Oxy . In this case, the computations are characterised by stretching direction aligned with the South-North direction, confirming that the loss of rotational invariance causes the dependence of the solution on the alignment of the reference framework. At this point, the model user may ask which is the correct computational grid alignment, and the answer is that all the grid alignments are equally wrong in decoupled-friction models. Flooding models, which incorporate a certain amount of theories, aim not only to computation but also to explanation (Clark et al. 2016). In particular, it is expected not only that a flooding model supplies acceptable results that are validated by available field data, but also that it supplies correct results for the correct reason (congruence with underlying physical laws). Satisfaction of this requirement supports the credibility of the model also in circumstances that are different from those that were chosen for validation. Unfortunately, the findings above are particularly discouraging because the friction decoupling introduces a loss of physical meaning that is hardly recognizable in topography-dominated applications (Bates and De Roo 2000, Horritt and Bates 2001) but it becomes evident in flat areas with truly two-dimensional applications. From the modeller’s point of view, the decrease of friction forces and the dependence on the reference framework may be acknowledged, but they re hardly quantified a priori because comparison with a coupled-friction model is needed in each practical application. On the other hand, the common user may be unaware of an error that increases the model epistemic uncertainty (Oberkampf et al. 2002). It is evident that the conceptual and computational error, which in the writers’ opinion is not compensated by the modest speed-up of calculations, can be eliminated by adopting a coupled-friction model. It must be noted that these adverse effects did not emerge in the literature of decoupled-friction flooding models until now, and this can be easily justified considering that decoupled-friction flooding models have been often validated first by means of some purely one-dimensional case study in controlled conditions and then by comparing the model results with field data available from historical floods (de Almeida et al. 2012, Sridharan et al. 2020). Of course, one-dimensional case studies are not able to shed light on two-dimensional modelling intricacies, while real-world case studies introduce interpretation difficulties originated by complicated topographies and/or land cover characteristics. This explains why a truly two-dimensional flow propagation problem characterized by uniform bed elevation and uniform friction coefficient was chosen in Section 3 to evaluate the effects introduced by friction vector decoupling . [Insert Figure 6 about here]
In the recent past, a class of algorithms for flooding simulation termed Cellular Automata (Thomas and Nicholas 2002, Rinaldi et al. 2007, Ghimire et al. 2013, Liu et al. 2015) has emerged for flood delineation applications. In these algorithms, the flow field is represented by means of a two-dimensional lattice, and the evolution of the flow storage in each lattice cell is driven by the state of the surrounding elements and by an appropriate rule of interaction between cells (Cai et al. 2014). he rule for the evolution of the cell storage is often simple and empirical but, in many cases, a more physical approach can be used to calculate flows between cells. An important adverse characteristic of Cellular Automaton algorithms with simple evolution rules is the fact that their results are strongly affected by the geometry of the underlying grid (Schönfisch 1995), because plane front waves appear where isotropic propagation is expected. For example, the binary spread on a rectangular grid with Neumann neighbourhood produces a diamond shape front while structured triangular grids produce a hexagonal shape (Ortigoza 2015). This means that Cellular Automata with simple rules are anisotropic at local scale, i.e. the signals propagate with different celerities in different directions, and this anisotropy reverberates at global scale (Schönfisch 1997) producing spurious geometric regularities. These effects can be quantified, and Weimar et al. (1992) found that the ratio between wave celerities in diagonal and horizontal directions is in Cellular Automata algorithms with square lattice and simple evolution rules. The continuous-state Cellular Automata paradigm (discrete variables in space and a rule for their evolution) is generic enough to embrace also classic numerical schemes for the solution of partial differential equations such as Finite Difference (Yang and Young 2006) and Finite Volume schemes. Not surprisingly, Caviedes-Voullième et al. (2018) have correctly stated that many Cellular Automata algorithms for flooding delineation correspond to Finite Volume discretizations of an appropriate system of differential equations. Nonetheless, these Authors have failed in individuating the correct correspondence between algorithms and differential equations. For example, the discrete model by Bates and De Roo (2000), which can be viewed as a Cellular Automata algorithm where the fluxes between cells are calculated by means of the one-dimensional Zero-Inertia model, is not a Finite Volume discretization of Eq. (25) as stated by Caviedes-Voullième et al. (2018), but it is a discretization of Eq. (30). Of course, there is a formal correspondence between flow differential models with decoupled friction and Cellular Automata with simplified exchange rules between cells. Not surprisingly, Figures 3 and 4 (right columns) show that decoupled-friction models produce ropagating waves whose isolines tend to be parallel to the grid axes. Strikingly, the ratio between wave celerities in diagonal and horizontal directions found by Weimar et al. (1992) coincides with the ratio ' f f S S for the bisectors of the quadrants in the flow propagation models with decoupled friction considered here.
In the present paper, it has been demonstrated that the decoupling of the friction vector introduces the loss of the mathematical model rotational invariance and the distortion of numerical solutions. In addition to these findings, we want to show how numerical choices may introduce unwanted decoupling effects even in the case that the original mathematical model is based on coupled friction. In Section 3, the discretization of the LInA numerical model by de Almeida and Bates (2013) has been constantly applied with relaxation parameter λ = 1, which ensures isotropy of the numerical model in the case of isotropic treatment of the friction. Nonetheless, fixed values of λ in the range [0.7, 0.9] (de Almeida et al. 2012, de Almeida and Bates 2013) or dynamically variable values with λ
1 (Sridharan et al. 2020) are recommended in the literature in order to ensure stability of computations in the case of flow propagation on dry bed with small friction. For this reason, the coupled-friction LInA model of Sub-section 3.1.1 with λ = 1 is contrasted in Figure 7 (left column) with the same model where the parameter λ = 0.7 is applied (right column). The inspection of the figure shows that the coupled-friction by de Almeida and Bates (2013) with λ < 1 causes the distortion of the propagating front, which is progressively converted into a straight line inclined by - /4 (right column of Figure 7). This effect is not present in the same model with λ = 1 (left column of Figure 7). It is evident that the expansion of the q x stencil along x nd the q y stencil along y for λ < 1 causes the partial decoupling of q x and q y , but this point needs further studies in the future. Also in this case, the propagating wave distortion cannot be individuated by one-dimensional case studies or complicated real-world problems, while it is well evidenced by the use of a truly two-dimensional case study in controlled conditions. [Insert Figure 7 about here]
5. Conclusions
Friction decoupling, i.e. the computation of friction vector components making separate use of the corresponding velocity components, is common in staggered grid models of the SWE simplifications (Zero-Inertia and Local Inertia Approximation), due to the programming simplicity and to the corresponding calculations speed-up. Despite the interest in understanding the modifications introduced by friction decoupling, scarce attention has been given in the literature to this topic. In the present paper, the effect of friction decoupling in flooding models has been studied from the theoretical and numerical point of view. First, it has been found that the decoupling of the friction vector causes the reduction of the computed friction force and rotation of the friction line of action. Second, it has been demonstrated that decoupled-friction models lack of rotational invariance, i.e. the results of the model depend on the alignment of the reference framework. These theoretical results have been confirmed by means of numerical experiments. In particular, it has been found that the friction decoupling causes the distortion of the propagating wave in a radially symmetric experiment, and that the distortion depends on the computational grid alignment. This issue is in agreement with the theory of Cellular Automata where oversimplified rules are used for the exchanges between cells. In addition to these results, it has been found that the celerity of the ropagating wave is overestimated in decoupled-friction models. Notably, calibration of the friction coefficient does not solve these issues. On the basis of the above findings, it is evident that the decoupling of friction vector causes a major loss of credibility of the corresponding mathematical and numerical models, because it is impossible to a priori predict how much the decoupled-friction solution will depart from the coupled-friction one. Despite the modest speed-up of decoupled-friction computations, the use of coupled-friction models is preferable in every case.
Acknowledgments
The research project “Sistema di supporto decisionale per il progetto di casse di espansione in linea in piccoli bacini costieri” was funded by the Italian Ministry for Environment, Land and Sea Protection through the funding program “Metodologie per la valutazione dell’efficacia sulla laminazione delle piene in piccoli bacini costieri di sistemi di casse d’espansione in linea realizzate con briglie con bocca tarata”.
Appendices
A. Rotation of the frame reference for the Shallow water Equations
In this Appendix it is demonstrated that the Eq. (14) is obtained from Eq. (10) after a rotation of the reference framework. The application of the first and the second of Eq. (13) to the first of Eq. (1) leads to A.1) cos sin sin cos 0 y yx x q qq qht X Y X Y . Recalling that cos sin 1 , the substitution of the first and the second of Eq. (12) into Eq. (A.1) leads to the first of Eq. (14). The application of the same procedure to the second of Eq. (1) leads to (A.2)
222 7 3222 7 3 cos 2sin 02
M XX X X Y M YY X Y Y n QQ Q Q Qgh zgh gt X h Y h X hn QQ Q Q Qgh zgh gt X h Y h Y h QQ . Similarly, the application of the rotation procedure to the third of Eq. (1) leads to (A.3)
222 7 3222 7 3 sin 2cos 02
M XX X X Y M YY X Y Y n QQ Q Q Qgh zgh gt X h Y h X hn QQ Q Q Qgh zgh gt X h Y h Y h QQ . If Eqs. (A.2) and (A.3) are summed after the multiplication by cos and sin , respectively, the second of Eq. (A.14) is obtained. Similarly, if Eqs. (A.2) and (A.3) are summed after the multiplication by -sin and cos , respectively, the third of Eq. (A.14) is obtained. The proof is complete. References rnold V.I. (1989) Mathematical methods of classical mechanics, 2 nd Ed., Springer-Verlag, New York. Artichowicz W., Gaşiorowski D. (2019) Computationally efficient solution of a 2D diffusive wave equation used for flood inundation problems, Water 11(10), 2195. Doi: 10.3390/w11102195. Bates P.D., De Roo A.P.J. (2000) A simple raster-based model for flood inundation simulation, Journal of Hydrology 236(1-2), 54-77. doi: 10.1016/S0022-1694(00)00278-X. Bates P.D., Horritt M.S., Fewtrell T.J. (2010) A simple inertial formulation of the shallow water equations for efficient two-dimensional flood inundation modelling, Journal of Hydrology 387(1-2), 33-45. doi: j.jhydrol.2010.03.027. Bouchut F. (2004) Nonlinear Stability of Finite Volume Methods for Hyperbolic Conservation Laws and Well-Balanced Schemes for Sources, Birkhäuser, Basel. Bradbrook K.F., Lane S.N., Waller S.G., Bates P.D. (2004) Two dimensional diffusion wave modeling of flood inundation using a simplified channel representation, International Journal of River Basin Management 2(3), 211-223. doi: 10.1080/15715124.2004.9635233. Cai X., Li Y., Guo X.-w., Wu W. (2014) Mathematical model for flood routing based on cellular automaton, Water Science and Engineering 7(2), 133-142. doi: 10.3882/j.issn.1674-2370.2014.02.002. Caviedes-Voullième D., Fernández-Pato J., Hinz C. (2018) Cellular Automata and Finite Volume solvers converge for 2D shallow flows modelling for hydrological modelling, Journal of Hydrology 563, 411-417. doi: j.jhydrol.2018.06.021. Clark M.P., Schaefli B., Schymanski S.J., Samaniego L., Luce C.H., Jackson B.M., Free J.E., Arnold J.R., Moore R.D., Istanbulloglu E., Ceola S. (2016) Improving the theoretical underpinnings of process ‐ based hydrologic models, Water Resources Research 52(3), 2350-2365. Doi: 10.1002/2015WR017910. Coulthard T.J., Neal J.C., Bates P.D., Ramirez J., de Almeida G.A.M., Hancock G.R. (2013) Integrating the LISFLOOD-FP 2D hydrodynamic model with CAESAR model: implications for odeling landscape evolution, Earth Surface processes and Landforms 38(15), 1897-1906. doi: 10.1002/esp.3478. Cozzolino L., Cimorelli L., Della Morte R., Pugliano G., Piscopo V., Pianese D. (2019) Flood propagation modeling with the Local Inertia Approximation: Theoretical and numerical analysis of its physical limitations, Advances in Water Resources 133, 103422. doi: 10.1016/j.advwatres.2019.103422. de Almeida G.A.M., Pates P.D., Freer J.E., Souvignet M. (2012) Improving the stability of a simple formulation of the shallow water equations for 2-D flood modeling, Water Resources Research 48(5), W05528. doi: 10.1029/2011WR011570. de Almeida G.A.M., Bates P. (2013) Applicability of the local inertial approximation of the shallow water equations to flood modeling, Water Resources Research 49(8), 4833-4844. doi: 10.1002/wrcr.20366. Delis A.I., Mathioudakis E.N. (2009) A finite volume method parallelization for the simulation of free surface shallow water flows, Mathematics and Computers in Simulation 79(11), 3339-3359. Doi: 10.1016/j.matcom.2009.05.010. Dottori F., Todini E. (2011) Developments of a flood inundation model based on the cellular automata approach: testing different methods to improve model performance, Physics and Chemistry of the Earth Parts A/B/C 36(7-8), 266-280. doi: 10.1016/j.pce.2011.02.004. Falter D., Vorogushyn S., Lhomme J., Apel H., Gouldby B., Merz B. (2013) Hydraulic model evaluation for large-scale flood risk assessments, Hydrological Processes 27(9), 1331-1340. doi: 10.1002/hyp.9553. Fewtrell T.J., Bates P.D., Horritt M., Hunter N.M. (2008) Evaluating the effect of scale in flood inundation modelling in urban environments, Hydrological processes 22(26), 5107-5118. doi: 10.1002/hyp/.7148. albiati G., Savi F. (1995) Evaluation of the comparative influence of soil hydraulic propoerties and roughness on overland flow at the local scale, Journal of Agricultural Engineering Research 61(3), 183-190. doi: 10.1006/jaer.1995.1045. Ghimire B., Chen A.S., Guidolin M., Keedwell E.C., Djordjević S., Savić D. (2013) Formulation of a fast 2D urban pluvial flood model using a cellular automata approach, Journal of Hydroinformatics 15(3), 676-686. doi: 10.2166/hydro.2012.245. Gioia G., Bombardelli F.A. (2002) Scaling and similarity in rough channel flows, Physical Review Letters 88(1), 014501. Doi: 10.1103/PhysRevLett.88.014501. Gualtieri P., De Felice S., Pasquino V., Pulci Doria G. (2018) Use of conventional flow resistance equations and a model for the Nikuradse roughness in vegetated flows at high submergence, Journal of Hydrology and Hydromechanics 66(1), 107-120. Doi: 10.1515/johh-2017-0028. Horritt M.S., Bates P.D. (2001) Predicting floodplain inundation: raster-based modelling versus the finite-element approach, Hydrological Processes 15(5), 825-842. doi: 10.1002/hyp.188. Horritt M.S., Bates P.D. (2002) Evaluation of 1D and 2D numerical models for predicting river flood inundation, Journal of Hydrology 268(1-4), 87-99. doi: 10.1016/S0022-1694(02)00121-X. Hromadka II T.V., Chintu Lai (1985) Solving the two-dimensional diffusion flow model, In: William R. Waldrop Ed., Hydraulics and Hydrology in the small computer age, ASCE, New York. Hunter N.M., Horritt M.S., Bates P.D., Wilson M.D., Werner M.G.F. (2005) An adaptive time step solution for raster-based storage cell modelling of floodplain inundation, Advances in Water Resources 28(9), 975-991. doi: 10.1016/j.advwatres.2005.03.007. Hunter N.M., Bates P.D., Horritt M.S., Wilson M.D. (2007) Simple spatially-distributed models for predicting flood inundation: a review, Geomorphology 90(3-4), 208-225. doi: j.geomorph.2006.10.021. eandro J., Chen A.S., Schumann A. (2014) A 2D parallel diffusive wave model for floodplain inundation with variable time step (P-DWave), Journal of Hydrology 517, 250-250. doi: j.jhydrol.2014.05.020. Liu L., Liu Y., Wang X., Yu D., Liu K., Huang H., Hu G. (2015) Developing an effective 2-D urban flood inundation model for city emergency management based on cellular automata, Natural Hazards and Earth System and Sciences 15, 381-391. doi: 10.5194/nhess-15-381-2015 Mateo C.M.R., Yamazaki D., Kim H., Champathong A., Vaze J., Oki T. (2017) Impacts of spatial resolution and representation of flow connectivity on large-scae simulation of floods, Hydrology and Earth System Sciences 21(10), 5143-5163. doi: 10.5194/hess-21-5143-2017. Mendicino G., Pedace J., Senatore A. (2015) Stability of an overland flow scheme in the frame work of a fully coupled eco-hydrological model based on the Macroscopic Cellular Automata approach, Communications in Nonlinear Science and Numerical Simulation 21(1-3), 128-146. doi: j.cnsns.2014.08.021. Neal J., Schumann G., Fewtrell T., Budimir M., Bates P., Mason D. (2011) Evaluating a new LISFLOOD-FP formulation with data from summer 2007 floods in Tewkesbury, UK, Journal of Flood Risk Management 4(2), 88-95. doi: 10.1111/j.1753-318X.2011.01093.x. Neal J., Villanueva I., Wright N., Willis T., Fewtrell T., Bates P. (2012) How much physical complexity is needed to model flood inundation?, Hydrological Processes 26(15), 2264-2282. doi: 10.1002/hyp.8339. Neal J., Dunne T., Sampson C., Smith A., Bates P. (2018) Optimisation of the two-dimensional hydraulic model LISFOOD-FP for CPU architecture, Environmental Modelling & Software 107, 148-157. Doi: 10.1016/j.envsoft.2018.05.011. Nguyen B.Y., Ichikawa Y., Ishidaira H. (2016) Estimation of inundation depth using flood extent information and hydrodynamic simulations, Hydrological Research Letters 10(1), 39-44. doi: 10.3178/hrl.10.39. berkampf W.L., DeLand S.M., Rutherford B.M., Diegert K.V., Alvin K.F. (2002) Error and uncertainty in modelling and simulation, Reliability Engineering and System Safety 75(3), 333-357. doi: 0.1016/S0951-8320(01)00120-X. Ortigoza G.M. (2015) Unstructured triangular cellular automata for modelling geographic spread, Applied Mathematics and Computation 258, 520-536. doi: j.amc.2015.01.116. Panday S., Huyakorn P.S. (2004) A fully coupled-physically-based spatially-distributed model for evaluating surface/subsurface flow, Advances in Water Resources 27(4), 361-382. doi: j.advwatres.2004.02.016. Rinaldi P.R., Dalponte D.D., Vénere M.J., Clausse A. (2007) Cellular automata algorithm for simulation of surface flows in large plains, Simulation Modelling Practice and Theory 15(3), 315-327. doi: j.simpat.2006.11.003. Savage J.T.S., Pianosi F., Bates P., Freer J., Wagener T. (2016) Quantifying the importance of spatial resolution and other factors through global sensitivity analysis of a flood inundation model, Water Resources Research 52(11), 9146.9163. doi: 10.1002/2015WR018198. Schönfisch B. (1995) Propagation of fronts in cellular automata, Physica D 80(4), 433-450. doi: 10.1016/0167-2789(94)00192-S. Schönfisch B. (1997) Anisotropy in cellular automata, BioSystems 41(1), 29-41. doi: 10.1016/S0303-2647(96)01664-4. Sridharan B., Gurivindapalli D., Kuiry S.N., Mali V.K., Devi N.N., Bates P.D., Sen D. (2020) Explicit expression of weighting factor for improved estimation of numerical flux in Local Inertial models, Water Resources Research 56(7), e2020WR027357. Doi: 10.1029/2020WR027357. Stelling G.S., Duinmeijer S.P.A. (2003) A staggered conservative scheme for every Froude number in rapidly varied shallow water flows, International Journal for Numerical Method in Fluids 43(12), 1329-1354. Doi: 10.1002/fld.537. Thomas R., Nicholas A.P. (2002) Simulation of braided river flow using a new cellular routing scheme, Geomorphology 43(3-4), 179-195. doi: 10.1016/S0169-555X(01)00128-3. oro E.F. (2001) Shock-capturing method for free-surface shallow flows, John Wiley & Sons, Chichester. Varra G., Pepe V., Cimorelli L., Della Morte R., Cozzolino L. (2020) On integral and differential porosity models for urban flooding simulation, Advances in Water Resources 136, 103455. Doi: 10.1016/j.advwatres.2019.103455. Viero D.P., Valipour M. (2017) Modeling anisotropy in free-surface overland and shallow inundation flows, Advances in Water Resources 104, 1-14. Doi: 10.1016/j.advwatres.2017.03.007. Wasantha Lal A.M. (1998) Performance comparison of overland flow algorithms, ASCE Journal of Hydraulic Engineering 124(4), 342-349. doi: 10.1061/(ASCE)0733-9429(1998)124:4(342). Weimar J., Tyson J.J., Watson L.T. (1992) Diffusion and wave propagation in cellular automaton models of excitable media, Physica D 55(3-4), 309-327. doi: 10.1016/0167-2789(92)90062-R Xanthopoulos T., Koutitas C. (1976) Numerical simulation of a two dimensional flood wave propagation due to dam failure, IAHR Journal of Hydraulic Research 14(4), 321-331. doi: 10.1080/00221687609499664. Yamazaki D., de Almeida G.A.M., Bates P.D. (2013) Improving computational efficiency in global river models by implementing the local inertial flow equation and a vector-based river network map, Water Resources Research 49(11), 7221-7235. doi: 10.1002/wrcr.20552. Yang X.-S., Young Y. (2006) Cellular Automata, PDEs, and Pattern Formation, In: Stephan Olariu and Albert Y. Zomaya Eds., Handbook of Bioinspired Algorithms and Applications, Chapman & Hall/CRC, Boca Raton. Yu D., Lane S.N. (2006) Urban fluvial modelling using a two-dimensional diffusion-wave treatment, part 1: mesh resolution effects, Hydrological processes 20(7), 1541-1565. doi: 10.1002/hyp.5935. Figures
Figure 1. Comparison of coupled and decoupled friction vectors: ratio between the moduli (a); absolute deviation angle (b). Figure 2. Staggered grid stencil. Figure 3. Propagation of a circular wave. Free-surface elevation with the Zero-Inertia model at time t = 60000 s. Left column, coupled friction: s = 3218.6 m (a); s = 1609.3 m (c); s = 804.65 m (e). Right column, decoupled friction: s = 3218.6 m (b); s = 1609.3 m (d); s = 804.65 m (f). Figure 4. Propagation of a circular wave. Free-surface elevation with the LInA model ( = 1) at time t = 60000 s. Left column, coupled friction: s = 3218.6 m (a); s = 1609.3 m (c); s = 804.65 m (e). Right column, decoupled friction: s = 3218.6 m (b); s = 1609.3 m (d); s = 804.65 m (f). Figure 5. Propagation of a circular wave. Free-surface elevation with the LInA model ( = 1) at time t = 60000 s. Thin line: coupled friction with n M = 0.03 s m -1/3 . Thick dashed line: decoupled friction with n M = 0.0325 s m -1/3 . Figure 6. Propagation of a circular wave. Free-surface elevation with the LInA model ( = 1) and decoupled friction at time t = 60000 s (refined grid): reference framework Oxy (a); reference framework
OXY rotated by = /4 with respect to the reference Oxy . Dashed line: stretching direction of the propagating wave. Figure 7. Propagation of a circular wave. Free-surface elevation with the coupled-friction LInA model at time t = 60000 s. Left column, = 1: s = 3218.6 m (a); s = 1609.3 m (c); s = 804.65 m (e). Right column, = 0.7: s = 3218.6 m (b); s = 1609.3 m (d); s = 804.65 m (f). Tables
Table 1. Resume of constant time steps and uniform grid sides used in the computations. Table 2. Run times for the refined grid. able 1. Resume of constant time steps and uniform grid sides used in the computations.
Grid Panel s (m) t (s) Coarse Upper: Fig. 3a,b; Fig. 4a,b; Fig. 7a,b 3218.6 50 Medium Central: Fig. 3c,d; Fig. 4c,d; Fig. 7c,d 1609.3 12.5 Refined Lower: Fig. 3e,f; Fig. 4e,f; Fig. 7e,f 805.65 3.125 able 2. Run times for the medium grid
Model Run time (ms)
Speedup
Zero-Inertia with coupled-friction 18626 1.25 Zero-Inertia with decoupled-friction 14867 LInA with coupled-friction 33493 1.17 LInA with decoupled-friction 28673 igure 1. Comparison of coupled and decoupled friction vectors: ratio between the moduli (a); absolute deviation angle (b). igure 2. Staggered grid stencil. igure 3. Propagation of a circular wave. Free-surface elevation with the Zero-Inertia model at time t = 60000 s. Left column, coupled friction: s = 3218.6 m (a); s = 1609.3 m (c); s = 804.65 m (e). Right column, decoupled friction: s = 3218.6 m (b); s = 1609.3 m (d); s = 804.65 m (f). igure 4. Propagation of a circular wave. Free-surface elevation with the LInA model ( = 1) at time t = 60000 s. Left column, coupled friction: s = 3218.6 m (a); s = 1609.3 m (c); s = 804.65 m (e). Right column, decoupled friction: s = 3218.6 m (b); s = 1609.3 m (d); s = 804.65 m (f). igure 5. Propagation of a circular wave. Free-surface elevation with the LInA model ( = 1) at time t = 60000 s. Thin line: coupled friction with n M = 0.03 s m -1/3 . Thick dashed line: decoupled friction with n M = 0.0325 s m -1/3 . igure 6. Propagation of a circular wave. Free-surface elevation with the LInA model ( = 1) and decoupled friction at time t = 60000 s (refined grid): reference framework Oxy (a); reference framework
OXY rotated by = /4 with respect to the reference Oxy . Dashed line: stretching direction of the propagating wave. igure 7. Propagation of a circular wave. Free-surface elevation with the coupled-friction LInA model at time t = 60000 s. Left column, = 1: s = 3218.6 m (a); s = 1609.3 m (c); s = 804.65 m (e). Right column, = 0.7: s = 3218.6 m (b); s = 1609.3 m (d);