Đề tài Hydraulic modeling of open channel flows over an arbitrary 3-D surface and its applications in amenity hydraulic engineering

Hydr. Eng., ASCE, Vol. 121(12), pp 900–905. 12. Yıldırım N. and Kocabas F. 1997. Closure to ‘Critical submergence for intakes in open channel flow. J. Hydr. Eng., ASCE, Vol. 123(6), pp 589–590. 13. Yıldırım N. and Kocabas F. 1998. Critical submergence for intakes in still-water reservoir. J. Hydr. Eng., ASCE, Vol. 124(1), pp. 103–104. 56 Chapter 5 UNSTEADY PLANE-2D ANALYSIS OF FLOWS WITH AIR-CORE VORTEX In most of the studies in the past, because of the complexity of 2D equation the variation of the flow quantities in the tangential direction at the same radius from the center of the vortex was neglected with the assumption of homogeneity of flow in circular direction. In this study, with the advantages of the 2D depth-averaged model, it is possible to estimate the variation in 2D flow. 5.1 Governing equations Coordinate setting: In order to apply the depth-averaged equations (3.26, 3.33 & 3.34) to simulate the 2D flows with air-core vortex, a new coordinate system is introduced as in Figure 5.1. Two coordinates ),( ηξ are set on the bottom and the other ( )ζ perpendicular with that surface. Parameters ( )0 1, , ,a b R R define the geometry of the bottom in which ( )1R is 57 radius of the intake, ( )0R refers to the curvature of the edge of the intake entrance, while a and b are radii of the flow approach area and the length of intake, respectively. Figure 5.1 Definition sketch of the new coordinates Governing equations: The original depth-averaged equations are used: Continuity equation: 01 000 = ∂ ∂ + ∂ ∂ + ∂ ∂ J N J M t h J ηξ (5.1) Momentum equations: ξ -direction ( )ξηηξηξξξηξξξηξ 0200020000 1 Γ+Γ+Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NNMMNM hJJ VM J UM J M t 000 000000 00 2 0 2 0 2 0 0 J dp J dp J G J h b h zzyyxx h zyx ρ τζ ρη ηξηξηξζ ρξ ξξξ ξξ − ∂ ∂++ − ∂ ∂++ −= ∫∫ (5.2) b 0R 1R a x y z η ξ ζ 58 η -direction ( )ηηηηηξηξηηξξηξ 0200020000 1 Γ+Γ+Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NNMMNM hJJ VN J UN J N t 000 000000 00 2 0 2 0 2 0 0 J dp J dp J G J h b h zzyyxx h zyx ρ τζ ρξ ηξηξηξζ ρη ηηη ηη − ∂ ∂++ − ∂ ∂++ −= ∫∫ (5.3) In this case, it can easily be seen that the coordinates setting on the bottom surface with two axis ),( ηξ are orthogonal system, therefore it was proved that: 0=Γ=Γ ξηξ ξ ξη 0=Γ=Γ ηηη η ξξ 0=Γ=Γ ζηξ ζ ξη 0000000 =++ zzyyxx ηξηξηξ and 0=ηG As a result, Equations (5.2-5.3) are reduced to ( )ξηηξξξηξ 02020000 1 Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NM hJJ VM J UM J M t 0 222 0 2 0 2 0 2 0 0 222 J hGNM J G J h bzyx ρ τ ξ ξξξ ξζζ µη ζ ξξ ξ −⎥⎦ ⎤⎢⎣ ⎡ −Γ+Γ ∂ ∂++ −= (5.4) ( )ηηξηξηηξ 000000 1 Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NMMN hJJ VN J UN J N t 0 222 0 2 0 2 0 2 0 222 J hGNM J bzyx ρ τ η ηηη ηζζ µη ζ ξξ −⎥⎦ ⎤⎢⎣ ⎡ −Γ+Γ ∂ ∂++ −= (5.5) These equations together with continuity Equation (5.1) are then applied to calculate 2D 59 water surface profile of flow with air-core vortex at a vertical intake in the following sections. j=2 j=3 i=2 j=1 i=1 Figure 5.2 Illustration of the computational grid 5.2 Numerical method Based on the generalized curvilinear coordinate system set on the bottom surface (Figure 5.1), the computation mesh is generated as shown in Figure 5.2. Then the set of equations (5.1, 5.4-5.5) is discretized by using the cell-centered staggered grid (Prayet and Taylor, 1983), in which the hydraulic variables are defined at different positions as shown in Figure 5.3 (Fletcher, 1991). The Finite Difference Method and the Finite Volume Method are used in dicretization process, in which the First-Upwind scheme is employed for convective terms while the Two-point Forward scheme is utilized for time derivatives. It follows that: Continuity equation: ⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + ∆ − + + ++ ++++ + ++ ++ 21,0 21, 21,10 21,121,21 1 21,21 21,210 11 ji ji ji ji n ji n ji ji J M J M t hh J ξ 60 Figure 5.3 Definition sketch of cell-centered staggered grid in the 2D calculation 01 ,210 ,21 1,210 1,21 =⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + + + ++ ++ ji ji ji ji J N J N η (5.6) Momentum equations: ξ -direction ⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + ∆ − +− ++−+− ++ +++++ + + + 21,210 21,121,21 21,210 21,21,2121, 1 21, 21,0 11 ji jbiji ji jaiji n ji n ji ji J MU J MU t MM J ξ [ ] =Γ+Γ+⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + + + − + ++ 21.0 2 0 2 21.0,0 ;21,, 1,0 21,1, 11 ji jiji djiji ji cjiji hVhU JJ MV J MV ξ ηη ξ ξξη [ −Γ++ ∆ −−= ++++ ++ + + + 21.21 2 21.21 21.0 2 0 2 0 2 0 21. 21. 21.0 21. 2 1 jiji ji zyx ji b ji ji ji M J G J h ζ ξξ ξ ξ ξξξ ξρ τ 21.21 2 21.2121.21 2 21.2121.21 2 21.21 +−+−+++++−+− Γ−Γ+Γ jijijijijiji NNM ζηη ζ ηη ζ ξξ ]2 21.2121.212 21.2121.21 +−+−++++ +− jijijiji hGhG ζζ (5.7) M, U M, U N, V N, V h i+1i j j+1 61 Figure 5.4 Illustration of the discretization scheme for the momentum equations N, V h i j-1/2 j+1/ h M, U i+1 j M, U M, U M, U M, U h i-1/2 j j+1 h N, V N, V N, V N, V i+1/i 62 η -direction ⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + ∆ − +− + ++++ + + + ji jbiji ji jaiji n ji n ji ji J NU J NU t NN J ,0 ,21, ,10 ,21,1,21 1 ,21 ,210 11 ξ [ ] =Γ+Γ+⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + + −+ +−+−+ ++ ++++ ji ji djiji ji cjiji J UVh J NV J NV ,21 021,210 1,2121,21 21,210 ,2121,211 η ηξ η ξηη [ 21,21 2 21,2121,21 2 21,21 ,210 2 0 2 0 2 0 ,21 2 1 −+−+++++ ++ Γ−Γ ++ ∆ −−= jijijiji ji zyx ji b MM J ζ ξξ ζ ξξ η ηηη ηρ τ 21,21 2 21,2121,21 2 21,2121,21 2 21,21 ++++−+−+++++ −Γ−Γ+ jijijijijiji hGNN ηη ζζ ηη ζ ηη ] 21,21 2 21,21 −+−+ + jiji hG ηη ζ (5.8) To close the above equation system, the boundary conditions were given with the discharge at the upstream end and Neumann condition at the downstream boundary. The computational domain is chosen from 1=j to 1+= JYj , therefore the hydraulic variables at 1+= JYj is used for boundary condition when calculating variables at 1=j and vice versa. 5.3 Results and discussions With different discharges ( )0M and tangential velocities ( )0V given at the outer-zone, some results of the surface profile were obtained as shown in Figures (5.5-5.7). Similarly as the phenomenon monitored in Chapter 4, the water surface and submergence of the flow with air-core vortex strongly depend on the flow condition such as discharge at the intake, the circulation (which in this simulation is defined from the velocity at boundary of outer zone) and also the geometry of the intake. When the discharge at the intake is increased, the submergence is also increased as shown in Figure 5.5 with the same circulation. If the discharge at the intake is maintained, when 63 0.0 0.1 0.2 0.3 0.4 0.5 z ( m ) -0.3 -0.1 0.0 0.2 0.3 x (m) -0.30 -0.15 0.00 0.15 30 y (m) 0.0 0.1 0.2 0.3 0.4 0.5 z ( m ) -0.3 -0.1 0.0 0.2 0.3 x (m) -0.30 -0.15 0.00 0.15 30 y (m) Figure 5.5 Water surface of flow with different discharges at the intake a) 0 0 0 10.35; 10.; 0.03 ; 0.05 ;QM V R m R m= = = = b) 0 0 0 10.5; 10.; 0.03 ; 0.05 ;QM V R m R m= = = = 0.0 0.1 0.2 0.3 0.4 0.5 z ( m ) -0.3 -0.1 0.0 0.2 0.3 x (m) -0.30 -0.15 0.00 0.15 30 y (m) 0.0 0.1 0.2 0.3 0.4 0.5 z ( m ) -0.3 -0.1 0.0 0.2 0.3 x (m) -0.30 -0.15 0.00 0.15 30 y (m) Figure 5.6 Water surface of flow with different velocity at the outer-zone boundary a) 0 0 0 110.; 0.5; 0.05 ; 0.05 ;V QM R m R m= = = = b) 0 0 0 115.; 0.5; 0.05 ; 0.05 ;V QM R m R m= = = = a) b) a) b) 64 velocity at boundary varies the water surface also varies. The higher circular velocity ( )0V , the higher submergence at the intake, or in other words, the circular velocity obstruct or impede the water drainage at the intake. Consequently to maintain the same discharge, water head would be increased with increased circular velocity (Figure 5.6). In order to investigate the effect of the intake geometry to water surface, the parameter 0R is used to control the shape of the join between intake and bottom. Figure (5.7) shows that the submergence increases when 0R increases while keeping other parameters constant. 0.0 0.1 0.2 0.3 0.4 0.5 z ( m ) -0.3 -0.1 0.0 0.2 0.3 x (m) -0.30 -0.15 0.00 0.15 30 y (m) 0.0 0.1 0.2 0.3 0.4 0.5 z ( m ) -0.3 -0.1 0.0 0.2 0.3 x (m) -0.30 -0.15 0.00 0.15 30 y (m) Figure 5.7 Water surface of flow with different shape of the intake a) 0 0 0 110.; 0.5; 0.05 ; 0.05 ;V QM R m R m= = = = b) 0 0 0 110.; 0.5; 0.03 ; 0.05 ;V QM R m R m= = = = 5.4 Summary Based on the simulation results, it can be concluded that, not only 1D water profile as described in Chapter 4 can be estimated by using analytical method, but also the 2D unsteady plane water surface can be reasonably simulated using the depth-averaged equations without any difficulty. Consequently, it also validates the broad perspective of the model’s applicability. a) b) 65 5.5 References 2 Peyret R. and Taylor T. D. 1983. Computational methods for fluid flow, Springer Ser. Comput. Phys. (Springer, Berlin, Heidelberg). 3 Fletcher C. A. J. 1991. Computational techniques for fluid dynamics. Springer Ser. Comput. Phys. (Springer Verlag). 66 Chapter 6 WATER SURFACE PROFILE ANALYSIS OF FLOWS OVER CIRCULAR SURFACE 6.1 Preliminary In a chute spillway, the stability of the free surface is strongly dependent on the geometry of the approach channel as well as that of spillway. If there is an abrupt drop of the channel bed, a free overfall or a vertical fully aerated drop may occur. This phenomenon has attracted considerable attentions of investigators (e.g., Davis et al. 1998; Dey 1998, 2002, 2003, 2005; Ahmad 2005; Guo 2005; Pal and Goel 2006) since the pioneering work of Rouse (1936) due to its practical engineering use in simple flow-measuring devices (e.g., Dey 1998; Guo 2005). But their approaches could not be applied in case of an abrupt change of channel’s bed with circular shape, accompanied with small discharge that cannot result in free overfall. In this chapter a mathematical model based on depth-averaged equations with a 67 generalized curvilinear coordinate system attached to the bottom surface is developed from the 2D depth-averaged model described in Chapter 3. The developed model was applied to calculate the water surface profile in open channel with circular surface at the end. To verify the model, an experiment was conducted to measure the water surface profile, where the model results of both steady and unsteady analysis were compared with the observations. 6.2 Hydraulic Experiment 6.2.1 Experimental Setup The experiments were conducted in an open-channel main flume consisted of two straight channels joined perpendicular to each other by a step ‘A’ at the end as described in Figure 6.1. The flume was 135cm long, 10cm wide and 20cm deep (Figure 6.1a) , made of a steel frame with glass in all side walls and bed which allowed the convenience for visual observations. In this study two detachable circular parts with radius of 15cm and 5cm respectively were used (Figure 6.1b). These parts were made such that, it is easily to be attached and detached to and from the main flume at step ‘A’, without altering the other parts of experiment facility to form a circular channel bottom at the joint (Figure 6.1c). These parts (B and C) were purposely constructed to investigate the effect of different curvatures to flows. For simplicity of discussion, from now onward, the curvature with 15cm and 5cm radius will be referred to “large case” and “small case” respectively. The water surface elevations were measured by the ultrasonic sensor Keyence UD-500 accompanied with the receiver and converter Keyence NR-2000 (Figure 6.2). For each case the measurements were focused in the region with the curvature of the bottom at four points 0o (point 1), 30o (point 2), 60o (point 3) and 90o (point 4) to the vertical direction (Figure 6.1d). 68 Figure 6.2 Experimental site Figure 6.1 Side view of the experimental facility (All units are in centimeter) b) Part B and part C R = 15 R = 5 c) Connection of main flume and parts B, C Point 1 Point 2 Point 3 Point 4 d) Sensors’ arrangement 135 100 20 20 20 75 20 20 a) Main flume Inflow Outflow A Sensor UD-500 Circular surface Digitizer NR-2000 69 The ultrasonic sensor Keyence UD-500 measures the distance from the sensor to the surface by calculating the reflected ultrasonic signal. Actually, the sensor is combined both the transmitter and receiver of the ultrasonic signals, and the output is voltage of the direct electric current, therefore it needs a receiver and digitizer to convert the analog signal into digital signal (Keyence NR-2000) for out-putting to the PC. Thus, the sensor is connected with the receiver and converter as in Figure 6.3. 6.2.2 Experimental procedures and results: Before using the sensor and converter, it is required to calibrate the sensor, and establish a correlation line between the voltage and the distance from the sensor to the surface. The calibration process was done by using the sensor to measure the known distances (Figure 6.4) and the correlation equation was found: tt VD *77.4035.52 −= (6.1) where: tD is the distance (in mm ) and tV is the measured voltage (in Voltage V ) at time t . At the first stage of experiment, the distance from the fixed sensor to the bed was measured in advance ( )bD and then the time series of distance from the sensor to the water surface were saved ( )sD , thereafter water depth was easily obtained by subtracting sD from bD . Water was drawn from a constant head tank being filled by a turbine pump to ensure the steady in-flow rate. The flow rate was controlled by a control valve and discharge measurements were carried out at the end of the flume for different flow conditions. At beginning, the discharge was set at very small value, then increased with small intervals. After each increase, the flows were kept for some time to get the equilibrium condition before taking measurements. Upon attainment of steady conditions at the 70 upstream end of the main flume, the water surface level in the region within circular bottom was measured by the sensor at four locations as in Figure 6.1d. The samples of time series results of water depth are shown in Figure 6.5. And the time-averaged values with the flow conditions are listed in Table 1. At the first stage, with relatively small discharge, the flows over the circular part were observed with very slightly fluctuation and considered to be stable (Figure 6.5a &6.5b). In the second stage, when the discharge was continuously increased, the oscillation occurs at the circular region and become significant with increasing flow rate while remaining other parameters as the same. With each subsequent increase of the discharge, the amplitudes of oscillation become larger (Figure 6.5c & 6.5d). The intensities of oscillation can be estimated by Equation (6.2) and the results are shown in Figure 6.6. ∑ ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − = n i h hh n 1 2 1 σ (6.2) From this figure, it is observed that, in both small- and large-case, the bigger intensities were observed in downstream end of the flow field, or in other words, the oscillation intensity increases from point 1 to point 4. 71 Figure 6.3 Schematic of sensor connection Figure 6.4 Sensor calibration Sensor UD-500 Digitizer NR - 2000 Transmitter & Receiver PC 72 0 1 2 3 4 5 Time (s) 0 5 10 15 D ep th (m m ) Point 1 Point 3 Point 2 Point 4 a) 0 1 2 3 4 5 Time (s) 0 5 10 15 D ep th (m m ) b) 0 1 2 3 4 5 Time (s) 0 5 10 15 D ep th (m m ) c) 0 1 2 3 4 5 Time (s) 0 5 10 15 D ep th (m m ) d) Figure 6.5 Time history of the free surface at four locations in different experiments: a) Exp1; b) Exp2; c) Exp3; d) Exp4; 73 Table 6.1 Experiment conditions Time averaged water depth (mm) Experiment Radius of circular part Discharge (l/s) Point1 Point2 Point3 Point4 Small-case Exp 1 5 cm 0.1704 6.2 3.8 2.2 1.4 Small-case Exp 2 5 cm 0.2441 8.6 4.8 3.1 1.7 Small-case Exp 3 5 cm 0.6350 11.1 7.8 5.0 3.5 Small-case Exp 4 5 cm 0.7754 14.4 10.5 7.8 5.6 Large-case Exp 5 15 cm 0.4330 10.8 5.1 2.2 1.8 Large-case Exp 6 15 cm 0.5421 12.7 5.7 3.9 3.1 Large-case Exp 7 15 cm 0.8295 15.8 8.6 5.1 3.8 Large-case Exp 8 15 cm 0.9936 18.6 10.8 7.3 5.1 0.00 0.05 0.10 0.15 0.20 0.25 0.30 Point1 Point2 Point3 Point4 σ Exp-1 Exp-2 Exp-3 Exp-4 Exp-5 Exp-6 Exp-7 Exp-8 Figure 6.6 The oscillation density at four locations in circular region 74 6.3 Steady analysis of water surface profile Firstly, the new generalized curvilinear coordinate is introduced. The coordinate system is generated based on the bottom surface with two axes ( )ηξ − attached to the surface and the other axis ( )ζ is a straight line perpendicular to it (Figure 6.7). In Figure 6.7, η is the straight axis normal to plane ( )xOz or is identical with y -axis. The surface is expressed mathematically by the equation 0),,( =zyxζ . The depth-averaged equations in this new coordinate system are used (Equations 3.26, 3.33 and 3.34): Continuity equation: 01 000 = ∂ ∂ + ∂ ∂ + ∂ ∂ J N J M t h J ηξ (6.3) Momentum equation in ξ -direction: ( )ξηηξηξξξηξξξηξ 0200020000 1 Γ+Γ+Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NNMMNM hJJ VM J UM J M t ∫∫ ∂∂ ++ −− ∂ ∂++ −= h zzyyxxb h zyx dp JJ dp J G J h 00 000000 000 2 0 2 0 2 0 0 ζ ρη ηξηξηξ ρ τζ ρξ ξξξ ξξ (6.4) Momentum equation in η -direction: ( )ηηηηηξηξηηξξηξ 0200020000 1 Γ+Γ+Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NNMMNM hJJ VN J UN J N t ∫∫ ∂∂ ++ −− ∂ ∂++ −= h zyxb h zzyyxx dp JJ dp J G J h 00 2 0 2 0 2 0 000 000000 0 ζ ρη ηηη ρ τζ ρξ ηξηξηξ ηη (6.5) 75 Figure 6.7 Curvilinear coordinates attached to the bottom ζ ξ W U ζe r ξe rηe r COVARIANT BASE VECTORS x z O 76 Neglecting transversal transport, then it follows that: 01 00 = ∂ ∂ + ∂ ∂ J M t h J ξ ç ç (6.6) 0 2 0 2 0 2 0 2 0 0 0 2 000 22 1 J hGM J G J hM hJJ UM J M t bzx ρ τ ξ ξξ ξ ξζζ ξξ ξξ ξξ −⎥⎦ ⎤⎢⎣ ⎡ −Γ ∂ ∂+ −=Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ (6.7) Using the system of equations (6.6-6.7), we can analyze the water surface profile of the flow over a circular surface as depicted in Figure 6.7. Consider the steady state, equations (6.6-6.7) can be recast as follows: Continuity equation: 000 00 const. J M 0 JQMQ J M d d =⇔==⇔=ξ (6.8) or h JQ UJQUh 0000 =⇔= (6.9) Momentum equation: ξξ ξξξ GJ hM hJJ UM d d 0 0 2 00 1 =Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ 0 2 0 2 0 2 0 2 0 22 J hGM d d J bzx ρ τ ξ ξξ ξζζ ξξ −⎥⎦ ⎤⎢⎣ ⎡ −Γ+− (6.10) Substitute Equations (6.8), (6.9) into (6.10), after some manipulations, the equation for water surface profile can be obtained in the form of Equation (6.11): ),( ),( 2 1 ξ ξ ξ hf hf d dh −= (6.11) with 20 0 0 00 2 0 2 0 2 034 2 0 2 0 1 222 ),( h d dJ d d JJQhGh d dGhf zxzx ⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ Γ+ Γ+ −+ + = ξξ ξξ ξ ξξξ ζξξ ζ ξξξζ 2 0 2 000 0 0 2 0 JfQhJd dJ JQ −⎥⎦ ⎤⎢⎣ ⎡ Γ+− ξξξξ (6.12) 77 and ( ) 2020320202 ),( JQhGhf zx ++= ζξξξ (6.13) The common method of analysis including singular point analysis similar as in Chapter 4 is applied to calculate the water surface profile in both the upstream and downstream direction from this point. The singular point is defined as the point at which both functions ),(1 ξhf and ),(2 ξhf in Equation (6.11) are equal zero. The definition of singular and an example of critical and quasi-normal depth lines can be found in Figure 6.8. Figures 6.9 to 6.12 show the comparison of water surface profile between the 1D model results under steady state condition and the experiment for Exp-01, Exp-02, Exp-05 and Exp-06. It is observed from the figures that the calculated results are consistent with the experimental data in both large and small cases. Note that during the experiment, because of limitation of sensors the measurements were only focused in the vicinity of the circular surface. 78 0.95 1.00 1.05 1.10 x (m) 0.50 0.55 z (m ) Bottom Critical line Quasi-normal depth line Water surface Figure 6.8 Illustration of computed water surface profile with quasi-normal and critical depth lines 79 0.90 0.95 1.00 1.05 1.10 x (m) 0.45 0.50 0.55 0.60 z (m ) Cal. Exp. Bottom Figure 6.9 Steady water surface profile with conditions of Exp1 0.90 0.95 1.00 1.05 1.10 x (m) 0.45 0.50 0.55 0.60 z (m ) Cal. Exp. Bottom Figure 6.10 Steady water surface profile with conditions of Exp2 80 0.90 1.00 1.10 1.20 x (m) 0.50 0.60 0.70 z (m ) Cal. Exp. Bottom Figure 6.11 Steady water surface profile with conditions of Exp5 0.90 1.00 1.10 1.20 x (m) 0.50 0.60 0.70 z (m ) Cal. Exp. Bottom Figure 6.12 Steady water surface profile with conditions of Exp6 81 6.4 Unsteady characteristics of the flows In order to study the unsteady characteristics of flows, the numerical method for unsteady analysis is developed by applying the common method used in the field of the numerical hydraulics. Figure 6.13 Illustration of staggered grid Equations (6.6-6.7) and the same curvilinear coordinate with previous steady analysis is also used, and for solving those equations numerically, the two-point forward scheme is applied for time integration, on the other hand the First Upwind scheme is employed for the convective term based on cell-centered staggered grid (Figure 6.13) for spatial discretization as in the following Equations (6.14-6.15). Continuity equation: 011 01 0 121 1 21 21 0 =⎥⎥⎦ ⎤ ⎢⎢⎣ ⎡ − ∆ + ∆ − + ++ + + + i n i i n i n i n i i J M J M t hh J ξ (6.14) Momentum equation: ⎥⎦ ⎤⎢⎣ ⎡ − ∆ + ∆ − +−−++ + 1-i 0 121 i 0 21i 0 1 J MU J MUJ t MM biiaiini n i ξ ξξ ξξ iiiii GhhU 0 2 =Γ+ ( ) [ ] 10 2 10 2 2 0 2 0 2 −− Γ−Γ ∆ + − iiii izx MM ζξξ ζ ξξξ ξξ ( ) [ ] i b iiii izx hGhG ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ −− ∆ + + −−+ ρ τ ξ ξξ ξζζ 2 211 2 21 2 0 2 0 2 (6.15) i-1 i-1/2 i i+1/2 i+1 : UM , : h 82 where: ⎪⎩ ⎪⎨ ⎧ < ≥ = ⎪⎩ ⎪⎨ ⎧ < ≥ = − − + + 0 if 1 0 if 0 ; 0 if 1 0 if 0 21 21 21 21 i i i i U U b U U a and 2 1 21 n i n in i UU U + = + + The bottom shear stress was evaluated using the following equation: 2UefUUfb ξ ξ ρ τ rr == (6.16) where f is the resistant coefficient. In this study, since all walls were made of smooth glass, the value of f = 0.01 was used. ξe r is covariant base vector as shown in Figure 6.7. The model has been applied in estimating the water surface profile of flows over a circular surface (Figure 6.7), based on the flow conditions identical with that of the experiments (Table 6.1). The inflow discharge is given at the upstream and the Neumann derivatives are provided at the downstream end to complete the boundary conditions. In order to determine the unsteady characteristics of the flow over circular surface, the same procedures as in experiment were applied for simulation. Firstly the small constant discharge was used as upstream boundary (Exp -1, Exp -2, Exp – 5 and Exp - 6). With constant initial depth, the equilibrium state of flow was reached after around 20 seconds. Then, comparisons were made between the numerical and experimental water surface profile as shown in Figure 6.14 – 6.15 for “small case” and Figure 6.16 – 6.17 for “large case”. Note that, with these relative small discharges the water surface profiles become stable as can be seen from Figures 6.14 – 6.17, and the good agreements between calculated and observed results can be observed from the figures implying the ability of the model in simulating the flows over a circular structure. These agreements are similar 83 with the steady analysis in section 6.3. In the second stage, the slightly increased upstream discharge were applied. Figure 6.18-6.21 shows water surface profile with different discharge in which the bold line is the bottom plane and the faint lines are the water surface. Note that the water surface profiles shown in the same figure are based on different time. It is observed that with increase of discharge the water surface oscillation started to appear (Figure 6.18, 6.20). Further increase of discharge resulted in increase of amplitude of variation (Fig. 6.19, 6.21). Moreover, in Exp - 3 and Exp - 4, the amplitude of fluctuations also increases toward downstream. The time history of water level was investigated to determine the frequency of fluctuation. The power spectrum analysis was done for the water depth data. Figures 6.22 – 6.23 show the measured oscillation spectra at point 3 and point 4 in Exp - 3. It can be easily seen from these experimental results that there are some significant contributions of the fluctuation with the periods approximately of 0.065~0.066s. These periods are consistent with the calculated results shown in Figure 6.24 and 6.25. In these graphs, the period of fluctuation in the calculated results is approximately 0.06s. The same picture is also monitored from the calculation for Exp – 4 in Figures 6.26-6.29. In observation, the significant periods are 0.065-0.069(s), while the computed one is 0.06(s). The response pattern in the calculated results is similar to the experimental data in term of frequency. It is furthermore confirmed by the comparison shown in Figures 6.30-6.33, in which the calculated results are for Exp-8. The observation periods are 0.075-0.087 (s) while the simulation one is 0.08 (s). But, it is clearly seen from these figures that the amplitude of oscillation for the calculated results is still larger than that of experiments. 84 0.90 0.95 1.00 1.05 1.10 x (m) 0.45 0.50 0.55 0.60 z (m ) Water surface Time averaged in Exp. Bottom Figure 6.14 Computed water surface profile in Exp1 0.90 0.95 1.00 1.05 1.10 x (m) 0.45 0.50 0.55 0.60 z (m ) Water surface Time averaged in Exp. Bottom Figure 6.15 Computed water surface profile in Exp2 85 0.90 0.95 1.00 1.05 1.10 1.15 1.20 x (m) 0.45 0.50 0.55 0.60 0.65 0.70 z (m ) Water surface Time averaged in Exp. Bottom Figure 6.16 Computed water surface profile in Exp5 0.90 0.95 1.00 1.05 1.10 1.15 1.20 x (m) 0.45 0.50 0.55 0.60 0.65 0.70 z (m ) Water surface Time averaged in Exp. Bottom Figure 6.17 Computed water surface profile in Exp6 86 0.90 0.95 1.00 1.05 1.10 x (m) 0.45 0.50 0.55 0.60 z (m ) Water surface at T=10.00s Time averaged in Exp. Bottom Water surface at T=10.15s Figure 6.18 Computed water surface profile in Exp3 0.90 0.95 1.00 1.05 1.10 x (m) 0.45 0.50 0.55 0.60 z (m ) Water surface at T=10.00s Time averaged in Exp. Bottom Water surface at T=10.15s Figure 6.19 Computed water surface profile in Exp4 87 0.90 0.95 1.00 1.05 1.10 1.15 1.20 x (m) 0.45 0.50 0.55 0.60 0.65 0.70 z (m ) Water surface Time averaged in Exp. Bottom Figure 6.20 Computed water surface profile in Exp7 0.90 0.95 1.00 1.05 1.10 1.15 1.20 x (m) 0.45 0.50 0.55 0.60 0.65 0.70 z (m ) Water surface at T=10.05s Time averaged in Exp. Bottom Water surface at T=10.00s Figure 6.21 Computed water surface profile in Exp8 88 10 15 20 25 30 Frequency (hz) 0.000 0.003 0.005 0.008 0.010 0.013 0.015 0.018 0.020 Po w er sp ec tru m (m m s)2 T = 0. 06 6s Figure 6.22 Power spectrum of water surface displacement at point 3 in Exp – 3 10 15 20 25 30 Frequency (hz) 0.000 0.003 0.005 0.008 0.010 0.013 0.015 0.018 0.020 Po w er sp ec tru m (m m s)2 T = 0. 06 5s Figure 6.23 Power spectrum of water surface displacement at point 4 in Exp – 3 89 0 0.2 0.4 0.6 0.8 1 Time (s) -0.003 -0.002 0.000 0.001 0.003 S ur fa ce di sp la ce m en t( m ) Cal. Exp. Figure 6.24 Comparison of calculated and experimental results at point 3 in Exp - 3 0 0.2 0.4 0.6 0.8 1 Time (s) -0.003 -0.002 0.000 0.001 0.003 S ur fa ce di sp la ce m en t( m ) Cal. Exp. Figure 6.25 Comparison of calculated and experimental results at point 4 in Exp - 3 90 10 15 20 25 30 Frequency (hz) 0.000 0.003 0.005 0.008 0.010 0.013 0.015 0.018 0.020 Po w er sp ec tru m (m m s)2 T = 0. 06 5s T = 0. 06 9s Figure 6.26 Power spectrum of water surface displacement at point 3 in Exp – 4 10 15 20 25 30 Frequency (hz) 0.000 0.010 0.020 0.030 0.040 0.050 0.060 Po w er sp ec tru m (m m s)2 T = 0. 06 5s T = 0. 06 7s Figure 6.27 Power spectrum of water surface displacement at point 4 in Exp – 4 91 0 0.2 0.4 0.6 0.8 1 1.2 1.4 Time (s) -0.003 -0.002 0.000 0.001 0.003 S ur fa ce di sp la ce m en t( m ) Cal. Exp. Figure 6.28 Comparison of calculated and experimental results at point 3 in Exp - 4 0 0.2 0.4 0.6 0.8 1 1.2 1.4 Time (s) -0.003 -0.002 0.000 0.001 0.003 S ur fa ce di sp la ce m en t( m ) Cal. Exp. Figure 6.29 Comparison of calculated and experimental results at point 4 in Exp – 4 92 10 15 20 25 30 Frequency (hz) 0.000 0.010 0.020 0.030 0.040 0.050 0.060 Po w er sp ec tru m (m m s)2 T = 0. 07 8s Figure 6.30 Power spectrum of water surface displacement at point 3 in Exp – 8 10 15 20 25 30 Frequency (hz) 0.000 0.010 0.020 0.030 0.040 0.050 0.060 Po w er sp ec tru m (m m s)2 T = 0. 08 1s T = 0. 07 5s T = 0. 08 7s Figure 6.31 Power spectrum of water surface displacement at point 4 in Exp – 8 93 0 0.2 0.4 0.6 0.8 1 1.2 Time (s) -0.003 -0.002 0.000 0.001 0.003 S ur fa ce di sp la ce m en t( m ) Cal. Exp. Figure 6.32 Comparison of calculated and experimental results at point 3 in Exp – 8 0 0.2 0.4 0.6 0.8 1 1.2 Time (s) -0.003 -0.002 0.000 0.001 0.003 S ur fa ce di sp la ce m en t( m ) Cal. Exp. Figure 6.33 Comparison of calculated and experimental results at point 4 in Exp – 8 94 6.5 2D simulation In this manuscript, in order to examine the 2D characteristics of the flow over circular surface, the 2D equations were applied with the same discretization scheme as described in Chapter 5. The similar boundary conditions to that used in the previous 1D simulation were applied, and the flows conditions were exactly as same with experiment as listed in table 6.1. After getting equilibrium state, the results of 2D simulation are shown in Figure 6.34 – 6.41. Generally, the 2D simulation is consistent with the 1D simulation. In Exp-1 to Exp-3 and Exp-5 to Exp7, the water surface is stable since upstream discharge is relatively small. And in Exp-4 and Exp-8, when the discharge is larger, some oscillation at surface can be observed. 6.6 Summary To investigate the oscillation of the flows over circular surface, the depth-averaged equations based on a generalized curvilinear coordinate attached to bottom were applied. A hydraulic experiment was conducted in laboratory. When the upstream discharge was relatively small, the water surface was stable, and with increased discharge, the oscillation of water surface occurred in the circular bottom region. The intensity of the oscillation increased toward downstream. Both steady and unsteady analysis with small discharges showed the good agreement with experimental results. Both 1D and 2D unsteady simulations showed that, the model has been able to reproduce the phenomena as well as the frequency of the oscillation. However, there are still some discrepancies in term of oscillation amplitude. 95 0.05 0.10 0.15 0.20 0.25 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.34 Carpet plot of water surface in 2D simulation of Exp-1 0.05 0.10 0.15 0.20 0.25 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.35 Carpet plot of water surface in 2D simulation of Exp-2 96 0.05 0.10 0.15 0.20 0.25 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.36 Carpet plot of water surface in 2D simulation of Exp-3 0.05 0.10 0.15 0.20 0.25 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.37 Carpet plot of water surface in 2D simulation of Exp-4 97 0.00 0.10 0.20 0.30 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.38 Carpet plot of water surface in 2D simulation of Exp-5 0.00 0.10 0.20 0.30 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.39 Carpet plot of water surface in 2D simulation of Exp-6 98 0.00 0.10 0.20 0.30 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.40 Carpet plot of water surface in 2D simulation of Exp-7 0.00 0.10 0.20 0.30 z (m ) 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 x (m ) 0.00 0.05 0.10 0.15y (m) Figure 6.41 Carpet plot of water surface in 2D simulation of Exp-8 99 6.7 References 1. Azmad Z. 2005. Flow measument using free overfall in inverted semi-circular channel. Flow. Meas. Instrum., Vol 16, pp 21-26. 2. Davis, A. C., Ellett, B. G. S. and Jacob, R. P. 1998. Flow measurement in sloping channels with rectangular free overfall. J. Hydr. Engrg., ASCE, Vol. 124 (7), pp. 760-763. 3. Dey, S. 1998. End depth in circular channels. J. Hydr. Engrg., ASCE, Vol. 124 (8), pp. 856-862. 4. Dey, S. 2002. Free overfall in circular channels with flat base: a method of open channel flow measurement. Flow. Meas. Instrum., Vol 13, pp 209-221. 5. Dey, S. 2003. Free overfall in inverted semicircular channels. J. Hydr. Engrg., ASCE, Vol. 129(6), pp. 438-447. 6. Dey S. 2005. End depth in U-shaped channels: a simplified approach. J. Hydr. Engrg., ASCE, Vol. 131(6), pp. 513-516. 7. Guo, Y. 2005. Numerical modeling of free overfall. J. Hydr. Engrg., ASCE, Vol. 131(2), pp. 134-138. 8. Pal M. and Goel A. 2006. Prediction of the end-depth ratio and discharge in semi-circular and circular shaped channels using support vector machines. Flow. Meas. Instrum., Vol 17, pp 49-57. 9. Rouse, H. 1936. Discharge characteristic of the free overfall. Civ. Engrg. ASCE, Vol. 6(4), pp. 257-260. 100 Chapter 7 MODEL REFINEMENT 7.1 Preliminary In Chapter 3, a 2D depth averaged model was introduced based on a generalized curvilinear coordinate system, in which the coordinate axis ζ was always perpendicular to the bottom surface. Therefore, the orthogonal conditions in Equation 3.10 and 3.11 were applied to obtain the simplified momentum equations as described in (3.33) and (3.34). Nevertheless, because axis ζ is always perpendicular to the bottom surface thus, for a concave bed there is a possibility that ζ will intersect some where above the surface. To avoid this difficulty, the above model should be applied in a convex topography as in Chapter 4-6 or we limit the applicability of the model in concave surfaces with a shallow depth therefore with these limitations the condition for ζ not to intersect within the computational domain is satisfied as shown in following figure (Figure 7.1). 101 In order to eliminate these limitations from the previuos model, we are now at the stage to derive a more general equation set that will not use the condition for ζ to be perpendicular to the bottom surface, (i.e. ζ could be an arbitrary axis). 7.2 Non-orthogonal coordinate system Coordinate system The new curvilinear coordinate system is set on the bottom surface similar with the previous application but the ζ -axis is arbitrary (as shown in Figure 7.2). The kinematic boundary condition is also derived at the free surface as in (3.23): ηξ ∂ ∂ + ∂ ∂ + ∂ ∂ = hVhU t hW sss (7.1) Depth-average equations: The transformed RAN equations in a new curvilinear coordinates are integrated through the depth from bottom to water surface with respect to ζ -axis, without applying the perpendicular relations, we obtain: Continuity equation ξ ζ Figure 7.1 Illustration for limitation of the model in concave topography 102 01 000 = ∂ ∂ + ∂ ∂ + ∂ ∂ J N J M t h J ηξ (7.2) Momentum equations: ( )ξηηξηξξξηξξξηξ 0200020000 1 Γ+Γ+Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NNMMNM hJJ VM J UM J M t 000 000000 00 2 0 2 0 2 0 0 J dp J dp J G J h b h zzyyxx h zyx ρ τζ ρη ηξηξηξζ ρξ ξξξ ξξ − ∂ ∂++ − ∂ ∂++ −= ∫∫ ( ) ∫∫ ∂∂+∂∂+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ++− hh b zzyyxx dJ d J p J 0000 000000 0 111 ζ ρ τ η ζ ρ τ ξρζξζξζξ ξηξξ (7.3) and ( )ηηηηηξηξηηξξηξ 0200020000 1 Γ+Γ+Γ+Γ+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ NNMMNM hJJ VN J UN J N t 000 000000 00 2 0 2 0 2 0 0 J dp J dp J G J h b h zzyyxx h zyx ρ τζ ρξ ηξηξηξζ ρη ηηη ηη − ∂ ∂++ − ∂ ∂++ −= ∫∫ ( ) ∫∫ ∂∂+∂∂+⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ++− hh b zzyyxx dJ d J p J 0000 000000 0 111 ζ ρ τ η ζ ρ τ ξρζηζηζη ηηηξ (7.4) In order to estimate the pressure distribution along the ζ -axis, the momentum equation in ζ -direction is recast with the assumptions of homogeneity in ζ and shear stress are negligible: ( ) ( ) ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ ++−=Γ+Γ+Γ+Γ ρζζζζ ζζ ηη ζ ηξ ζ ξη ζ ξξ pVVUUVU J zyx 22222 J 1 J G1 ( ) ( ) ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ ++−⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ ++− ρξξζξζξζρηηζηζηζ p J p J zzyyxxzzyyxx 11 (7.5) or ( ) ( ) ( ) ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +++⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ +++⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ ++ ρξξζξζξζρηηζηζηζρζζζζ ppp zzyyxxzzyyxxzyx 222 103 ( )ζηηζηξζξηζξξζ Γ+Γ+Γ+Γ−= 22G VVUUVU (7.6) To solve equation (7.6) by iteration method, the first estimation of pressure term is given as: ( ) ( ) ( )ζηηζηξζξηζξξζρζζζζ Γ+Γ+Γ+Γ−=⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ ++ 22 0 2 0 2 0 2 0 G VVUUVU p zyx (7.7) hence, taking integration from ζ to h gives ( ) ( ) ( )[ ]( )ζζζζρ ζηηζηξζξηζξξζ −Γ+Γ+Γ+Γ−++= hVVUUVUp zyx 22 2 0 2 0 2 0 0 G1 ( ) ( )ζ−= hf p0 Consequently, ( ) ( ) ( ) ( ) ( )( ) ( )ξζξξξζρξ ∂ ∂ ⋅− ∂ ∂ = ∂ ∂ + ∂ ∂ −=⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ 000 00 p pp p fhfhf f hp (7.8) and similarly, ( ) ( ) ( ) ( ) ( )( ) ( ) η ζ ηηη ζ ρη ∂ ∂ ⋅− ∂ ∂ = ∂ ∂ + ∂ ∂ −=⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ ∂ ∂ 000 00 p pp p fhfhf f hp (7.9) Substituting Equations (7.8) and (7.9) into (7.7) and taking integration we attain the equation for the new estimation of pressure distribution as follows ( ) ( )( )( ) ( ) +⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − ∂ ∂ −− ∂ ∂ ++ 22 220 0 0 ζ ξζξξζξζξζ hfhhf ppzzyyxx ( ) ( )( )( ) ( ) ( ) ( ) ρ ζζζζ η ζ η ηζηζηζ 1 0 222 220 0 0 22 phfhhf zyx p pzzyyxx ++−⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − ∂ ∂ −− ∂ ∂ ++ ( )[ ]( )ζζηηζηξζξηζξξζ −Γ+Γ+Γ+Γ−= hVVUUVU 22G as a result, 104 ( ) ( ) ρ ζζζ 1 0 222 p zyx ++ ( )[ ]( )ζζηηζηξζξηζξξζ −Γ+Γ+Γ+Γ−−= hVVUUVU 22G ( ) ( )( )( ) ( ) ⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − ∂ ∂ −− ∂ ∂ +++ 22 220 0 0 ζ ξζξξζξζξζ hfhhf ppzzyyxx ( ) ( )( )( ) ( ) ⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ − ∂ ∂ −− ∂ ∂ +++ 22 220 0 0 ζ η ζ η ηζηζηζ hfhhf ppzzyyxx (7.10) Taking integration Equation (7.10) gives ( ) ( ) ( )[ ] 2 2 22 0 1 0 222 hGVVUUVUdp h zyx ζζ ηη ζ ηξ ζ ξη ζ ξξζρζζζ −Γ+Γ+Γ+Γ=++ ∫ ( ) ( )( ) ( ) ⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ∂ ∂ − ∂ ∂ +++ 32 302 0 0 hfhhf ppzzyyxx ξξξζξζξζ ( ) ( )( ) ( ) ⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ∂ ∂ − ∂ ∂ +++ 32 302 0 0 hfhhf ppzzyyxx ηη ηζηζηζ (7.11) Equation (7.11) describe the pressure distribution along the ζ -axis. It is easily seen that, two last terms are added in the pressure distribution and if ζ normal to the bottom surface, the distribution reduced as the first term - a term representative for the hydrostatic pressure and effects of centrifugal force due to bottom curvature shown in Chapter 3. The pressure at the bottom is then obtained: ( ) ( )[ ] +−Γ+Γ+Γ+Γ++=⎟⎟⎠ ⎞ ⎜⎜⎝ ⎛ hGVVUUVUp zyxb ζζ ηη ζ ηξ ζ ξη ζ ξξζζζρ 22 222 1 ( ) ( )( ) ( ) ⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ∂ ∂ − ∂ ∂ +++ ξξξζξζξζ 02 0 0 2 p pzzyyxx fhhfh ( ) ( )( ) ( ) ⎪⎭ ⎪⎬ ⎫ ⎪⎩ ⎪⎨ ⎧ ∂ ∂ − ∂ ∂ +++ ηη ηζηζηζ 02 0 0 2 p pzzyyxx fhhfh (7.12) 105 The contravariant components of gravitational vector are derived similar as in Chapter 4: ( )ηζζηξ yxyxgJG −−= (7.13a) ( )ξζζξη yxyxgJG −−= (7.13b) ( )ξηηξξ yxyxgJG −−= (7.13c) 7.3 Application In order to verify the validity and the applicability of the new general equations, in this first version of manuscript, the dam-break flows in a horizontal and sloping channel bed are investigated. In this application, the 1-D flow is considered using the new coordinate system as described in Figure 7.2. The parameters 1α and 2α define the slope of the channel and the angle of the ζ -axis with the channel bottom (also means the angle between two new axes). By changing the value 1α , it is easily to set the slope of channel and 01 =α defines for the case of horizontal bed. Figure 7.2 Definition sketch of new generalized coordinate system z ζ x ξ 2α 1α 106 In this research, the second and third terms in equation (7.6) are neglected, it gives ( ) ( )2 2 2 2 20 0 0 Gx y z p U UV VU Vζ ζ ζ ζ ζξξ ξη ηξ ηηζ ζ ζ ζ ρ⎛ ⎞∂+ + = − Γ + Γ + Γ + Γ⎜ ⎟∂ ⎝ ⎠ (7.14) Consider 1D flow (i.e. 0V = ) and the calculation gives 0ζξξΓ = then pressure distribution in (7.14) becomes hydrostatic: ( ) ( ) ( )2 2 2 2 2 20 0 0 0 0 0 G 1 G x y z x y z p p h ζ ζ ζζ ρ ρζ ζ ζ ζ ζ ζ ⎛ ⎞ ⎛ ⎞∂ = ⇒ = −⎜ ⎟ ⎜ ⎟∂ + + + +⎝ ⎠ ⎝ ⎠ hence, the pressure at the bottom is obtained: ( )2 2 20 0 0 G b x y z p hζ ρ ζ ζ ζ ⎛ ⎞ =⎜ ⎟ + +⎝ ⎠ (7.15) and ( ) 2 2 2 2 0 0 0 0 1 G 2 h x y z p d hζζ ρ ζ ζ ζ ⎛ ⎞ =⎜ ⎟ + +⎝ ⎠∫ (7.16) Substituting (7.15-7.16) into (7.3) and neglecting the effect of internal stresses gives 2 2 22 0 0 00 0 0 0 0 0 0 h x y zM pM UM h G d t J J J h J J ξ ξξ ξ ξ ξ ξ ζξ ξ ρ + +Γ⎛ ⎞ ⎛ ⎞∂ ∂ ∂ + + = −⎜ ⎟ ⎜ ⎟∂ ∂ ∂⎝ ⎠ ⎝ ⎠ ∫ ( )0 0 0 0 0 0 0 0 1 b x x y y z z b p J J ξτξ ζ ξ ζ ξ ζ ρ ρ ⎛ ⎞ − + + −⎜ ⎟⎝ ⎠ or ( ) 2 2 22 0 0 00 2 2 2 2 0 0 0 0 0 0 0 0 1 G 2 x y z x y z U hM UM h G h t J J J J J ξ ξξ ξ ζξ ξ ξ ξ ξ ζ ζ ζ ⎡ ⎤+ +Γ⎛ ⎞ ⎛ ⎞∂ ∂ ∂ ⎢ ⎥+ + = −⎜ ⎟ ⎜ ⎟ ⎢ ⎥∂ ∂ ∂ + +⎝ ⎠ ⎝ ⎠ ⎣ ⎦ ( ) ( ) 0 0 0 0 0 0 2 2 2 0 00 0 0 1 G x x y y z z b x y z h J J ξ ζξ ζ ξ ζ ξ ζ τ ρζ ζ ζ + + − − + + (7.17) Equation (7.17) is solved in plugging with continuity equation (7.2) and the geometric 107 variables are calculated from the mesh depicted in Figure 7.2, to investigate the dam-break flow in horizontal and slopping channels. Figure 7.3 and 7.4 show the dam-break flow water surface at different times in a dried and wetted bed for horizontal channel. The bores are visually observed, and they were compared with the conventional depth-averaged model (i.e., Kuipers and Vreugdenhill equations). The good agreement was achieved in the bore of both cases: dried and wetted bed , but there are still some bias at the upstream of the dam, that might be caused by the difference in boundary conditions at upstream end due to the difference of the grid. It infers that, the new model can be applied in such cases that the third axis is not necessary always perpendicular to the bottom surface. It also was emphasized in Figure 7.5-7.6, which show the dam-break flows in channels with slop 030 . 108 2.0 4.0 6.0 8.0 x (m) 0.0 2.0 z (m ) Figure 7.3 Calculated water surface profile at different time steps of dam break flow in a dried-bed horizontal channel: mHini 0.1= ; 01 30=α ; 0 2 60=α 2.0 4.0 6.0 8.0 x (m) 0.0 2.0 z (m ) Figure 7.4 Calculated water surface profile at different time steps of dam break flow in a wetted-bed horizontal channel: mHinimHini downup 5.0;0.1 == ; 0 1 0=α ; 0 2 60=α T = 0.00s T = 0.20s T = 0.40s T = 0.60s Channel bed T = 0.00s T = 0.20s T = 0.40s T = 0.60s Channel bed 109 2 4 6 8 10 12 14 16 18 x (m) 0.0 0.5 1.0 1.5 z (m ) Figure 7.5 Comparison of water surface profile for dried horizontal channel at different times: T = 0.4s; 1.0s; and 1.6s; 2 4 6 8 10 12 14 16 18 x (m) 0.0 0.5 1.0 1.5 z (m ) Figure 7.6 Comparison of water surface profile for wetted horizontal channel at different times: T = 0.4s; 0.6s; and 0.8s; Initial depth Conventional depth-averaged model New model Initial depth Conventional depth-averaged model New model 110 2.0 4.0 6.0 8.0 10.0 x (m) 0.0 2.0 4.0 6.0 z (m ) Figure 7.7 Calculated water surface profile at different time steps of dam break flow in a dried-bed sloping channel: mHini 0.1= ; 01 30=α ; 0 2 60=α 2.0 4.0 6.0 8.0 x (m) 0.0 2.0 4.0 6.0 z (m ) Figure 7.8 Calculated water surface profile at different time steps of dam break flow in a dried-bed horizontal channel: mHinimHini downup 5.0;0.1 == ; 0 1 30=α ; 0 2 60=α T = 0.0s T = 0.2s T = 0.4s T = 0.6s Channel bed T = 0.0s T = 0.2s T = 0.4s T = 0.6s Channel bed 111 Chapter 8 CONCLUSIONS In this research, a depth-averaged model of open channel flows in a generalized curvilinear coordinate system over an arbitrary 3D surface was developed. This model is the inception for a new class of the depth-averaged models classified by the criteria of coordinate system. This work focused on the derivation of governing equations with examples of application in hydraulic amenity. The major results and contributions of this research are summarized as follows. 1. Mathematical model - The original RAN equations were firstly transformed into the new generalized curvilinear coordinate system, in which two axes lay on an arbitrary surface and the other perpendicular to that surface. The depth-averaged continuity and momentum equation were derived by integrating the transformed RAN equation and substituting the kinematic boundary condition at the free surface. - The pressure distribution was obtained by integrating the momentum equation along the perpendicular axis as a combination of hydrostatic pressure and the effect of 112 centrifugal force caused by bottom curvature. 2. Flow with air-core vortex at a vertical intake - Firstly, 1D analytical analysis was applied to calculate the water surface profile of the flows with air-core vortex at a vertical intake. The use of curvilinear coordinate avoids invalidation of the Kelvin’s theorem in the core of vortex, therefore the conservation of circulation can be applied for whole flow. The comparisons of calculated data using the model and an empirical fomula show that the proposed model yields reliable results in predicting the critical submergence of the intake without any limitation of Froude number, a problem that most of existing models cannot escape. - Secondly, 2D unsteady analysis was applied to simulate the water surface profile and also showed the applicability of the model in computing the flow into intake with air-entrainment. 3. Flows over circular surface - A hydraulic experiment was conducted in the laboratory to measure the surface profile of the flows over circular surface. It was pointed out that: with small discharge, the free surface remained stable but some oscillations occurred when the discharge was increased. The intensity of the oscillation was observed to increase toward the downstream direction of the flow. - Both 1D steady and unsteady analysis were applied to estimate the water surface profile of the flows with small discharge. Good agreements in comparison with measured data were achieved for both approaches. Results from 1D unsteady model demonstrate the fluctuation at the region with circular bottom with large discharge. Spectrum analysis was applied to investigate the characteristics of the oscillation. The comparison with observed data showed the fairly good agreement in term of 113 frequency or period but some discrepancies in magnitude of fluctuation were still observed. 2D simulation was conducted, and also could reproduce the phenomena as well. 4. Model refinement - In order to apply the model for the case with highly concave topography, the new depth-averaged equations were derived without using the perpendicular conditions. The improved model therefore has the ability in more general geometry of the bottom surface - A simple example was presented to verify the model with dam-break flows in horizontal and sloping channels.