10
pages

A Shallow Water model for the numerical simulation of overland flow on surfaces with ridges and furrows Ulrich Razafison?,a, Stephane Cordierb, Olivier Delestrec, Frederic Darbouxd, Carine Lucasb, Franc¸ois Jamesb aUniversite de Franche-Comte, Laboratoire de Mathematiques, CNRS UMR 6623, 16 route de Gray, 25030 Besanc¸on Cedex, France bUniversite d'Orleans, Laboratoire MAPMO, CNRS UMR 6628, Federation Denis Poisson, B. P. 6759, 45067 Orleans Cedex 2, France cUniversite de Nice - Sophia Antipolis, Laboratoire de Mathematiques J. A. Dieudonne, UMR 6621 CNRS UNSA, Parc Valrose, 06108 Nice Cedex 2, France dINRA, UR 0272 Science du sol, Centre de recherche d'Orleans, CS 40001, F-45075 Orleans Cedex 2, France Abstract We introduce a new Shallow Water model for the numerical simulation of overland flow with furrow effects without representing them explicitly. The model is obtained by adding an anisotropic friction term that takes into account these effects to the classical Shallow Water equations. We validate the model with numerical tests, and we compare it with the classical Shallow Water model where the furrows are explicitly and precisely described. Key words: Overland flow, Shallow Water equations, Furrows, Friction 1. Introduction During rainfall, overland flow on cultivated lands induces problems at the watershed scale for soil conservation (decreases soil thickness by erosion and causes nutrient loss), infrastruc- tures (flooding and destruction of roads and buildings), preser- vation of water quality (drinking water) and sustainability of aquatic ecosystems (chemical pollution).

- sinus furrows
- furrows
- friction term
- velocity over
- muscl reconstruction
- friction law
- horizontal flow velocity
- overland flow
- better predict

Voir plus
Voir moins

Vous aimerez aussi

on surfaces with ridges and furrows

∗,a b c d b bUlrich Razaﬁson , Ste´phane Cordier , Olivier Delestre , Fre´de´ric Darboux , Carine Lucas , Franc¸ois James

aUniversite´ de Franche-Comte´, Laboratoire de Mathe´matiques, CNRS UMR 6623, 16 route de Gray, 25030 Besanc¸on Cedex, France

bUniversite´ d’Orle´ans, Laboratoire MAPMO, CNRS UMR 6628, Fe´de´ration Denis Poisson, B. P. 6759, 45067 Orle´ans Cedex 2, France

cUniversite´ de Nice - Sophia Antipolis, Laboratoire de Mathe´matiques J. A. Dieudonne´, UMR 6621 CNRS UNSA, Parc Valrose, 06108 Nice Cedex 2, France

dINRA, UR 0272 Science du sol, Centre de recherche d’Orle´ans, CS 40001, F-45075 Orle´ans Cedex 2, France

Abstract

We introduce a new Shallow Water model for the numerical simulation of overland ﬂow with furrow eﬀects without representing

them explicitly. The model is obtained by adding an anisotropic friction term that takes into account these eﬀects to the classical

Shallow Water equations.

We validate the model with numerical tests, and we compare it with the classical Shallow Water model where the furrows are

explicitly and precisely described.

Key words: Overland ﬂow, Shallow Water equations, Furrows, Friction

1. Introduction lands. This interaction can be seen as the interaction between

three types of roughness. The topography is the roughness of

During rainfall, overland ﬂow on cultivated lands induces the Earth and is described on Digital Elevation Maps with a

problems at the watershed scale for soil conservation (decreases horizontal resolution larger than one meter and commonly of

soil thickness by erosion and causes nutrient loss), infrastruc- ten meters and more. Furrows are the roughness due to agri-

tures (ﬂooding and destruction of roads and buildings), preser- cultural practices and create a strong directional heterogeneity

vation of water quality (drinking water) and sustainability of inside a ﬁeld. They are characterized by their wavelength (one

aquatic ecosystems (chemical pollution). to a few decimeters), their amplitude (a few centimeters to one

These troubles can be prevented by improving watershed decimeter) and their direction. Finally, the random roughness

management in connection with overland ﬂow. Thus, the water due to soil aggregates and clods is homogeneous in space and

ﬂux at the outlet not only must be simulated well but also must has an amplitude of a few millimeters to about one decimeter.

predict the spatial distribution of the water ﬂux and velocity To our knowledge, most of the research on the interaction be-

over the whole watershed well. However, current hydrological tween roughness and ﬂow have been dedicated to topography

models predict overland ﬂow within small watersheds ineﬃ- [5, 6] or to random roughness [7, 8, 9].

ciently [1, 2, 3]. In agricultural watersheds, one of the main Few studies have addressed furrows, and most are concerned

diﬃculties is that ﬂow directions are controlled not only by the with the storage capacity of the furrows, i.e., the amount of

topography but also by ditches along the ﬁeld boundaries and water stored in the puddles created by the furrows (e.g., [10]).

by ridges and furrows created by tillage operations inside the These studies do not consider the water ﬂowing on the soil sur-

ﬁelds. The ﬂow pattern is clearly the result of the interaction faces but rather the water stored in puddles. The few studies

between these objects [4], but the way they interact remains considering both overland ﬂow and the furrows-topography in-

mostly unspeciﬁed. Therefore, this interaction must be better teraction are empirical studies [4, 11]. They lead to empirical

understood to better predict the spatial and temporal distribu- laws giving an on/oﬀ prediction: the predicted ﬂow direction is

tions of overland ﬂow and to improve the decisions made by either the direction of the topographic slope or the furrow direc-

watershed managers. tion, while water can ﬂow in both directions at the same time

In this paper, we focus on the interaction between topogra- in reality. Moreover, these laws are limited by their empirical

phy and furrows, a feature encountered in almost all cultivated basis.

To be of practical use, a model accounting for the eﬀects of

∗ furrows on overland ﬂow direction must not require an explicitCorresponding author. Tel: +33 381666397; fax: +33 381666623

Email addresses: ulrich.razafison@math.cnrs.fr (Ulrich representation of the furrows: such a requirement would require

Razaﬁson), stephane.cordier@univ-orleans.fr (Ste´phane Cordier), the use of a digital topographic map with a horizontal resolution

olivierdelestre41@yahoo.fr (Olivier Delestre),

of about one centimeter for the whole watershed, and a small

frederic.darboux@orleans.inra.fr (Fre´de´ric Darboux),

watershed covers approximately one square kilometer. Suchcarine.lucas@univ-orleans.fr (Carine Lucas),

francois.james@univ-orleans.fr (Franc¸ois James) digital maps are not available and, even if available, require too

Preprint submitted to European Journal of Mechanics-B/Fluids December 2, 2011free surface

many computational resources.

The purpose of this study is to propose a model that can ac-

count for the eﬀects of the furrows on overland ﬂow. Numerical

results are presented. The model is a ﬁrst step in an attempt to h(t,x)

u(t,x)

predict overland ﬂow directions controlled by furrows and to-

pography without representing the furrows explicitly. Indeed,

only average amplitude, wavelength and direction are used to

characterize the furrows. In this paper, the furrow direction is

Z(x)

kept perpendicular to the slope. Our model is based on the Shal-

low Water equations that are widely used to describe ﬂows in

rivers and ocean and overland, among other applications. Figure 1: Notations for a 1D Shallow Water ﬂow.

The outline of the paper is as follows. In the next section,

we ﬁrst present the Shallow Water model. Then, we propose a

topographynew model in which we add a new friction term to account for

the eﬀects of the furrows on overland ﬂow. Section 3 describes z

the numerical scheme used to solve the model, and in Section 4,

0.02we present and discuss the numerical results that we obtain with 0

-0.02

our model. Conclusions are outlined in Section 5. -0.04

-0.06

-0.08

-0.1

-0.122. Mathematical models -0.14

-0.16

-0.18

The starting point is the 2D classic Shallow Water sys-

0.2

tem [12] in a bounded domainΩ: 0.15

0 0.1

0.5 ∂h ∂(hu) ∂(hv) x ∈ [0;ℓ] 1 1.5 0.05 + + = R, 2 2.5 ∂t ∂x ∂y y ∈ [0; L] 3 0 3.5 2 ∂(hu) ∂(hu ) ∂(huv) ∂h + + + gh ∂t ∂x ∂y ∂x Figure 2: An example of topography with furrows. ∂Z 2 −1/3 (2.1) +gh + gk h |u|u= 0, ∂x 2 ∂(hv) ∂(huv) ∂(hv ) ∂h 1. The direction of the ﬂow is parallel to the length of the + + + gh ∂t ∂x ∂y ∂y domainΩ with respect to y (pseudo-1D case) and, conse- ∂Z 2 −1/3 quently, perpendicular to the furrows. +gh + gk h |u|v= 0. ∂y 2. We only consider ﬂuvial ﬂows, which means that |u| <p

gh.For t > 0 and x = (x, y) ∈ Ω, the unknowns are the water

height h= h(t, x) and the horizontal ﬂow velocity u= u(t, x)= 3. Inﬁltration and soil erosion are not taken into account.

T(u(t, x), v(t, x)) . Furthermore, Z(x) describes the bottom to-

Under such assumptions, the furrows overﬂow at the same time

pography of the domain; therefore, h+ Z is the level of the wa-

during rainfall events or one after the other during inﬂow from

ter surface (Figure 1). In equations (2.1), g is the acceleration

upstream.

due to the gravity and R is the rainfall intensity. Several studies

We propose a model that takes into account the eﬀects of the

have shown a derivation of the Shallow Water system originat-

furrows without explicitly representing them in the topography

ing from the free surface Navier-Stokes equations [13, 14, 15].

Z. In other words, we want to ﬁnd an equivalent model to the

For the friction term, we choose the Manning law with

Shallow Water system onΩ that can be used at a macroscopic

k as the Manning coeﬃcient. We also denote q(t, x) =

scale, i.e., on a topography that is only an inclined plane. WeT(q (t, x), q (t, x)) = h(t, x)u(t, x) as the water ﬂux.x y want to force the ﬂow to slow down when its depth is smaller

Now, we consider a rectangular domainΩ= ℓ × L and a to-

than the value corresponding to the water height that can be

pography Z with furrows. We suppose that the topography is an

trapped in the furrows. The idea of this article is to model this

inclined plane with sinus furrows and that the geometry of the

eﬀect of the furrows through an additional friction term that

furrows is known through their amplitude and their wavelength.

forces the ﬂow to slow down for a low water depth. To that end,

Note that realistic (measured) furrows are shaped slightly dif-

we ﬁrst introducehh i the average height of the water trappedF

ferently due to the existence of random roughness. However,

in the furrows. This value is given by

random roughness, because it is isotropic, does not aﬀect ﬂow

direction at the scale of the furrows. We also suppose that the hh i= V/(L ×ℓ) (m), (2.2)F F

furrows are perpendicular to the length ofΩ with respect to y.

An example of such topography is illustrated in Figure 2. where V is the volume of trapped water in a furrow, L is itsF

Next, we complement the problem with the following as- wavelength (see Figure 3) and ℓ is the length of the domainΩ

sumptions. (with respect to x). Note that the value ofhh i only depends onF

2

Trapped water (volume V) Finally, we propose the following new Shallow Water modelz

with a “furrows-friction” coeﬃcient:

h(t,y) ∂h ∂(hu) ∂(hv)Z(y) + + = R, ∂t ∂x ∂y 2 ∂(hu) ∂(hu ) ∂(huv) ∂h + + + gh ∂t ∂x ∂y ∂xy yy +L0 0 F ∂Z 2 −1/3 (2.4) +gh + gk h |u|u= 0, ∂xFigure 3: Water trapped in a furrow. 2 ∂(hv) ∂(huv) ∂(hv ) ∂h + + + gh ∂t ∂x ∂y ∂y ∂Z 2 −1/3 +gh + gk h |u|v+ K(h)hv= 0.

0.25 ∂y

0.20

Remark 2.2. 1. Note that, because the furrows are perpen-

dicular to the slope, the additional friction law K(h)hv only

0.15 appears in the third equation of (2.4), and therefore, it only

acts on the ﬂow in the y-axis direction. This assumption is

not restrictive: in general, the direction of the furrows is0.10

constant on each agricultural ﬁeld, and, if necessary, we

apply a rotation to obtain the equations for an arbitrary di-

0.05 rection of the furrows.

2. The form of the new friction law is chosen arbitrarily. The

α β

0.00 general form of friction laws is Kh |u| u whereα andβ are

0.000 0.005 0.010 0.015 0.020 0.025 0.030

positive real numbers. For example, we obtain the Man-h

ning law for (α,β) = (−1/3, 1) and Darcy-Weisbach law

for (α,β) = (1, 1). Note that these laws are empirical and

Figure 4: Shape of the friction term K(h) forhh i= 0.01 m.F are obtained considering stationary ﬂows, and their valid-

ity is still discussed among hydrologists (e.g., [16]).

For the numerical experiments presented in Section 4 withthree parameters: the slope of the domain, the furrows average

system (2.4), we chose (α,β) = (1, 0) but we could haveamplitude and the furrows average wavelength.

made another choice, should we change the value of K or0Next, we consider the following additional friction coeﬃ-

the form of K(h).cient: !

−h+hh iF

K(h)= K exp , (2.3)0 Because shallow water ﬂows can also be described by the

Chh iF

so-called multi-layer Shallow Water system (e.g., [17, 18, 19]

where C is a characteristic constant increasing function of the for a derivation and numerical studies), we can propose another

small random variations of the height of the furrows, and K is0 approach based on multi-layer models to take into account the

a coeﬃcient we determine in the following. eﬀects of the furrows on overland ﬂows.

In Figure 4, the general shape of K(h) is plotted for hh i =F In this work, we introduce the following two-layer like

0.01 m. We clearly see that K(h) is large for h ≤ hh i. Thus,F model:

when the water height h is lower than the average height of the

if h(t, x) ≤ hh i, then u(t, x)= 0 and h(t, x)= Rt,furrowshh i, the ﬂow slows as a result of K(h),. FF

(2.5) if h(t, x) > hh i, FRemark 2.1. 1. In (2.3), if hh i tends to 0, then K(h) alsoF then solve (2.1) with an inclined plane topography.

tends to 0 for any h > 0. In other words, the additional

friction coeﬃcient disappears when there are no furrows. In (2.5), the lower layer corresponds to the ﬁlling up of the fur-

2. If C tends to 0, then we obtain the empirical models rows; note that the upper layer is active only when the furrows

that are usually used. These models consist in giving overﬂow. The initial conditions for the upper layer are then

an on/oﬀ prediction of the furrows-topography interaction ˆ ˆu(0, x)= 0 and h(0, x)= h−hh i, where h is the water height atF

(see [4, 11]); more precisely, while the critical water height the overﬂow time. Note that in one dimension, this model pro-

is not attained, there is no ﬂow, and after this threshold, the vides more satisfactory results than the model (2.4) (Section 4).

furrows are not taken into account. However, its extension to more complex two-dimensional prob-

The new Shallow Water model we introduce here (2.4) can lems requires modeling the coupling between the two layers

be seen as an improvement of these models. carefully, and it is more diﬃcult than extending the model (2.4).

3

K !

3. Numerical scheme

where c < c are given by c = inf inf λ (U) , c =1 2 1 j 2

U=U ,U j=1,2l r

In this section, we explain the numerical scheme we used in p sup sup λ (U) and where λ (U) = u − gh, λ (U) = j 1 2our numerical simulations. The Shallow Water system is not

U=U ,U j=1,2l rpas easy to solve. In hydrology, the Mac Cormack scheme is ′u+ gh are the eigenvalues of the Jabobian matrix F (U). A

usually used for overland ﬂow simulation (e.g., [20, 21]). How-

study compared diﬀerent numerical ﬂuxes for (3.6) and showed

ever, it is not well adapted to this system because of several

that, in the framework of overland ﬂow, the HLL ﬂux is the

problems, such as the preservation of the positivity of the water

best compromise between the accuracy of the approximation

height and of the steady states or the behavior at the wet/dry in-

and the computational cost [27]. To have a second order ac-

terface. To that end, we used well-balanced schemes [22] based

curacy scheme, we use the modiﬁed ENO reconstruction [28]

on the hydrostatic reconstruction [23, 24]. This ﬁnite volume

deﬁned as follows

scheme has shown to be adapted to overland ﬂow simulation at

small scales [25, 26, 27]. Δx Δx

h = h − D h , h = h + D h ,i−1/2+ i enom i i+1/2− i enom iTo make this presentation simpler, we describe the numeri- 2 2

cal scheme on the classical the one-dimensional Shallow Water h Δxi+1/2−

u = u − D u ,i−1/2+ i enom imodel with variable topography and Manning’s friction law : h 2i

h Δxi−1/2+ u = u + D u∂h ∂(hu) i+1/2− i enom i h 2+ = R i ∂t ∂x

! (3.6) 2 with, for a spatially discretized function s, ∂(hu) ∂ h ∂Z 2 2 −1/3 + hu + g = −gh − gk h |u|u.

∂t ∂x 2 ∂x

D s = minmod(D s , 2θ D s )enom i eno i enom mm i

The model (3.6) can be written into a conservative form

where

∂U ∂F(U) + = S (U)+ S (U), (3.7)0 f min(x, y) if x, y> 0∂t ∂x

minmod(x, y)= max(x, y) if x, y6 0where 0 otherwise, ! ! hu h h 2 U = = , F(U)= ,h 2 s − s Δxhu q i i−1hu + g 2D s = minmod +θ D s ,eno i eno i−1/22 Δx 2! ! s − s Δx R i+1 i 20 −θ D s , eno i+1/2 S (U)= ∂Z and S (U)= .0 f 2 −1/3 Δx 2 −gk h |u|u !−gh

∂x s − 2s + s s − 2s + si+1 i i−1 i+2 i+1 i2D s = minmod , ,i+1/2 2 2System (3.7) is discretized using the ﬁnite volume method for Δx Δx

s − s s − si i−1 i+1 ihyperbolic conservation laws. We introduce a space-time grid D s = minmod ,mm i

Δx Δxin which the space and the time steps are, respectivelyΔx and

nΔt. We set x = iΔx, t = nΔt and C = x , x . We with θ ,θ ∈ [0, 1]. Note that for θ = 0, this reconstruc-eno enom enoi i i−1/2 i+1/2

n ndenote by U the approximation of the average of U(t , x) over tion is exactly the usual MUSCL reconstruction.i

the cell C , namely, To take into account the topography while preserving the steadyi

state of a lake at rest, i.e.,Z

1n nU ≃ U(t , x)dx.i Δx h+ Z = cst and u= 0,Ci

Considering only the homogeneous part of (3.7), then the ﬁnite we use a hydrostatic reconstruction described elsewhere [23,

volume scheme is of the form 24]. First, we need to deﬁne the reconstructed values Zi+1/2−

and Z that can be deduced from the reconstructed valuesi−1/2+Δtn+1 n n nU − U + (F − F )= 0,i i i+1/2 i−1/2 h , h and the following reconstruction of Z+ h:i−1/2− i−1/2+Δx

n n n Δxwhere F = F (U , U ) is the HLL numerical ﬂux (e.g.,i+1/2 i i+1 (Z+ h) = Z + h − D (Z + h )i−1/2+ i i enom i i

[24]) through the interface between C and C . Note that the 2i i+1

HLL ﬂux is deﬁned by

and

Δx

F(U ) if 0 < c ,l 1 (Z+ h) = Z + h + D (Z + h ). i+1/2− i i enom i i 2 c F(U )− c F(U ) c c2 l 1 r 1 2 Next the hydrostatic reconstruction requires the following new + (U − U )r l

F (U , U )= c − c c − cl r 2 1 2 1 values to be deﬁned: if c < 0 < c , 1 2 h = max(0, h + Z − max(Z , Z )),F(U ) if c < 0, i+1/2 l i+1/2− i+1/2− i+1/2− i+1/2+r 2

4 n+2 h h = max(0, h + Z − max(Z , Z )), i−1/2 r i+1/2+ i+1/2+ i+1/2− i+1/2+ in+2 e • Compute U = by solving the Manning fric- i ! ! n+1 n+2h eu

i ih hi+1/2 l i−1/2 r

U = , U = . tion termi+1/2 l i−1/2 r h u h ui+1/2 l i+1/2− i−1/2 r i−1/2+ n+2h i ∗∗ n+2 ∗∗ = S f (U ).The positive parts in the deﬁnitions of h and h (we eu − u Δti+1/2 l i−1/2 r i i i ∗∗ hihave max(0,·)) ensure the positivity of the water height. There- Δt

fore, the scheme can be written into the form

n+2• Compute U using the Heun method deﬁned by (3.8).

i

Δtn+1 n n n nU − U + (F − F − Fc )= 0,i i i+1/2 l i−1/2 r i We now describe the treatment of the boundary conditions atΔx

the inﬂow and the outﬂow. We denote by b subscript the values

where on the (ﬁctive) boundary cell and by the index “in” the values in

the ﬁrst cell inside the domain. The normal n is equal to -1 on

0 the left boundary (x= 0) and 1 on the right boundary (x= ℓ).n n n gF = F (U , U )+ , n 2 n 2i−1/2 r i−1/2 l i−1/2 r (h ) − (h )i−1/2+ i−1/2 r2 • A solid wall is modeled by imposing u = −u , h = h ,b in b in

on the condition that the topography be extended horizon-

0 n n n tally on the ﬁctive cells. g F = F (U , U )+ , n 2 n 2i+1/2 l i+1/2 l i+1/2 r (h ) − (h )i+1/2− i+1/2 l2 • At the inﬂow boundary, we impose the discharge q sat-b

isfying nq < 0. Because we only consider ﬂuvial ﬂows,and b

the water height h is computed using Riemann invariantsb p

0 (e.g., [29, 30]). More precisely, assume that c = gh.n gFc = . n n n ni − h + h Z − Z It is well known that for the Shallow Water system (3.6),i−1/2+ i+1/2− i+1/2− i−1/2+2

the quantity u ∓ 2c is constant along the characteristic

dxn = u∓ c. Thus we haveThe term Fc is added to obtain a well-balanced and consistent dti

scheme [23]. Now, to have a second order scheme in time, we

u + n× 2c = u + n× 2c . (3.9)b b in inuse the Heun method,

Multiplying (3.9) by h , we obtainn+1 bn neU = U +ΔtΦ(U )i i i √ 3/2

n+2 n+1 n+1 −n× 2 gh + (u + n× 2c )h − q = 0.e e e in in b bU = U +ΔtΦ(U ) bi i i (3.8)

n n+2eU + U Newton method is used to solve the last equation and ton+1 i i

U = ,i get h .2 b

• At the outﬂow boundary, we always impose the waterwhere

1 height h and we use Riemann invariants to compute then n n n bΦ(U )= − (F − F − Fc ).i i+1/2 l i−1/2 r iΔx discharge. We obtain

Concerning the Manning friction term, we introduce a semi-

u = u + n× 2(c − c )b in in bimplicit treatment of this term [21, 29]. Then, the scheme is

modiﬁed as follows and we can easily deduce q = h u .b b b

• Solve the Shallow Water system All the previous steps are usual for the resolution of the Shal-

low Water system (2.1). Note that, for our new 2D model (2.4),

∗ n nU = U +ΔtΦ(U ). the additional friction term is treated in an explicitly. As thei i i

obtained results are satisfactory, we did not try to use a semi-

n+1 implicit treatment of the additional friction term. Keeping the h in+1 e • Compute U = by solving the Manning fric- explicit treatment also allows the numerical model to be ex-i n+1 n+1h eui i tended easily to more complex two-dimensional ﬂows.

tion term

∗n+1 h h i 4. Numerical results i ∗ n+1 ∗ n+1 ∗ = S f (U ) ≡ eq |q | . eu − u Δt i i i 2 i i ∗ −gk In Section 2, we introduced several models: the usual Shal- h n n+1i 4/3h (h )Δt i i low Water system (2.1), the Shallow Water system with an ad-

ditional friction coeﬃcient that represents the furrows when we

• Solve the Shallow Water system consider plane topography (2.4) and a two-layer model (2.5).

In this section, we present several results obtained with these

∗∗ n+1 n+1e eU = U +ΔtΦ(U ).i i i three models in order to show the ability of our new model (2.4)

5 0.05to approximate the exact solution. Namely, we consider two

topography

water heighttypes of test cases: in the ﬁrst one, we only take into account

rainfalls, and in the second one, the water comes from up- 0

stream. For these two numerical experiments, the “exact” (or

reference) solution is that from the Shallow Water system (2.1)

-0.05

with a precise description of a topography with furrows. The

domain Ω we consider here is Ω = ℓ × L, where ℓ = 0.2 m

-0.1

and L = 4 m (Figure 2). We assume that the plane topography

has a constant slope of 5%. The amplitude of the furrows is

-0.15

0.01 m (0.02 m peak-to-peak), and their wavelength is 0.1 m.

1/3 −1We choose a friction coeﬃcient k = 0.04 m s . For the

-0.2following computations, we use a time stepΔt = 0.001 s (this

time step is imposed by the resolution of (2.1), because we need

a small space step to get a good representation of the furrows). -0.25

0 0.5 1 1.5 2 2.5 3 3.5 4

y (m)We note that all the numerical results are obtained using C++

software for the resolution of the Shallow Water system, and a

Figure 5: Side-view of the water height at the ﬁnal time for the rainfall test case.

new library for the new friction coeﬃcient.

0.05

topography

water height

4.1. Reference solutions and description of the test cases

0

This paragraph describes how the solutions of the Shallow

Water model (2.1) are computed when the geometry of the fur- -0.05

rows is known explicitly. These solutions will be considered

here as reference solutions. According to the parameters given

-0.1

above, the topography is modeled by the equation:

-0.15

Z(x, y)= −0.05 y+ 0.01 cos(20π y). (4.10)

-0.2

The space steps (with respect to x and y) are equal to 0.01 m,

which means that each furrow is described by 200 cells. We

-0.25assume that the domain is initially empty, i.e., u(0, x) = 0 and 0 0.5 1 1.5 2 2.5 3 3.5 4

y (m)h(0, x)= 0.

Let h, u, and q denote these reference solutions at the small

Figure 6: Side-view of the water height at the ﬁnal time for the inﬂow test case.

scale.

4.2. Numerical comparisons of the models

4.1.1. Rainfall test case

We perform numerical tests on the new model (2.4). TheIn this case, we impose rainfall on the whole domain with

−4 −1 furrows are removed from the topography deﬁned by (4.10).a constant permanent rain intensity R = 8 × 10 m s . The

−3 2 −1 Thus, the topography is now reduced to an inclined plane withrain discharge is then Q = 3.2 × 10 m s . The ﬁnal timeR

the same general slope:is T = 22.5 s. Note that because we are interested in the ef-

fects of the furrows, we focus on the transitional stage of the

Z(x, y)= −0.05 y. (4.11)ﬂow. Therefore, the ﬁnal time T is chosen such that the outﬂow

discharge is approximately equal to half of the rain discharge.

The space step with respect to y is set equal to the wave-We assume here that the upstream boundary is a solid wall. We

length of the furrows, i.e., 0.1 m. The initial conditions remainshow the side-view of the water height at the ﬁnal time in Fig-

unchanged: u(0, x) = 0 and h(0, x) = 0. The capital lettersure 5 .

(H, U and Q) denote the solutions of (2.4) at the large scale.

To improve the comparison in the rainfall test, we also com-

pute the solution (h,u,q) of the two-layer system (2.5). The dis-4.1.2. Inﬂow test case

cretization parameters and the topography we use are the sameWe also consider a permanent inﬂow from upstream. We set

−2 2 −1 as for the system (2.4).Q = 3.132× 10 m s as discharge on the inﬂow boundary.I

To compare the three models (2.1), (2.4) and (2.5), we con-The ﬁnal time is T = 27.75 s. As for the rainfall test case, the

sider the water height (h, H and h, respectively) and the dis-ﬁnal time was chosen such that the outﬂow discharge at T was

charge (q, Q and q, respectively) at the outﬂow. For this pur-approximately equal to half of the inﬂow discharge. We show

i

nthe side-view of the water height at the ﬁnal time in Figure 6. pose, we ﬁrst introduce h as the average of the reference water

6

Z, Z+ h (m) Z, Z+ h (m)nheight h contained in the furrow i at the time t . We also con- which corresponds to K = 0.02 and C = 0.4. The correspond-0

n n Q −2sider H , the water height in the furrow i at time t , which is ing discharge error is e ≃ 5.8× 10 . We notice that the newi

2computed with the model (2.4) for a given K and a given C, model (2.4) allows the L error on the water height to diminish0

nand H in the case K = 0 (i.e., the Shallow Water system by a factor of 7 with respect to the case K = 0, which indicates0 0 0i

on the coarser grid with plane topography and without the new that the furrow eﬀects are taken into account satisfactorily.

Hfriction term). Next, we denote by e the relative water height We now report the results obtained with the two-layer

error deﬁned by model (2.5), for the rainfall test. For diﬀerent Manning coef-

hﬁcients k, Figure 8 shows the water height error e (obtained 1/2N XX 2 i