Simulation of Dam-Break Problem using Random Choice Method

Computers & Fluids 111 (2015) 187–196

Contents lists available at ScienceDirect

Computers & Fluids j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c o m p fl u i d

Simulation of Dam-Break Problem using Random Choice Method Himanshu Gupta, L.P. Singh ⇑ Department of Mathematical Sciences, Indian Institute of Technology (Banaras Hindu University), Varanasi 221005, India

a r t i c l e

i n f o

Article history: Received 20 June 2014 Received in revised form 24 September 2014 Accepted 4 February 2015 Available online 11 February 2015 Keywords: Dam-Break Problem Shallow water Random Choice Method FORCE method Hybrid method Riemann Problem

a b s t r a c t In the present study we investigate the application of Glimm’s Random Choice Method to the one dimensional homogenous shallow water equations of the Dam-Break Problem and extend the same to two-dimensional case using operator splitting procedure. A comparison is drawn between the Random Choice Method and the FORCE (First Order Centered) scheme. It was observed that the Random Choice Method gives sharper shock resolution as compared to FORCE method which causes the smearing of discontinuities. On further analysis of the two-dimensional case via operator splitting procedure, we observe that the result is highly erroneous and haphazard due to the presence of shocks oblique to the cell interfaces. We resolve this error by applying a Hybrid method to the two dimensional case; the Glimm’s method is applied to the continuous parts of the flow while the FORCE method is applied at large pressure jumps. This increases the accuracy of the results and provides good shock capture as well. Ó 2015 Elsevier Ltd. All rights reserved.

1. Introduction The hyperbolic quasi-linear system of equations governing the flow of an ideal incompressible fluid on almost horizontal free surface under the influence of gravity, based on the assumption that the depth of fluid is small as compared to the width of the flow, is known as shallow water equations [23]. They govern a large number of physical phenomenon in nature including the shallow coastal waters flows, tidal flows of ocean, river floods, shallow lakes, Dam break flows, etc. The Dam-Break Problem [1,2,3,11,12,19] is an example of the typical Riemann Problem of Fluid Dynamics (Shock Tube Problem) [13,20]. Just like the original shock tube, Dam-Break flows too have an initial discontinuity due to the existence of two separate states on either side of the Dam Wall, one representing the high level of water contained in the Dam while the other the lower surface water level. The problem is initiated with the failure of the Dam at t = 0 and owing to the resemblance with the shock tube, the solution of the Dam Break depicts characteristic shocks in pressure and velocity values (discontinuities) at some points while the flow is mostly smooth at the others. The non-linear nature of the shallow water equations inhibits the solution by analytical methods and in general, specific numerical difference methods are used to obtain the best results. There are several well established methods to handle shock discontinuities in shallow water models [4,6,9,14,22], including the ⇑ Corresponding author. E-mail address: [email protected] (L.P. Singh). http://dx.doi.org/10.1016/j.compfluid.2015.02.001 0045-7930/Ó 2015 Elsevier Ltd. All rights reserved.

exact Riemann solutions by Han [12] and the most prevalent upstream type ones which use either Flux-difference Splitting (FDS) or Flux Vector Splitting (FVS) [23] method. One such FVS scheme was presented by Baghlani [3]. Another category of schemes include those which utilize the Godunov [10,16,17] type Flux and are known as Godunov type methods. Toro used this flux to solve the two dimensional equations by using a Weighted Average Flux [21] scheme and then again to solve dry-bed problems [20]. More recently the finite volume methods particularly those of MUSCLE [23] type have gained popularity to solve the multi-dimensional cases. Zhou [24] used the MUSCLE scheme in combination with a Surface Gradient Method (SGM) to model the flow. Generally, the approximate solutions should meet the following criteria: i. The approximate solution must be reasonably accurate as compared with the exact one in terms of wave position, speed and shape both in the continuous and discontinuous parts. ii. The discontinuities (both linear and non-linear) in the flow should be depicted as sharply as possible by the approximate solution with little or no smearing across the cells. In this work, we apply the Random Choice Method of Glimm [7,8] and its Godunov type extension known as the FORCE [23,2] method to one dimensional and two dimensional Dam-Break Problems, and the results are compared. We also apply a Hybrid method involving a combination of RCM and FORCE method for

188

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

the two dimensional flow. Further, the numerical results for 1-D Dam Break flow have been validated with the analytical solution [18]. Simulation of Partial Dam Break Flow [5] have also been presented to test the capability of the code. 2. Governing equations and methodology The two-dimensional homogenous shallow water equations in Cartesian coordinate system is given by [23],

@U @F @G þ þ ¼ 0; @t @x @y

h i Ii ¼ xi1 ; xiþ1 ; 2

2

i;j;k

in which,

hv

huv

hv huv 2

2

hv þ gh =2

nþ1

 Dt  n F i1;j  F niþ1;j ; 2 2 Dx

2

2

2

nþ1

 Dt  nþ12 nþ1 Gi;j1  Gi;jþ21 ; Dy 2 2

2

2

2

2

! Dxi;j;k Dyi;j;k Dzi;j;k ; ; ; Sn;x Sn;y Sn;z i;j;k i;j;k i;j;k

ð3Þ

ð4Þ

n;d

n Sn;d i;j;k ¼ jdi;j;k j þ ai;j;k :

2

2

2 2

ð6Þ

ð7Þ

n;d

Here di;j;k is the velocity component in d direction at time n and ani;j;k pffiffiffiffiffiffi is the wave celerity at time n in the cell Ii;j;k given by gh where g is the acceleration due to gravity and h is height of water level. 3. Random Choice Method This method was first proposed by Glimm [8] in 1965 to solve a class of non-linear hyperbolic systems and was first applied successfully by Chorin [7] to solve the Euler equations of Gas Dynamics. This method provides fairly good results when the solution exhibits a pattern of shock and slip lines with very little computational effort as compared to most classical methods which generally employ high CPU time. Random Choice Method has two prerequisites: (i) solution of local Riemann Problem (ii) a random sampling procedure to select the states to be assigned to the next level. The interesting feature of this method is its unconditional stability and its ability to neglect all waves with speed greater than Dx=Dt, a feature which is difficult to achieve in traditional conservative methods, due to the fact that RCM picks up a state randomly from the local solutions and assigns it to the next state. The sampling range for our case is given in Fig. 3. The method for x direction is shown below, which can then be extended in y direction, for the two dimensional case, using dimensional splitting by changing the co-ordinates in computation procedure. Given an initial state Uðx; t n Þ at time t ¼ tn , we need to

ð5Þ

As for all Godunov type methods we discretize the domain ½0; L; ½0; K in x and y directions respectively into N, M grid cells

Fig. 1. Cross sectional view of water level.

Dx ¼ xiþ1  xi1 ¼ Dy ¼ yiþ1  yi1 ¼

sents the speed of the fastest wave at time n traveling in d direction where d = (x, y, z). The simplest choice is often

where G is the numerical flux at the cell interface j + 1/2 and j  1/2 given by

    Gniþ1;j ¼ G U niþ1 ð0Þ Gni1;j ¼ G U ni1 ð0Þ :

and

where C cfl is the CFL (Courant Friedrichs–Lewy) stability criteria

where U niþ1 ð0Þ and U ni1 ð0Þ represents the solution of local Riemann 2  2   Problem RP U ni ; U niþ1 and RP U ni1 ; U ni respectively along the line x=t ¼ 0. Likewise, for the y direction, the system is described by

U nþ1 ¼ U i;j 2 þ i;j

size

which varies from [0, 1] according to the method chosen, Sn;d i;j;k repre-

where U is the vector of conserved variables and F is the numerical flux at the cell interface i + 1/2 and i  1/2 given by 2

of

direction

7 5

ð2Þ

    F niþ1;j ¼ F U niþ1 ð0Þ F ni1;j ¼ F U ni1 ð0Þ ;

2

x

3

,h is height of water level, u and v are the velocity components in x and y direction respectively and g is the acceleration due to gravity. Due to homogenous form, the bed slope component and surface friction component are reduced to zero and are therefore not accounted for. U is the vector of conserved variables, F is the flux in x direction and G is the flux in y direction. A cross sectional view of the water level is shown in Fig. 1. This non-linear non homogenous equation is of hyperbolic type. We solve this equation by utilizing dimensional splitting in x and y directions i.e. first solving it in the x-direction using an explicit conservative method of Godunov type (for non-linear hyperbolic systems) which is described by

U i;j 2 ¼ U ni;j þ

for

K=M; j ¼ 1; . . . ; M in y direction respectively.Dx is the step size in x-direction, Dy is the step size in y direction as shown in Fig. 2, the superscript n; n þ 12 ; n þ 1 in (2)-(5) represent the advancement to the next time interval respectively. The time step Dt is dependent on (i) intercell flux of the particular method used (ii) the wave speed(s) present and (iii) the grid dimensions. It is defined as:

Dt ¼ C cfl  min

3 2 3 2 h hu 6 7 6 2 7 6 2 U ¼ 4 hu 5; F ¼ 4 hu þ gh =2 5; G ¼ 4

2

L=N; i ¼ 1; . . . ; N

ð1Þ

2

h i Ji ¼ yi1 ; yiþ1

Fig. 2. The grid discretization.

189

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

find the solution at time t ¼ t nþ1 . The Random Choice Method assumes a piecewise constant distribution of data across the entire cell Ii . In Godunov type methods, this is generally achieved by taking the integral averages across the cell, however for the Random Choice Method this is not required and we can safely assume that the value at xi is constant throughout the entire cell Ii . 3.1. Non-staggered cell division The method for non-staggered grid is divided into two steps: i. The constant states U ni1 ; U ni , U niþ1 define local Riemann Prob    lems RP U ni ; U niþ1 , RP U ni1 ; U ni the solution of which are 1

1

b nþ12 ðx; tÞ and U b nþ12 ðx; tÞ respectively. obtained as U iþ i 2

2

ii. Random sample these solutions at time t ¼ Dt nþ1 to find the solution for the cell Ii . The chosen point to randomly sample data depends upon the choice of hn which lies in the interval ½0; 1. The solution is thus obtained as:

U nþ1 i

8 n n 1 > < U i12 ðh Dx=DtÞ; if 0 6 h 6 2 n ¼ U iþ1 ððh  1ÞDx=DtÞ; if 1 < hn 6 1 : 2 > 2 :

ð8Þ

3.2. Staggered cell division For staggered grid, the solution for non-linear systems is obtained in two steps as: i. The constant states U ni1 ; U ni , U niþ1 define local Riemann Prob    lems abbreviated as RP U ni ; U niþ1 , RP U ni1 ; U ni the solution of 1

1

b nþ12 ðx; tÞ and U b nþ12 ðx; tÞ respectively. which are obtained as U iþ i 2

3.3. Solution of Riemann Problem The Riemann Problem for the shallow water equations abbreviated as RP is a problem involving two separate states says U ni and U niþ1 . At the point i þ 12, there appears to be a discontinuity because to its left lies the state U ni1 and the state U niþ1 to the right which are abbreviated as U L and U R respectively. The solution to   this RP U ni ; U niþ1 proceeds in the form of three characteristic waves (as the Jacobian Matrix of the flux vector F has three eigenvalues) originating from the point i þ 12 and traveling with speeds pffiffiffiffiffiffi u þ a; u; u  a respectively, with a ¼ gh (g is the acceleration due to gravity and h is height of water level). Out of these, the left and right waves can be either shock or rarefaction waves whereas the middle one is always a contact discontinuity. Thus, the entire solution plane seems to be divided into three separate regions as shown in Fig. 4. The region U L is to the left of the characteristic u-a, while U R lies to the right of the characteristic u + a. The region in between U L and U R is depicted by the star region with U L and U R the regions to the left and right of the contact discontinuity respectively. The solution to the given RP may lie in any of the above four regions depending upon the location of the point at which the solution is required. This solution to the Riemann Problem can either be calculated using the time consuming exact mathematical method or the same may be estimated by using approximate Riemann Solvers which are computationally more efficient and the solution provided gives better results with Godunov type methods.   In this case also, the given Riemann Problems, RP U ni ; U niþ1 ,     nþ1 nþ1 RP U ni1 ; U ni and RP U i12 ; U iþ12 are solved using the HLL approxi2

2

mate Riemann Solver which does not account for the contact wave present in the solution. However in the context of shallow water equations this does not account for much difference in solution and hence the HLL solver can be used.

2

We random sample these solutions in the range Dx at time 1

t ¼ Dtnþ2 , and obtain the solution as

  1 nþ1 b nþ12 hn Dx; Dtnþ12 ; U i12 ¼ U i 2

2

  1 nþ1 b nþ12 hn Dx; Dt nþ12 : U iþ12 ¼ U iþ 2

ð9Þ

2

 nþ1 nþ1  ii. Solve the Riemann Problem RP U i12 ; U iþ12 to find the solu2

2

b nþ1 ðx; tÞ and random sample it at time t ¼ Dtnþ1 to tion U i obtain the final solution as

b nþ1 ðhn Dx; Dt nþ1 Þ: U nþ1 ¼U i i

ð10Þ

Fig. 4. The Riemann Problem.

Table 1 Van der Corput sequence for k1 ¼ 2 and k2 ¼ 1.

Fig. 3. Sampling process for Random Choice Method.

n

m

a0

a1

1 2 3 4 5 6 7 8 9 10

0 1 1 2 2 2 2 3 3 3

1 0 1 0 1 0 1 0 1 0

1 1 0 0 1 1 0 0 1

a2

1 1 1 1 0 0 0

a3

1 1 1

A0

A1

0 1 0

1 1 0 0 1 1 0 0 1

A2

1 1 1 1 0 0 0

A3

hn

1 1 1

0.5000 0.2500 0.7500 0.1250 0.6250 0.3750 0.8750 0.0625 0.5625 0.3125

190

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

The HLL approximate Riemann Solver provides a solution   b nþ1 U ðx; tÞ to the Riemann Problem RP U ni ; U niþ1 according to followi1 2

ing equation

8 n > < U i1 b nþ1 U ðx; tÞ ¼ U hll i12 > : n Ui

if if if

9 > = x Si1 6 t 6 Si ; > ; x P Si t x t

6 Si1

ð11Þ

Si U ni Si1 U ni1 þF ni1 F ni , fU ni ; U ni1 g is the solution at time Si Si1 fF ni ; F ni1 g is the flux vector at t ¼ Dt n . Si1 represents the

where U hll ¼

t ¼ Dt n ; speed of discontinuity at the left side of Riemann Problem while Si represents the same at the right side of the Riemann Problem. For the shallow water equations, they are given by

 qffiffiffiffiffiffi qffiffiffiffiffiffiffiffiffiffiffi ^ ; ^  gh Si1 ¼ min ui1  ghi1 ; u  qffiffiffiffiffiffi qffiffiffiffiffiffiffi ^ ; ^  gh Si ¼ min ui  ghi ; u

Table 3 The initial conditions for two dimensional problem. X length Y length (breadth) X grid points Y grid points CFL factor Centre of Dam (X, Y) Radius of Circular Dam Water level inside Dam Water level outside Dam X component of initial velocity(Inside) X component of initial velocity(Outside) Y component of initial velocity(Inside) Y component of initial velocity(Outside) Van der Corput sequence

50 m 50 m 500 500 0.40 (25 m, 25 m) 10 m 10 m 3m 0 0 0 0 Binary(2, 1)

ð12Þ

where ui , ui1 represents the velocity and hi ; hi1 the height at U ni ; ^ are the Roe [23] averaged values of velo^; h U n respectively and u i1

city and height taken as

^ ¼ hi1 þ hi h 2

^¼ and u

pffiffiffiffiffiffiffiffiffi pffiffiffiffi ui1 hi1 þ ui hi pffiffiffiffiffiffiffiffiffi pffiffiffiffi : hi1 þ hi

ð13Þ

3.4. Time step and random number sequence For both staggered and non-staggered versions, the time steps 1

t ¼ Dtnþ2 and t ¼ Dt nþ1 need not be the same but must satisfy the stability criteria. Here hn represents a random number for the nth time step i.e. at a particular time step, t ¼ Dtn the sampling along all the cells Ii is done using the same random number hn : This choice for the same value of hn per time step is necessary because it is a well-established fact that more random the generation of sequence fhni g, the worse are the results of Random Choice Method unless the data set is similar or near similar. As a first step towards improvement, Chorin [7] has suggested the use of only one random number per time level rather than a new number per time interval and per cell. In a second step, Collela [8] advocates the use of pseudo random numbers of Van der Corput [8,23] type for best results.

Fig. 5. Comparison between FORCE method and Random Choice Method (t = 0.6 s).

3.5. Van der Corput numbers The numbers of Van der Corput sequence are well suited for use in Random Choice and related methods due to their low-discrepancy

Table 2 Initial conditions for one dimensional case. X length Y length (breadth) X grid points Y grid points CFL factor Left side water level height Right side water level height Separation distance X component of initial velocity(Left) X component of initial velocity(Right) Y component of initial velocity(Left) Y component of initial velocity(Right) Van der Corput sequence

50 m 50 m 500 500 0.40 10 m 3 m (Wet Bed), 0.001 m (Dry Bed) 25 m 0 0 0 0 Binary(2, 1)

Fig. 6. Comparison between the FORCE method and Random Choice Method (t = 1.8 s).

and uniform distribution over the range ½0; 1 with mean close to 0.5. For generating a sequence of Van der Corput numbers, two relatively prime or prime integers, ðk1 ; k2 Þ are chosen such that k1 > k2 > 0: The ðk1 ; k2 Þ Van der Corput sequence fhn g is then generated in accordance with the following.

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

Fig. 7. Comparison between the FORCE method and Random Choice Method (t = 0.6 s).

hn ¼

m X ðiþ1Þ Ai k 1 ; i¼0

9 > > > > > > =

Ai ¼ k2 ai ðmod ki Þ; > > m > X > i > > n¼ ai k1 : ;

191

Fig. 8. Comparison between the FORCE method and Random Choice Method (t = 1.8 s).

ð14Þ

i¼0

For binary choice of ðk1 ; k2 Þ i.e. k1 ¼ 2 and k2 ¼ 1 the sequence is depicted in Table 1. 4. The FORCE method The FORCE (First order Centered) Scheme is basically a re-interpretation of the original Random Choice Method on a staggered grid. It works out by replacing the stochastic versions of the solutions (Random Choice Method) by a deterministic version through integral averages of the Riemann Problems. b 1 ðx; tÞ the solution of a local We begin by denoting by U

Fig. 9. FORCE method (t = 0.5 s).

iþ2

b 1 ðx; tÞ the solution of a local Riemann Problem RPðU ni ; U niþ1 Þ; U i 2

Riemann Problem RPðU ni1 ; U ni Þ. 1

Then, setting Dt nþ2 ¼ 12 Dt; nþ1

U i12 ¼ 2

1 Dx

  1 b nþ12 x; Dt dx: U i2 2 12Dx

ð15Þ

  1 b nþ12 x; Dt dx: U iþ2 2 12Dx

ð16Þ

Z

1Dx 2

and nþ1

U iþ12 ¼ 2

1 Dx

Z

1Dx 2

Applying the integral form of conservation laws to (15) and (16) nþ1

U i12 ¼ 2

nþ1

U iþ12 ¼ 2

1 n Dt n ðU þ U ni Þ þ ðF  F ni Þ; 2 i1 2Dx i1

ð17Þ

1 n Dt n ðU þ U niþ1 Þ þ ðF  F niþ1 Þ: 2 i 2 Dx i

ð18Þ

b i ðx; tÞ the solution of a local Riemann Next, denoting by U  nþ1 nþ1  Problem RP U i12 ; U iþ12 ; the solution U nþ1 is obtained by taking 2

b i ðx; tÞ at local time t ¼ 1 Dt. the integral average of U 2

U nþ1 i

1 ¼ Dx

Z

1D x 2

12Dx

Fig. 10. Random Choice Method (t = 0.5 s).

2

b n ðx; DtÞdx: U i

The integral form of (19) using conservation laws is

ð19Þ

U ni ¼

  1  nþ12 Dt  nþ12 nþ1 nþ1 U 1 þ U iþ12 þ F 1  F iþ12 ; 2 i2 2Dx i2 2 2

ð20Þ

192

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

Fig. 11. Comparison between two FORCE method and Random Choice Method.

Fig. 13. Random Choice Method (t = 1.1 s).

Fig. 12. FORCE method (t = 1.1 s).

Fig. 14. Comparison between the FORCE method and Random Choice Method.

where

 nþ1  nþ1 F iþ12 ¼ F U iþ12  F RI iþ1 : 2

2

2

ð21Þ

F RI iþ1 is the Richtmyer [23]. Flux calculated as 2

F RI iþ12

  ¼ F U iþ1 ; 2

where U iþ1 and U i1 are given by Eqs. (18) and (17) respectively. 2

2

The Eq. (20) of FORCE Scheme can be written in conservative form

U nþ1 ¼ U ni þ i

 Dt  n F i1  F niþ1 ; 2 2 Dx

ð22Þ

with intercell numerical flux:

F force iþ1 ¼ 2



 1 nþ12 1  n F iþ1 þ F i þ F niþ1 2 2 2 1 Dx n ðU  U niþ1 Þ and F FORCE þ i12 4 Dt i

 1 nþ12 1 n 1 Dx  n ¼ F 1 þ ðF þ F ni Þ þ U  U ni : 2 i2 2 i1 4 Dt i1

Fig. 15. FORCE method (t = 0.5 s).

ð23Þ

The above Flux is the arithmetic mean of the Richtmyer [23] and Lax–Freidrichs [23] method:

F FORCE ¼

ðF RI þ F LF Þ : 2

ð24Þ

5. Geometry and boundary conditions For both one dimensional and two dimensional cases, the mesh size is set at 500  500, giving a total of 250,000 mesh elements (uniform, square shaped) with the x dimension domain length as 50 m, the y dimension domain length as 50 m and the spatial step sizes (Dx and Dy) as 0.1 in both directions. The boundary

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

Fig. 16. Hybrid method (t = 0.5 s).

193

Fig. 19. Hybrid method (t = 1.1 s).

Fig. 20. Comparison between the FORCE method and Hybrid method.

Fig. 17. Comparison between the FORCE method and Hybrid method.

Fig. 18. FORCE method (t = 1.1 s).

conditions are set in such a way that the height of water level, x component of velocity and the y component of velocity at the boundaries take the value same as the cells to their left, right, top or bottom, accordingly as the implied boundary is the right, left, bottom or top one.

Fig. 21. Comparison between analytical solution and FORCE method (t = 1.2 s).

5.1. One dimensional case This problem features an upstream reservoir filled with water at height 10 m separated from the downstream which in turn has the

194

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

Fig. 25. The contour plot for the Partial Dam-Break Flow at time t = 1.2 s (FORCE method). Fig. 22. Comparison between analytical solution and RCM method (t = 1.2 s).

5.2. Two dimensional case The two dimensional case is of the Circular Dam Break Problem in which a circular tank of radius 10 m is located in the center of the 50  50 m area and is filled with water at height 10 m. The surroundings are filled with water at height 3 m. At t ¼ 0 s, the Dam is broken and the flow is simulated. The initial conditions and values are given in Table 3. 6. Results 6.1. One dimensional simulation (Wet Bed)

Fig. 23. Comparison between analytical solution and FORCE method (t = 1.2 s).

The computed results are presented through graphs in Figs. 5 and 6. From the above figures it may be noted that the FORCE method is dissipative as compared to Random Choice Method, as the solution obtained by FORCE method has smearing over many cells near the shock discontinuities and the desired sharpness is not obtained. Random Choice Method on the other hand gives better resolution of the shock characteristics by giving the desired sharpness to the graph owing to the random approach of the method. 6.2. One dimensional simulation (Dry Bed) The computed results are presented through graphs in Figs. 7 and 8. It is observed that in case of dry bed too, the solution obtained by Random Choice Method is much more sharper as compared to the one obtained by FORCE method. 6.3. Two dimensional Dam Break Problem The graphs presented in Figs. 9–14 clearly depict that the Glimm’s Random Choice Method generates more significant errors as well as inaccurate results and non smooth solutions everywhere when applied to the two dimensional case. The FORCE method on the other hand gives less significant errors but smearing of values near shock discontinuities leads to a smoother curve. According to Collela [8] this is largely due to the fact that Random Choice Method gives incorrect results where significantly large jumps occur in pressure and velocity values across oblique to the mesh directions, which is the case for the Circular Dam-Break Problem.

Fig. 24. Comparison between analytical solution and RCM method (t = 1.2 s).

7. Correction of results (Hybrid method) water level of 3 m. At t ¼ 0 s, the Dam is broken and the flow is simulated. The initial conditions and values are given in Table 2.

In order to correct the obnoxious behavior of the Random Choice Method for the two dimensional case, we employ a method

195

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

Firstly, we estimate the maximum and minimum height (pressure) across three cell interfaces; before, after and on the current cell as shown below max

¼ maxðhi1 ; hi ; hiþ1 Þ

min

¼ minðhi1 ; hi ; hiþ1 Þ

hi hi

n

n

n

n

n

n

ð25Þ

:

The respective method (Random Choice or FORCE) is then chosen as

ðU nþ1 Þhybrid i

Fig. 26. The contour plot for the Partial Dam-Break Flow at time t = 1.2 s (RCM).

¼

8 < ðU nþ1 Þ i

:

FORCE

if

ðU nþ1 ÞRCM i

hmax i hmin i

9 > C0 =

otherwise

;

:

ð26Þ

The choice of C 0 depends on several factors and it depicts the strength of the smallest shock observed across the cell interface, which will be treated as a discontinuity and hence applied with the FORCE method. The choice of C 0 needs to be optimum as incorrect values may lead to unnecessary smoothing of shocks or erroneous results. 7.1. Two dimensional Dam Break Problem revisited using the Hybrid method The results obtained by using the Hybrid method are presented in Figs. 15–20. It may be noted here that the extreme haphazardness visible in Random Choice Method (Figs. 11–16) is moderated to a great extent by using the Hybrid method. However, the sharpness across shocks which is a characteristic feature of Random Choice Method is removed and the solutions now seem to have increased smoothness. This is accounted due to the use of FORCE method across large pressure jumps in Eq. (26). 8. Validation of results 8.1. One dimensional problem

Fig. 27. The contour plot for the Partial Dam-Break Flow at time t = 1.2 s (Hybrid method).

in which we apply a Godunov type method(FORCE in current case) across the cells which show large changes in the pressure values for both x and y sweeps of dimensional splitting. For the remaining cells we proceed with the normal Random Choice Method as applied previously. There is no fixed method for identifying the cells on which the FORCE method shall be applied and on which RCM is applied. We have used the following criteria derived from the original wok of Collela [8] and that of Miniati [15].

Table 4 The initial conditions and values for Partial Dam Break Flow problem. X length Y length (breadth) X grid points Y grid points CFL factor Left side water level height Right side water level height Breach width Separation distance X component of initial velocity(Left) X component of initial velocity(Right) Y component of initial velocity(Left) Y component of initial velocity(Right) Van der Corput sequence

50 m 50 m 500 500 0.40 10 m 3 m(Wet Bed) 18.7 m (from Centre of Dam Wall) 25 m 0 0 0 0 Binary(2, 1)

To provide a validation to our methods, we compare our numerical solutions against the analytical solution provided by Stoker [18] for both the wet bed and dry bed. It is evident from Figs. 21–24 that the numerical solutions obtained here are in close agreement with the analytical solutions provided by Stoker [18]. 8.2. Two dimensional problem The code used for generating solution of two dimensional problems is applied on another problem in the same domain, the Partial Dam Break Flow problem, to provide a wider view of its capacities. In this problem, the one dimensional Dam is broken partially i.e. only a certain fixed part of the separating wall is broken at time t = 0 and the flow is simulated. The computed results for FORCE method, RCM and Hybrid method is presented in Figs. 25–27 which is found to be in close agreement with the results presented in Baghlani [3]. The initial conditions and values are given in Table 4. 9. Conclusion The present paper investigates the application of Glimm’s Random Choice Method to the one dimensional homogenous shallow water equations of the Dam-Break Problem. Further the method is extended to two-dimensional case using operator splitting procedure. We observe that the Random Choice Method gives very sharp shock capture and accurate results for the one dimensional case. However for two dimensional problem, the approach of splitting in x and y direction is not a good choice with Random Choice Method probably due to the jump in pressure (height) values

196

H. Gupta, L.P. Singh / Computers & Fluids 111 (2015) 187–196

across oblique cell interfaces. The FORCE method on the other hand is dissipative and provides unwanted smoothness in the solutions, even at the shock discontinuities for both one dimensional and two dimensional cases. The Hybrid method attempts to combine the properties of FORCE method (smoothness) with the original Random Choice Method (sharpness) to obtain better solutions. This in turn is heavily dependent on several criteria and one has to compensate for added smoothness by reduction in the sharpness level of the solution at shock discontinuities. An incorrect choice of the same may lead to extreme smoothness (dissipative) or erroneous sharp values (close to original Random Choice Method). References [1] Ahmada MF, Mamatb M, Wan Nika WB, Kartono A. Numerical method for dam break problem by using Godunov approach. Appl Math Comp Intell 2013;2(1). [2] Amiri SM, Talebbeydokhti N, Baghlani A. A two-dimensional well-balanced numerical model for shallow water equations. Sci Iran 2013;20(1):97–107. [3] Baghlani A. Simulation of dam-break problem by a robust flux-vector splitting. Sci Iran 2011;18(5):1061–8. [4] Einfeldt Bernd. On Godunov type method for gas dynamics. SIAM J Numer Anal 1988;25(2):298–318. [5] Chaudhary M Hanif. Open-channel flow, vol. 2. Colombia: Springer; 2008. [6] Chippada S, Dawson CN, Martinez ML, Wheeler MF. A Godunov type finite volume method for the system of shallow water equations. Comput Methods Appl Mech Engrg 1998;151(1–2):105–29. [7] Chorin Alexendre Joel. Random choice solution of hyperbolic systems. J Comput Phys 1976;22:517–27. [8] Collela Phillip. Glimm’s method for gas dynamics. SIAM J Sci Stat Comput 1982;3(1):76–110. [9] Concus Paul, Proskurowski Wlodzimierz. Numerical solution of nonlinear hyperbolic equations by the random choice method. J Comput Phys 1979;30(2):153–8.

[10] Pan Cun-hong, Dai Shi-qiang, Chen Sen-mei. Numerical simulation for 2d shallow water equations by using Godunov-type scheme with unstructured mesh. J Hydrodyn 2006;18(4):1–2. [11] Gustafson John L. Reevaluating Amdahl’s law. Communications of the ACM 0001–0782/88/0500-0532, 1988. [12] Han EE, Warnecke Gerald. The exact Riemann solutions to shallow water equations. Preprint, 2011–027. [13] Lencina Ignacio Valcárcel. Comparison between 1D and 2D models to analyze the dam break wave using the FEM and shallow water equations. TRITA-LWR Master Thesis, 2007, (ISSN 1651–064X, LWR-EX-07-17). [14] Marshall Guillermo, Mendez Rau L. Computational aspects of the random choice method. J Comput Phys 1981;39:1–21. [15] Miniati Francesco. Glimm–Godunov’s method for cosmic ray hydrodynamics. J Comput Phys 2007;227:776–96. [16] Silva Luis Moura E, Buyyaz Rajkumar. Parallel programming models and paradigms. In: Buyya R, editor. High performance cluster computing: architectures and systems, vol. 2. NJ, USA: Prentice Hall PTR; 1999. [17] Skoula ZD, Borthwick AGL, Moutzouris CI. Godunov-type solution of the shallow water equations on adaptive unstructured triangular grids. Int J Comput Fluid Dyn 2006;20(9):621–36. [18] Stoker JJ. Water waves mathematical theory with applications, vol. IV. New York: Interscience Publishers Inc.; 1957. [19] Sun Xian-He, Chen Yong. Reevaluating Amdahl’s law in the multicore era. J Parallel Distrib Comput 2010;70:183–5. [20] Toro EF. The dry-bed problem in shallow-water flows. College of Aeronautics Reports, 1990. [21] Toro EF. Riemann problems and the WAF method for solving the twodimensional shallow water. Philos Trans: Phys Sci Eng 1992;338(1649):43–68. [22] Toro EF. Shock-capturing methods for free-surface shallow flows. Wiley and Sons. Ltd.; 2001. [23] Toro Eleuterio F. Reimann solvers and numerical methods for fluid dynamics. 3rd ed. Springer; 2009. [24] Zhou JG, Causon DM, Ingram DM, Mingham CG. Numerical solutions of shallow water equations with discontinuous bed topography. Int J Numer Methods Fluids 2001:1–2.