(College of Water Resources & Environment, Hohai University, Nanjing, 210098)

Abstract: In tidal areas, natural land boundary is complex and underwater topography varies acutely due to influence of upstream runoff and outer tide. The simulation and forecast of water current and mass transportation play an important role in practical engineering. According to the situation of irregular natural boundaries in tidal region, unstructured triangular grid arrangement is applied to suit for complex condition. A finite difference method with alternating directional implicit scheme for triangular grid is established in this paper. The model has been applied in calculation of flow and concentration fields for Nantong reach of the Yangtze River. It is satisfied that the calculated values are in agreement with observed data.

 

Keywords: tidal areas; triangular grid; finite difference method; numerical simulation

Up to now, the reports for solving 2D depth-integrated equations of flow and mass transportation in tidal region have been introduced in large numbers. From the beginning of Leendertse (1970), the finite difference method (FDM) has been widely applied by Falconer (1989), Vincenz Casulli (1988,1990), He Shaoling (1985) etc due to its simple scheme, high efficiency and good stability. However, the rectangular mesh layout in FDM can not well fit for complex land boundary conditions in tidal areas. In order to get high computational precision, the scale of rectangular meshes must be reduced, so that the time-step must be cut down, as a result, the efficiency of the model becomes low.

For the finite element method (FEM), the arbitrarily triangular and quadrilateral grid layout suited for natural boundaries can be used. Grotkop (1973), Gray (1977), Shi Jinsong (1987), Tan Weiyan (1985), Geng Zhaoquan (1984) etc have obtained successful calculation result, but FEM has some disadvantages such as complex node number, slow convergence speed and more limited time-step.

According to irregular boundaries of tidal region, a body-fitted technique is presented. The body-fitted grid can be approximately classified into structured and unstructured grid.

Several researchers such as Thompson (1980), W. Rodi (1989), Wijbenga (1985), Hausser (1986), Jothi Shankar (1997), Binliang (1996), HUA Zulin (1994,1999) etc have reported some transformation from practical irregular water region plane into regular one by mathematical conversion (algebraical or differential conversion). This method has a good precision suited for real boundaries, but it needs a large amount of generation of curvilinear grids. At present, unstructured grid (e.g. triangular grid) has been applied to the numerical simulation of water current and mass transportation in tidal areas, for example, by Liu Zhi (1987), Dou Xiping (1995), Chen Yang (1998), K. Anastasiou (1997) etc. Because of the main use of explicit scheme, the stability and time-step can be limited for these models.

Based on the triangular grid, a finite difference method model with alternating directional implicit scheme is established in this paper. The model has been applied in simulation of flow and concentration fields for Nantong reach of Yangtze River, and satisfactory results of calculation have been obtained.

For a widely open water areas, the vertical velocity component w is smaller than other two directional velocity components u and v respectively, i.e. w<<u and w<<v. In addition, the vertical acceleration component of the flow is much smaller than gravitational acceleration. Consequently, the hydrostatic pressure approximation (or shallow water approximation) can be introduced. Integrating 3d Reynolds equations along water depth, the 2D basic equations are obtained as follows

Continuity equation:

                        (1)

Momentum equations:

                            (2)

                          (3)

Mass transportation equation:

                         (4)

where, x, y and t --- variables of space in the Cartesian coordinate system and time variable respectively; u, v --- depth-averaged velocity components in the directions of the x- and y-axes respectively; ξ--- tidal level; H --- water depth; f --- geostrophic inertial force; ρ--- density of water; t _sx,t _sy --- shear stresses in the directions of the x- and y- axes of wind; t _{bx,t by}--- shear stresses of the river bottom in the directions of the x- and y- axes; C --- concentration; D_x,D_y --- depth-averaged diffusion coefficients in the directions of the x- and y- axes; Sc --- source and sink; n --- depth-averaged turbulent viscosity in the directions of the x- and y- axes.

The formulas of calculation for some parameters in the model are as follows

(1) Shear stresses of wind:

where r _a --- density of air; C_f --- drag coefficient of the interface between water and air，depending on the velocity of wind; w_x, w_y --- components of the wind vector in the directions of the x- and y- axes respectively.

(2) Sheer stress of bottom:

where ，n --- manning roughness coefficient.

(3) Turbulent viscosity: semi-experiential formula is used, ，L_m --- mixing length, selected as L_m =0.1H in this paper.

(1) Closed boundaries (land boundaries)_: At the closed boundaries the normal velocity component equals zero, i.e. , and ;

(2) open boundaries (water boundaries)_: For the boundaries of the flow field, ; or ,.

For the concentration boundaries,

at the inflow boundaries:

at the outflow boundaries: （s --- streamline）.

The finite difference method model with alternative directional implicit scheme for 2D open channel flow was given early by Leendertse (1970). The model for rectangular mesh has been applied widely due to its easy solution and good stability. In this paper, Leendertse’s alternative directional implicit (ADI) scheme is used for the discretization of basic equations in the triangular meshes. The process of discretization can be divided into two steps as following

(1) The first time-step ( from to ())

(2) The second time-step ( from () to ())

then, the partial differential equations can be discretized

(1) The first time-step ( from to ())

(2) The second time-step ( from () to ())

The Nantong reach located in the tidal current of the Yangtze River is affected by runoff and tide. The underwater topography shown in Fig.1 is quite undulatory with a lot of shoals and submerged reefs. The current is in an irregular semi-diurnal tide. The average period of tide is nearly 12.4 hours, in which flood and ebb take about 4.15 hours and 8.25 hours respectively. The average tide range is 2.68 meters. The mean annual maximum and minimum tide levels are 5.66 meters and 0.93 meters respectively.

The length of computational region is about 26.80 km, and the average width is 9.60 km approximately. The distribution with 9408 triangular grids shows in Fig.2, in which the densely scattered grids are distributed in the area for acutely undulatory topography. The calculational model can be verified using observed data for the spring and neap tides in January of 1999, and the calculational results are as follows:

The comparison of calculated velocity values with observed data at the observed point A1 (see Fig.1) for the spring and neap tides are shown in Fig.3, and both relations are quite in agreement. The flow fields of maximum flood and ebb during spring tide are shown in Fig.4 and Fig.5 respectively.

Fig.1 The underwater topography of Nantong reach

Fig.2 The distribution of triangular meshes

Fig.3 The comparison of calculation velocity values with observed data

Fig.4 The flow field of maximum flood (spring tide)

Fig.5 The flow field of maximum ebb (spring tide)

 

The drogue route observations have been made in the calculated region during the periods of flood and ebb at twice, i.e. in January 13, 1999 (neap tide) and January 17, 1999 (spring tide) respectively. The drogue route is determined using Lagrangian method and compared with observed data in the flow field, as shown in Fig.6. The fact of the small difference indicates that the hydrodynamic simulation results in the tidal areas are believable.

Fig.6  A comparison between calculation and observation for the drogue route

The fluorescent dye tracer 5 kg in weight was put in the flow near the waste discharge of the calculated areas during the period of ebb, then the diffusion cloud of the tracer was live followed, furthermore the maximum concentration of the diffusion cloud centre was observed. In general, the comparison between observed and predicted values at twice ebbs (i.e. in January 18 and 19, 1999), as shown in Fig.7, is satisfactory.

Fig.7 The relations of the maximum concentration for the diffusion cloud centre with time

Fig.8 shows the concentration distribution of waste discharge for the maximum ebb. It can be seen that comparison of predicted values with observed values is in a good agreement.

Fig.8   The concentration field of waste discharge (maximum ebb)

The triangular mesh arrangement is complete suitable for complex and irregular boundaries in tidal areas. The triangular grid method is rather simpler and more convenient than the body-fitted method in the distribution of meshes. In addition, the size and distribution of the triangular grids can be flexible selected according to the practical condition. Therefore, the method can be applied for the numerical simulation of flow and mass transportation problems having complex and irregular boundaries in tidal region.

(1) Based on the triangular meshes, the finite difference method model with alternating directional implicit scheme has been established in this paper. This model still retains some advantage of FDM for rectangular mesh such as good stability, easy solution, fast convergence.

(2) By means of the typical calculational case of Nantong reach of yangtze River, the flow and concentration fields affected by tide, runoff, underwater topography and etc. can be calculated successfully. The computed values can be verified by considering several factors from different point of view. The comparison shows that the calculated values are in agreement with observed data satisfactorily.

[1] J. J. Leendertse, “A water quality simulation model for well-mixed estuaries and coastal seas”, Vol.1, Principle of Computation, The Rand Corporation RM-6230, 1970.

[2] R. A. Falconer, P. Goodwin and R.G.S. Matthew(Eds.), Hydraulic and Environmental Modelling of Coastal , Estuarine and River Water, Gower Technical, Aldershot,1989.

[3] Vincenzo Casulli, “Semi-implicit finite difference method for the two-dimensional shallow water equation”, J. Comput. Phys., 1990, 86: 56-74.

[4] He Shaoling, Lin Bingnan, “Fractional steps method for 2d tidal flow computation”, Journal of Oceanography (In Chinese), 1984, 6(2): 260-271.

[5] Vincenzo Casulli, Filippo Notamicola, “An Eulerian-Lagrangian method for tidal current computation”, In: Computer Modelling in Ocean Engineering (Schrefler & Zienkiewocz, eds), Ballkema, Rotterdam,1988, 237-244.

[6] G. Grotkop, “Finite element analysis of long period water waves”, Comp. Meth. In Appl. Mech. And Eng., 1973, (2): 147-160.

[7] W. G. Gray, “An efficient finite element scheme for two-dimensional surface computation”, Proceeding of the first Int. Conf. On Finite Elements in Water Resources, Pentech Press, 1977, 433-449.

[8] Shi Jinsong, “Conservational finite element method for two dimensional unsteady flow”, Journal of Hydraulic Engineering (In Chinese), 1987, (10): 41-48.

[9] Tan Weiyan, Zhao Dihua, “Finite element method and program packet for two-dimensional subcritical unsteady shallow flow”, Journal of Hydraulic Engineering (In Chinese), 1984, (10): 1-13.

[10] Geng Zhaoquan, “Explicit up-wind finite element model for two-dimensional unsteady flow”, The Proceedings of first National Hydrodynamics Conference (In Chinese), China Ocean Press; Beijing, 1984.

[11] J. F. Thompson, “Numerical solution of flow problems using body-fitted coordinate system”, in W. Kollmann(ed.), Computational Fluid Dynamics, Hemisphere, New York,1980.

[12] W. Rodi, S. Majumdar and B. Schonung, “Finite volume methods for two-dimensional incompressible flows with complex boundaries”, Comput. Methods Appl. Mech. Eng., 1989, 75,369-392.

[13] J. Hausser, H.G. Paap and D. Eppel, “Boundary conformed coordinate system for fluid flow problems”, Numerical Grid Generation in Computational Fluid Dynamics, Ed. C. Taylor and J. Hausser, Pineridge Press, Swansea, 1986.

[14] J.H.A. Wijbenga, “Determination of flow pattern in river with curvilinear coordinate”, Congress IAHR, Melbourne, Australia, 1985.

[15] Jothi Shankar N, Hin-Fatt Cheong and Chun-Tat Chan, “Boundary fitted grid models for tidal motions in Singapore coastal waters”, J. of Hydraulic Research, 35(1), 1997, 3-19.

[16] Binliang Lin, Simon N, Chandler-Wilde, “A depth-integrated 2D coastal and estuarine model with conformal boundary-fitted mesh generation”, International Journal for Numerical Methods in Fluids, Vol.23, 1996,819-846.

[17] HUA Zulin, Bian Hua, “Selection of adjusting factors in the generation of curvilinear grid”, Journal of Hohai University(Natural Sciences) (In Chinese), Vol.27, No.2 ,1999,40-44.

[18] HUA Zulin, “Numerical Simulation of flow field in tidal estuary using body-fitted coordinate system”, China Ocean Engineering, Vol.8, No.4, 1994, 447-455.

[19] Liu Zhi, Lin Bingnan and He Shaoling, “The application of irregular triangular meshes in 2-D unsteady flow”, Journal of Hydraulic Engineering (In Chinese),1987, No.9, 25-33.

[20] Dou Xiping, Li Tilai, “Application of triangular grid generation method in tidal mathematical models of coastal engineering”, Journal of Nanjing Hydraulic Research Institute(In Chinese), 1995, No.3, 65-69.

[21] Chen Yang “Two-dimensional mathematical model research and application for navigable flow condition in dam region”, Tianjin Research Institute of Water Transport Engineering (In Chinese), 1998, Rep. No. 95-06-01-04.

[22] K. Anastasiou and C. T. Chan, “Solution of the 2D shallow water equations using the finite volume method on unstructured triangular meshes”, International Journal for Numerical Methods in Fluids, Vol.24, 1997, 1225-1245.

[23] A. S. Sens and G. D. Mortchelewicz, “Implicit schemes for unsteady Euler equations on triangular meshes”, International Journal for Numerical methods in fluids, Vol.18, 1994, 647-668.