{'fraction_non_alphanumeric': 0.04343060287364835, 'fraction_numerical': 0.023285439194193455, 'mean_word_length': 3.8548827844096074, 'pattern_counts': {'":': 0, '<': 6, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 22, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Angle dependence in coupling conditions for shallow water equations at canal junctions

In this paper we propose a numerical Riemann problem solver at the junction of one dimensional shallow-water canal networks. The junction conditions take into account the angles with which the channels intersect and include the possibility of canals with different sections. The solver is illustrated with several numerical tests which underline the importance of the angle dependence to obtain reliable solutions.


Introduction
The shallow water equations are a model to describe free surface water flows. They are a non-linear hyperbolic system of PDEs consisting of a mass and momentum balance. They are used to describe flows in artificial canals and water channels with applications for instance to environmental problems. In water management problems, these equations are often used as a fundamental tool to describe the dynamics of networks of canals or of the branching of rivers. Networks occur in different type of configurations. The most straightforward treatment from a numerical point of view consists in considering the network as a two dimensional domain covered with an unstructured grid [10,11]. However from a computational point of view it is much more efficient to consider the network as a set of one dimensional canals coupled through junctions.
The main difficulty mathematically is the definition of the coupling conditions at the junction between the adjoining channels. We mention the following reviews for one-dimensional flows on networks [4,13]. The coupling condition can be seen as a Riemann problem involving a constant state for each of the adjoining canals. Riemann problems at a junction are widely discussed in literature, see [14,8,19,31,29]. To 1 close the problem one completes the Riemann problem with physical conservation properties across the junction, see [27,26].
From a numerical point of view, one has to couple one dimensional numerical solvers in the 1D channels with an approximate junction Riemann solver, see [1,3,6,24,32].
We consider a junction of three canals and we assume by convention to have one incoming canal which ends at the junction and two outgoing canals which start at the junction. To solve the junction problem we need to find the three states (mass and discharge) facing each of the three one dimensional channels at the junction for a total of six unknowns. One imposes mass conservation at the junction which yields one equation, then one formulates a left-half Riemann problem for the incoming canal and a right-half Riemann problem for each of the two outgoing canals. Under subcritical flow assumptions, we obtain therefore three more conditions. Thus, two remaining equations have to be specified in order to define the junction model. In some works the set of equations is completed by assuming the continuity of water levels [8,19,31,5] or the continuity of energy [27,24]. However, none of these works use a condition that takes into account the geometry and especially the angles formed by the channels in the fork. For an attempt to include an angle dependency in the solver, see [7,25,16]. Different approaches covert the junction with two-dimensional elements and project the computed 2D solution along the one dimensional channels [1].
These studies have several applications, such as optimization [29,22,20,21] or see [2] for a nice application to the modeling of a particular wave energy converter, the so-called oscillating water column.
In this article, we propose new coupling conditions at the junction that depend on the angles with which the channels intersect at the junction allowing also for channels with different sections, [12,28]. Away from the junction we assume the solution to be 1D, while we describe the junction as a 2D region where coupling occurs between the branches. We then consider the triangle formed by the intersection points of the walls of the three channels and, to this two dimensional domain, we apply conservation of mass and of the two components of momentum. We obtain three non linear equations which include a dependence on the angles, for the six unknowns at the junction to be coupled with the three equations of the characteristic curves. Extending this study to a network of canals with several nodes is straightforward. This work therefor extends the results of [7] by considering branches with different sections.
We prove that the conservation condition based on the continuity of the energy across the junction occurs in a particular configuration of our setting. We also prove the existence of the solution of our junction Riemann problem in a few particular cases.
Validation of numerical schemes obtained in this way is carried out comparing the numerical 1D solution with the junction, with a fully 2D solver, see [16,18,23,17]. We compare our numerical solver with a fully 2D solver for shallow water equations showing that the numerical approximation improves as the width of the 2D channels is reduced.
The paper is organized as follows. In Section 2, we concentrate on the solution of the Riemann problem for shallow-water equations. We then present the junction geometry in Section 3, defining our coupling conditions in Section 3.1. We discuss extensions for special configurations in Section 3.2, this includes the case of a single channel with varying cross-section. In Section 4 we merge the relations at the junction with the numerical approximation of shallow water equations along the channels. Section 4.3 is devoted to a discussion of the existence of the numerical solution in a few cases. We end in Section 5 with the numerical tests.

The shallow water equations and its standard Riemann problem
Let us first recall the shallow water or Saint Venant equations, and some of theirs properties that will be useful in the following.
2.1. The shallow water equations. The 1D shallow water equations, introduced by Saint-Venant in [9] and derived in [15] from Navier-Stokes incompressible equations with a free moving boundary, describe the water propagation in a canal with rectangular cross-section and constant slope as follows: the water height, v(x, t) the water velocity at time t and location x along the canal, g the gravity constant, S 0 the bed slope function and S f the friction slope function. The first equation comes from mass conservation and the second one from momentum balance. For the purpose of this work, we assume a steady state friction on all canals and we assume horizontal canals with zero slope. Thus, the source term is zero.
We set q = hv (the quantity hv is often called the discharge in shallow water theory, since it measures the flow rate of water past a point) and we reformulate system (1) in vector form as For smooth solutions, system (2) can equivalently be written in the quasilinear form with eigenvalues Note that in general λ 1 and λ 2 can be of either sign. When the velocity v = q/h of the fluid is smaller than the speed √ gh of the gravity waves, that is |v| < √ gh, the flow is said to be fluvial or subcritical and then one has Hence, under the subcritical condition (7), there are two waves propagating in opposite directions. The left and right characteristics are associated to λ 1 and λ 2 respectively. The ratio F r = |v|/ √ gh is called the Froude number and the flow is subcritical iff F r < 1.

2.2.
The standard Riemann problem for shallow-water equations. Here we are in particular interested in the solution of the Riemann problem: where U (x, 0) = (h(x, 0), q(x, 0)) is the initial condition and U l = (h l , q l ) (resp. U r = (h r , q r )) is the initial constant state on the left (resp. on the right) of the interface x = 0. The characteristic fields of the shallow water equations are genuinely nonlinear and so the Riemann problem always consists of two waves, each of which is either a shock or a rarefaction. Under the subcritical flow condition (7), there will be one left (with negative speed) and one right (with positive speed) going wave. In the sequel the left and right going waves are denoted by l-wave and r-wave, respectively. The solution to this Riemann problem consists of the l-wave and the r-wave separated by an intermediate stateÛ = (ĥ,q). We remark that the solution at the interface x = 0 coincides withÛ , which is the intersection point of the two functions φ l and φ r defined by which return the physically correctĥ andv intermediate values connecting the left and right states with an entropic solution.

Angle dependent conditions at the junction
In this work, a junction is defined as the intersection of three channels. We assume that 1D shallow water equations hold on each canal of the network and we aim at deriving coupling conditions at the junction. These conditions enable to compute the intermediate states at the junction for the Riemann problem under consideration.
3.1. Definition of the coupling conditions at the junction. The channels will be labeled 1, 2 and 3 respectively, where channel 1 is assumed to be parallel to the x axis. We fix the origin of the reference system in the point where the three channels intersect. Let θ and φ be the angles that channel 3 and 2 respectively form with the x axis. We will assume that θ ≥ 0, while φ ≤ 0, obtaining the geometry in Fig. 1 on the left. This is the one dimensional set up.
Further, we will suppose that the channels can have different widths. Let then 2s j , j = 1, 2, 3 be the width of each channel. Therefore, we can think that the 1D setup is the core of a two dimensional junction, as shown in Fig. 1 on the right.
Let I k , k = 1, 2, 3, be the interface separating the k-th channel from the junction. Let U * k , k = 1, 2, 3, denote the state variable in channel k at the side of I k facing the channel obtained with the 1D solver used in the canal, while U k , k = 1, 2, 3, is the state variable at the side of I k facing the junction. The purpose of the junction Riemann solver is to compute U k given U * k . Since each state consists of the couple (h, v), we need to find 6 unknowns at the junction. Three conditions are obtained finding the intermediate states of the one dimensional Riemann problem defined at each interface I k and in order to compute the three other missing data we shift to the 2D setting of Fig.1 on the right. We consider the triangle formed by the intersection points of the walls of the three channels, and to this two dimensional figure we apply conservation of mass and of the two components of momentum, which gives us the 3 missing equations. Once the three states U k , k = 1, 2, 3, at the junction have been computed, we have at each interface I k the left and right states which are needed to compute the numerical flux at the boundary interfaces of the channels.
3.1.1. Junction conditions coming from the Riemann solver. Let us begin with the 3 equations coming from the Riemann solver. We emphasize that, by convention, the given configuration fixes channel 1 as entering the junction and channels 2 and 3 as leaving the junction.
Three relations are obtained matching the unknowns U k at the junction with the data U * k coming from the three channels through equations (9) and (10). More precisely, . We note again that this construction requires a fluvial regime, in which only one wave exits the junction towards each of the three adjoining channels.
3.1.2. Junction conditions coming from mass and momentum conservation. We now derive the 3 supplementary equations coming from conservation of mass and of the two components of momentum. For that purpose, we come back to the 2D configuration of the junction and we use the following notations • h denotes the height of water in the 2D configuration, • v = (v x , v y ), denotes the 2D velocity in the 2D junction domain, see Fig.1, , v y,k ), k = 1, 2, 3, denotes the average discharge on the edge of the junction triangle corresponding to channel k, see Fig.1. We first recall the shallow-water equations in 2D, composed of the mass conservation equation and of the momentum conservation equation: In the following, we call T the triangle formed by the intersection points of the walls of the three channels and its boundary ∂T is composed of three edges, denoted by e k , k = 1, 2, 3, see Fig.1.
Mass conservation across the triangle T with boundary ∂T yields where n is the outer normal of ∂T , while the conservation of the two components of momentum gives the two relations (14) ∂T ∂T v y q + 1 2 Decomposing ∂T as the sum of the three edges e k , k = 1, 2, 3, the three conditions at the junction, given by mass conservation and the two components of momentum conservation, can then be written as: where k is the length of the edge e k of the triangle, n k is the outer normal to e k and q k is the average of q on the side e k of the triangle. To specify all quantities appearing in system (16), we need to compute the normals n k to the sides of the triangle and their lengths k . To fix notation, refer to Fig. 2.
To begin with, we need to give the coordinates of the intersection points of the walls, namely points P 12 , P 13 and P 23 that are displayed on Fig. 2. Let us recall that θ and φ are the angles of canals 2 and 3 with the x-axis, while 2s k is the section of canal k.
The equations for the straight lines composing the 1D skeleton of Fig.1 written in parametric form are with t k ∈ R, k = 1, 2, 3. Then, to obtain the walls of the channels, i.e, to construct the 2D setting of Fig. 2, we just need to write the equations of the two straight lines parallel to the axis y k at the center of the channel, and at a distance ±s k from the axis, for each channel. The walls of the three channels are Fig. 2 across which the 2D interaction occurs is obtained intersecting the straight lines defining the walls of the channels. More precisely, P 13 , is the intersection of y + 1 with y + 3 , P 12 is defined by the intersection of y − 1 with y − 2 , and the last point P 23 lies at the intersection of y − 3 and y + 2 . We obtain, If θ = 0, the system has a solution only provided s 1 = s 3 , and the two straight lines actually coincide. In this case we define P 13 = (0, s 1 ). Analogously, If φ = 0, we must have s 1 = s 2 , and we pick P 12 = (0, −s 1 ). With this approach, we cannot treat the case in which both φ = θ = 0, unless we consider the two channels y 2 and y 3 superposed one on top of the other. We will see in the next subsection how to extend the construction also to the case φ = θ = 0. Finally, The quantity sin(φ − θ) can be zero either for φ = θ = 0, in which case the two pipes coincide, or when φ = −π/2 and θ = π/2. Now you have solutions only for s 3 = s 2 , which means that y − 3 and y + 2 coincide, and we fix the intersection point to P 23 = (s 3 , 0).
Once the points P 13 , P 12 , P 23 are defined, we can compute all quantities n 1 , n 2 , n 3 , 1 , 2 and 3 depending on the geometry appearing in (16). The length of the sides is and the normals are Remark 1. The construction is well defined as long as the triangle formed by P 13 , P 12 , P 23 is non degenerate. We say that the triangle is degenerate when the three points lie on the same straight line. Straightforward computations show that this occurs when det(n 1 , n 3 ) = 0 which is equivalent to the particular combination In the frame of reference we have chosen, the discharge in the three canals can be written as where q k = q k .
Let v k = q k /h k be the velocity along the k−th channel. Then the conservation laws (16) across the junction can be written as where we used the fact that the axis of channel 1 is parallel to the x axis.
3.1.3. Solutions for the whole system of equations at the junction. Combining the three equations (24) with the three equations (11) coupling the states in the junction with the 1D channels, we find a system of 6 non linear equations at the junction, whose solution is given by the three intermediate states Remark 2. If we consider a stationary solution of (1) such that the velocity is null and the height is constant in space, i.e. h * k =h and v * k = 0, k = 1, 2, 3 then h k =h and v k = 0 k = 1, 2, 3 is a trivial solution of system (24)-(11) since n 1 + n 2 + n 3 = 0. This means that the coupling condition at the junction preserves the lake at rest stationary solution on the whole network. (24), we find a system of three non linear equations in the three unknowns h 1 , h 2 and h 3 which gives the solution at the junction. Once the parameters s 1 , s 2 , s 3 , θ, φ and U * 1 , U * 2 , U * 3 are fixed, these three equations define three hypersurfaces whose zeros surfaces can be plot in the plane h 1 , h 2 and h 3 ( Figure  3). The intersection of this surfaces is the required solution. In Figure  3-right we show an example.
3.2. Special cases and extensions. In this section, we consider three particular cases.
We start with the simplified case in which the canals are orthogonal to the sides of the triangle. In that case, the junction is defined uniquely by the three sections and the equations (24) simplify loosing the dependency on the angles. In this case, see Fig. 4, it is easy to see that the angles θ and φ defining the skeleton of the junction coincide with the angles labelled θ and φ internal to the triangle in Fig. 4, and the length of the sides coincides with the amplitude of the channels, namely l k = 2s k , k = 1, 2, 3.
Then, it is straightforward to see that the sections depend on the angles through the following relations (25) Since in the present case, q k is parallel to n k , equation (16a) becomes Equation (16b)-(16c), corresponding to the conservation of momentum at the junction in 2D give:  Using the identities in (25), we can rewrite (27) and (28) as Therefore, conservation of mass and of the two components of momentum at the junction in this particular case yield Note that in this case, the junction conditions do not depend on the angles with which the canals intersect. We thus recover the equal energy condition at the junction used by several authors, see [24] and references there in. This condition derives from the 2D momentum conservation at the junction, but we stress that it holds only for the particular case in which the channels are orthogonal to the sides of the triangle defining the junction. Moreover, tedious but straightforward calculations show that equations (25) imply the geometry in Fig. 4. Since the conservation condition can be multiplied by a constant without changing the result, we see that for each pair of angles φ and θ there exists a one parameter set of sections λ(s 1 , s 2 , s 3 ), with λ > 0 for which momentum conservation coincides with energy conservation. We now consider the cases θ = φ = 0 and θ = −φ = π/2 which where excluded in the generic case described in section 3. We call : • T-junction: θ = −φ = π/2. • Straight channel: θ = φ = 0.

T-junction.
In the case of a T-junction, for which the angles are equal to θ = −φ = π/2 and s 2 = s 3 , the points P 12 , P 13 and P 23 can be defined as in Fig.5-(A), namely Equations (24) reduce to When s 2 = s 3 , the point P 23 can be defined as in Fig.5-(B), then 3.2.2. Straight channel. Now, we consider the case when θ = φ = 0. There is a natural way to define points P 12 and P 13 , see Fig.6. By symmetry, the y-coordinate of P 23 should be set to 0 but the x-coordinate is undetermined. We fix P 23,x = s 1 , thus and equations (24) reduce to: Note that, with this configuration we can extend our construction to the case of a single channel with a varying cross section.

Numerical scheme for shallow-water equations complemented with junction conditions
In this section, we couple a standard final volume scheme for the shallow water equations along each channel with the numerical flux consistent with our junction conditions (11)-(24). 4.1. One dimensional finite volume scheme. For the sake of simplicity, we will suppose that each canal has the same length, discretized with a uniform grid. Then, the computational domain in each canal is defined by the finite interval [0, L], which is divided in M equal cells, of amplitude ∆x = L/M . The cell centers are given by x j = (j − 1 2 )∆x, j = 1, . . . , M , and the cell average of the numerical solution at time t in the j-cell is defined as (30) U j (t) = 1 ∆x with appropriate boundary conditions for U 0 (·) in channel 1 and U M +1 (·) on the two outgoing channels 2 and 3. The system is evolved until the final time T , with time step ∆t. We denote by U n j the approximate value for the average of U in cell j at the discrete time t n = n∆t. Hence, the finite volume approximation of system (2) can be written under the form with F (·, ·) a proper numerical flux. We apply a Godunov type numerical flux [30,33], computing the exact intermediate state for the Riemann problem between the two cells defined by U l = U n j−1 and U r = U n j for j = 1, . . . , M . From now on we shall add a canal index and we use the notation U n j,k , j = 1, . . . , M , to indicate the numerical solution along canal k = 1, 2, 3 computed at time t n , in x j,k . The time step ∆t is fixed to satisfy the stability condition (33) ∆t ≤ ∆x max Let U k = (h k , v k ), k = 1, 2, 3 be the solution to system (11)-(24) with U * 1 = U n M,1 , U * 2 = U n 1,2 , U * 3 = U n 1,3 where U n M,1 , U n 1,2 and U n 1,3 are the values computed by the 1D scheme along the channels. Then, at the junction points we imposeF n M +1/2,1 = f (U 1 ),F n 1/2,2 = f (U 2 ) and F n 1/2,3 = f (U 3 ).

4.3.
Solving the non-linear system at the junction. Now, let us study the solutions to the nonlinear system at the junction. We recall the notations used in Sec.3 and denote by h * k and v * k the approximate values of h and v near the junction at channel k = 1, 2, 3 given by the 1D numerical scheme, see also Sec.4 for their exact definition. Let Ω be the open set of admissible states, Ω = {h k ∈ R + * , v k ∈ R, |v k | < √ gh k , k = 1, 2, 3}. The approximate values h k and v k are then obtained solving the non-linear system (24) with v k given by (11). Let us denote We can rewrite the system under the following form: where Ψ : Ω × Ω → R 6 . In general, existence and uniqueness results for solutions to non linear systems are difficult to prove. Assume that we have solved the system up to t = t n , this gives the solution X n, * which faces the junction. Suppose that we have found a solution X n such that Ψ(X n ; X n, * ) = 0. If we can prove that Det DΨ (X n ; X n, * ) = 0, where DΨ denotes the Jacobian with respect to the first argument, then there exists a unique X = X(X * ), for X * − X n, * < , with small enough, such that Ψ(X; X * ) = 0. Therefore, if the flow is smooth, one can find ∆t small enough such that X n+1, * − X n, * < and the implicit function theorem guarantees that there exists a unique solution X n+1 such that Ψ(X n+1 ; X n+1, * ) = 0. So, the procedure can be iterated provided one can prove at each step that Det DΨ (X n ; X n, * ) = 0.
In the particular case when all waves are rarefactions, the relations (11) become (36) and it is clear that the Jacobian DΨ does not depend on the data X * . Thus, starting from a set of data X * and a solution X such that Ψ(X; X * ) = 0, once one can prove that (37) Det DΨ (X; X * ) = 0, for all X * , X ∈ Ω × Ω, the solution exists at each time step. Note that, for the steady solution h = const. and v = 0 one has X * = (h, h, h, 0, 0, 0) T and Ψ(X; X * ) = 0 for X = X * , thus there exists at least one case for which Ψ(X; X * ) = 0.

Example 1.
In the case when the canals are orthogonal to the sides of the triangle, we have Tedious algebra gives, We can therefore conclude that since we are in the sub-critical case, for which we have Det DΨ (X; X * ) < 0, for all X * , X ∈ Ω, which implies condition (37). Thus, starting from a point for which Ψ(X; X * ) = 0 we can prolong the solution for all of times. This coincides with the case in which one assumes the continuity of the energy.
Example 2. Case with vanishing velocities. Now, consider the case given by system(24)- (11), with only rarefaction waves. In order to simplify the expressions arising in the computations, we introduce the following notations Since explicit computations are too difficult, we restrict ourselves to the particular case when solutions with vanishing velocities v n 1 = v n 2 = v n 3 = 0 at the junction exist, that is to say X n = h n 1 h n 2 h n 3 0 0 0 . We already note that n 1 + n 2 + h, h, 0, 0, 0) is a trivial solution of system (24)- (11).
We prove that in this case Det DΨ (X; X * ) = 0. In fact, We notice that this expression is the product of two terms: the first one, The second term depends only on s 1 , s 2 , s 3 , θ, φ, that is to say on the triangle geometry. This second term vanishes iff the triangle is degenerate, see Remark 1.

Numerical tests
In order to evaluate the effectiveness of our method we perform various tests consisting of a subcritical wave propagating across junctions of different geometries. For that purpose, we will use the numerical scheme presented in the previous section. First of all we check the numerical convergence of the scheme on the whole network under grid refinement. Then, we compare numerical simulations on a network with θ = φ = 0, with simulations on a single canal, to show the consistency of the junction conditions with the traditional one-dimensional Riemann solver. Subsequently, we increase θ and φ in order to enhance the influence of the angles in our junction conditions. We also investigate numerically the case when water flows out canals 1 and 3 and pours into canal 2. Finally, we compare the dynamics of the 1D solver with the numerical solution obtained with a two-dimensional code. 5.1. Convergence of the 1D numerical scheme. In this test case we check the numerical convergence of the 1D scheme (31) coupled with system (24)-(11) at the junction under grid refinement. We set the geometry parameters s 1 = s 2 = s 3 = 1, θ = −φ = π/6, L 1 = L 2 = L 3 = 3.5, where L k is the length of channel k. The initial data in the three channels is We expect the formation of a rarefaction wave propagating backwards on channel 1 and a shock crossing the junction and travelling with positive speed along the two outgoing canals. In Figure 7 we show the solution obtained after the water wave has reached the junction. The dynamic involves the three canals and we observe the convergence of the numerical solution under grid refinement. The number of grid points on each channel is N = 12, 24, 48, 96. We observe that the solution of system (24)-(11) at the junction does not depend on the grid parameter N , proving the consistency of the Riemann solver. We observe the formation of a stationary shock at the junction.

5.2.
Comparison with the solution on a single canal. Here we study the influence of our junction conditions involving the angles θ and φ on the solution. We set s 1 = 1 and s 2 = s 3 = 1/2, L 1 = L 2 = L 3 = 5 and the initial data are the same as in (38). In all tests we fix the grid parameters ∆x = 0.01 on each channel and ∆t to satisfy (33).
We first consider θ = φ = 0 which corresponds to a single channel. As expected for θ = φ = 0 and s 1 = s 2 + s 3 the solution of our algorithm coincides with the solution computed on a single canal, i.e. without the junction, see blue line and black dashed line on Figure 8. Next, we change θ and φ to study the influence of the angles on the dynamics. Specifically, in Figure 8 the angles vary symmetrically with θ = −φ and θ = 0, π/12, π/6, π/3. We observe that the symmetry of the configuration is preserved and that the solution varies monotonically increasing the angles and moving away from the single channel profile. Note that to see the dependence of the solution on the angles, it is essential to include an angle dependence in the junction condition.
In Figure 9 we study the influence of a non symmetric variation of the angles: we fix θ = π/8 and vary φ = −kπ/8 with k = 0, 1, 2, 3. We observe that the symmetry of the solutions of the two outgoing channels is lost and that water meets more resistance as channel 2 becomes more and more bent. So, the water level decreases in channel 2 and increases in channel 3. The solution on channel 1 does not change because the total lumen of the outgoing channels remains the same.
Finally, we compare the 1D solution fixing the two angles and varying the channel sizes. Specifically, we fix θ = −φ = π/4, s 1 = s 2 = 1 and consider s 3 = 0.5, 1, 1.5, 2. In Figure 10 we observe that, as the section of channel 3 increases, the water height decreases in canals 1 and 3 and increases in canal 2. The water velocity in channels 2 and 3 follows the same behavior, while the velocity in channel 1 increases. However, the dynamics in channel 2 does not vary significantly.

Merging canals.
In this section we display numerical results for a 2-to-1 or merging junction for which the water flows from channels 1 and 3 towards channel 2. We set as initial data Then, we fix s 1 = s 2 = s 3 = 1, θ = π/3. In Figure 11, we compare the solutions obtained for φ = −π/3, −π/6, −π/12, 0. As the angle φ widens, the water height increases in channel 1 and decreases in channel 3. This asymmetry explains why the dynamic in channel 2 is almost unaffected by the angle variation, the amount of water entering remains almost constant.

5.4.
Comparison of 1D and 2D solutions. We compare our 1D solver with the 2D shallow water solution (12). The numerical solution of (12) has been computed by the free and open source ToolBox FullSWOF2D (Full Shallow-Water equations for Overland Flow in 2D), which is a C++ code for simulations in two dimensions [10]. We compute the solution on the rectangle in Figure  12 with (v x , v y ) · n = 0 on the boundary. To obtain the 2D geometry of the junction we are interested in, we use a bottom topography which is z in = −1 within the dashed region Ω and z out = 0 in the complement Ω C . We choose the initial water height so that h + z in is less than z out . In this way, the water flows only inside the dashed region while Ω C is seen as a dry state region.
To compare the results of the 2D solution with the 1D code, the values of h, hv x , hv y of the 2D solution are sampled on the straight lines at the middle of the channels (black solid lines in Figure 12) of lengths L k , k = 1, 2, 3. We compare the 1D velocity with the 2D velocity norm v 2 x + v 2 y . The mesh of the two-dimensional domain contains about 26 × 10 4 grid nodes.
We fix the canal lengths as L 1 = L 2 = L 3 = 5 and the simulation final time to T = 1.5. We set s = s 1 = s 2 = s 3 , and we compare the 1D solution obtained with junction conditions (24)- (11) to the 2D solution with s = 1, 0.5, 0.25. As s decreases, the 2D configuration becomes closer to the 1D network. Recall that, in the 1D case, if s = s 1 = s 2 = s 3 , the solution of the shallow water equations complemented with junction conditions (24)- (11) does not depend on the value of s.
We consider two different cases, a diverging junction, that it to say the case when water flows from the single channel 1 towards the two channels 2 and 3 and a merging junction, when the water flows from the two channels 1 and 3 towards the single channel 2. Diverging junction. We fix θ = π/3, φ = −π/12 and initial states as in (38) for the 1D configuration and for the 2D system such that h(x, y, t = 0) = 1.5 (x, y) ∈ Ω ∩ {0 ≤ x ≤ L 1 /2} 1 otherwise.
In Figure 13, the 1D solution with junction conditions (24)-(11) is presented in solid red, the dotted black curve is the 1D solution assuming equal energy at the junction, as in [24]. We have three 2D solutions which are displayed in dashed blue, magenta and green for s = 1, 0.5, 0.25 respectively. Decreasing s we see that the 2D solutions converge towards the 1D wave front, both in the two shocks and the receding rarefaction. The only differences are observed in the flat intermediate states at the junction. The junction solver proposed in this work seems more accurate than the one with equal energy condition at the junction. Some wiggles appearing in the 2D solution might be due to numerical artefacts at the interface between dry and wet states of the 2D code.
Merging junction. We fix θ = 5π/12, φ = −3π/8 and initial data as in (39) for the 1D configuration while for the 2D system as (40) with h = 1.5 in the right half of channel 3 too, see Figure 12 on the right. As before, in Figure 14, we superpose the 1D and 2D curves. Again, our 1D solution is displayed in solid red and the dotted black curve represents the solution with equal energy condition at the junction, [24]. The 2D solutions are drawn in dashed blue, magenta and green for s = 1, 0.5, 0.25 respectively. Decreasing s, we see that the 2D solutions converge towards the 1D water front, as expected.

Conclusions
In this paper we have presented a numerical solver for one dimensional channels in a network. The solver is based on a finite volume scheme in each canal and coupling conditions at the junction are obtained with a single 2D element at the junction across which mass and the two components of momentum are conserved. This approach allows to take into account quite general geometries including the dependence on the angles with which the canals intersect at the junction and the sections of the canals. In this framework we can also include the construction of solvers for shallow-water problems along a channel with varying section.
The solver is based on the assumption that the flow across the junction is fluvial. Future work on this topic will be concentrated on the case of torrential flows and on the dependence of the bottom topography. (a) water height