Modeling the Impact of Wind Drag Coefficient on Wind-Driven Currents in Lake Taihu, China

Abstract

The wind drag coefficient, $C_d$, has a great influence on the numerical results obtained from shallow lakes. To analyze the modeling impacts of $C_d$ on wind-driven currents, a series of numerical simulations of Lake Taihu were conducted at three grid resolutions (800 m × 800 m, 400 m × 400 m, and 100 m × 100 m) using the empirical formulae of Flather (F76), Large and Pond (LP81), Large and Yeager (LY04), Andreas (A12), and Gao (G20). The G20 formula produced the optimum results of all the formulae for both the water level and velocity simulations; however, the grid resolution was found to have a significant influence on simulation in G20 cases. Thus, the G20 formula is only recommended when using a high-resolution grid to meet the accuracy requirements of analyzing wind-driven currents in the numerical modeling of Lake Taihu. A combination of the A12 formula and a coarse grid is preferred when taking computational efficiency into consideration.

1. Introduction

Lake Taihu, the third largest freshwater lake in China with a surface area of 2388 ${\mathrm{km}}^2$ and an average water depth of 1.9 m, is a typical fetch-limited shallow lake with wind-driven currents (Tsanis 1989 [1], Wu 2006 [2], Li 2011 [3]). Its basin accounts for 7.0% of China’s water resources and supports 18.8% of the national GDP; therefore, it is vital that the water quality of Lake Taihu is maintained. Eutrophication has been the most serious environmental issue of the around-Taihu area over the past two decades. As nutrients are released into the lake’s water body from both intenal and extenal resources and varying both spatially and temporally, one of the accepted research methods is to construct high resolution numerical models of Lake Taihu, taking into account the effects of sediment release and the external loads from the tributaries (Mao 2008 [4]). The movement of wind-driven currents is considered to be the main factor affecting the water quality of large shallow lakes due to sediment redistribution, contaminant transport, and nutrient resuspension (Li 2017 [5], Chu 2019 [6]). In numerical models, the influence of wind-driven current is traditionally imposed through the wind stress term, expressed as $\tau =\rho C_d{U}_{10}^2$, where $\tau$ is the wind stress; $\rho$ is the air density; ${\mathrm{kg/m}}^3$; $U_{10}$ is the wind speed at 10 m height measured over the free surface, m/s; $C_d$ is the wind stress drag coefficient, representing the momentum flux (Kranenburg 2010 [7], Li 2012 [8], Wróbel-Niedźwiecka 2019 [9], Zhang 2023 [10]). $C_d$ is invariably the most important physical parameter for determining the accuracy of wind-driven models.

As the collection of observed data has increased over the past hundred years, many empirical parameterization methods used to obtain the value of $C_d$ have been formulated, and they can be divided into several groups as follows:

1. Ignore changes in the roughness of free surfaces and assume that the value of $C_d$ is constant in both time and space (Morales-Marin 2017 [11], Xu 2019 [12], Oliveira 2019 [13]);

2. Assess a variety of observed data to establish a linear relationship between $C_d$ and the wind speed (Sheppard 1972 [14], Flather 1976 [15], Large and Pond 1981 [16], Geernaert 1987 [17], Ataktürk and Katsaros 1999 [18]). These kinds of parameterization methods can be written in the general form given by Guan and Xie (2004) [19];

3. Consider the effects of surface tension and swell at low wind speeds to obtain the nonlinear relationship between the $C_d$ value and the wind speed (Yelland and Taylor 1996 [20], Large and Yeager 2004 [21], Andreas 2012 [22]);

4. Take both wind speed and sea state parameters into consideration. Based on the theory of Janssen (1991) [23], some $C_d$ parameterizations have been conducted using an ocean wave model (Tolman and Chalikov 1996 [24]). These mathematical expressions of $C_d$ are often used in numerical simulations of storm surges. Multiple nonlinear parameterization methods of $C_d$ are indicated, and they are related not only to wind speed but also to wind blowing fetch and water depth (Mueller and Veron 2008 [25], Gao 2020 [26]).

In the first three parameterization methods, the value of $C_d$ is either constant or is only related to wind speed. However, many studies have suggested that wind stress in ocean conditions is related not only to the wind speed, but also to the sea state (Mastenbroek 1993 [27], Smith 1992 [28]). Smith (1992) [28] analyzed data measured in the North Sea and found that the $C_d$ value was strongly related to the wave age. By analyzing the observed data, Ebuchi (1990) [29] demonstrated that the $C_d$ value is negatively correlated with wave age; wave-induced stress in mature wind waves is larger than that of young wind waves. When determining the wind-induced current in fetch-limited shallow lakes, it was found that the state variables affecting the $C_d$ value are different from those in the ocean (Gao 2019 [30]). Elfouhail (1997) [31] studied the influence of water state parameters on wind stress and found that the effect of wind blowing fetch on wind stress value should not be ignored. It was found that when the wind speed is constant, the $C_d$ value has a positive correlation with wind blowing fetch and a negative correlation with water depth. Most of the parameterization methods used to calculate $C_d$ are based on a marine environment; previous studies have concentrated more on the effects that these parameterization methods have on wind-induced currents in this context (Chu 2019 [6], Fang 2018 [32], Moon 2004 [33]). As such, the different $C_d$ parameterization methods lack applicability in the context of wind-induced currents in fetch-limited shallow lakes. This study aims to investigate the influence of various $C_d$ parameterization methods on wind-induced current modeling in a typical fetch-limited shallow lake—Taihu Lake. For this purpose, a high-resolution unstructured model is developed based on the TELEMAC3D module and evaluated against the data obtained from Taihu Lake. Five $C_d$ parameterization methods are chosen for their importance in the field for testing in Taihu Lake:

1. The linear relationship developed by Flather (1976) [15];

2. The linear relationship developed by Large and Pond (1981) [16];

3. The nonlinear relationship developed by Large and Yeager (2004) [21];

4. The nonlinear relationship developed by Andreas (2012) [22];

5. The nonlinear relationship considering the influence of wind blowing fetch and water depth developed by Gao (2020) [26].

The details of the five parameterization methods are described in Section 2.2. The parameterization method developed by Gao (2020) [26] requires the spatial distribution of wind blowing fetch and water depth to be known, and the grid resolution may have a great influence on the numerical results; therefore, additional cases with different grid resolutions (800 m × 800 m, 400 m × 400 m, and 100 m × 100 m) were set to study the effect of each resolution on the numerical results for the five parameterization methods. Descriptions of the methodology of this study are given in Section 2. Section 3 introduces the numerical setting and model configuration. Section 4 describes the results of the model in response to the different $C_d$ parameterization methods and grid resolutions. A summary and our conclusions are given in Section 5.

2. Methodology

2.1. Description of 3D Numerical Model

We chose to build our numerical model of wind-induced current and waves based on the OPEN TELEMAC-MASCARET system https://www.opentelemac.org (accessed on 5 December 2022), which was developed by the Company Electricité de France (EDF), Paris, France. The TELEMAC3D module is a three-dimensional hydrodynamic model which solves the Navier–Stokes equations using finite element techniques and uses the sigma coordinate for vertical discretization. The unstructured mesh used by the Telemac3D module is especially suitable for numerical simulations of shallow lakes possessing a complex geometry and dynamic physical processes (Costi 2018 [34], Kopmann and Markofsky 2000 [35]). The non-hydrostatic Navier–Stokes equations used in TELEMAC3D can be written as follows:

(1)$\nabla ·\mathit{u}=0$

(2)${\displaystyle \frac{\partial \mathit{u}}{\partial t}}+\left(\mathit{u}·\nabla \right)\mathit{u}=-{\displaystyle \frac{1}{\rho }}\nabla p+{\displaystyle \frac{1}{\rho }}\nabla ·(\mu \nabla \mathit{u})+\mathit{G}+\mathit{F}$

(3)$p=p_{air}+\rho g(Z_s-z)+\rho g{\int }_z^{z_s}dz+p_d$

where t is time; u represents velocity; $\rho$ and ${\rho }_{air}$ are the water density and air density, respectively; p and $p_d$ are the total pressure and dynamic pressure, respectively; G is the gravity vector; g is the acceleration due to gravity; $\mu$ is the dynamic molecular viscosity; F is the source terms denoting the bottom friction and wind force; p is the atmospheric pressure; $Z_s$ and z are free surface elevation and depth per layer, respectively. Equations are spatially discretized using the finite element technique. In case the influence of wind is taken into account, wind forcing appears in the equations through the condition at the surface:

(4)$\mu {\displaystyle \frac{\partial {\mathit{u}}_h}{\partial n}}=p_{air}+C_d{U}_{10}{\mathit{U}}_{\mathbf{10}}$

where $u_h$ is the horizontal velocity at the free surface; $C_d$ is the wind drag coefficients; $U_{10}$ is the wind velocity measured 10 m above the free surface. The k-$ϵ$ model is employed to compute the horizontal and vertical turbulence; the turbulent kinetic energy k and the dissipation of turbulent kinetic energy $ϵ$ at the free surface consider the influence of wind and are computed as follows:

(5)$k{|}_{z=\eta }={\left({\displaystyle \frac{C_d{\rho }_{air}{\mathit{U}}_{10}^2}{\rho \sqrt{C_{\mu }}}}\right)}^2$

(6)$\epsilon {|}_{z=\eta }={\displaystyle \frac{{k|}_{z=\eta }}{\alpha h}}$

where $C_{\mu }$ is the Prandtl–Kolmogorov constant, $\alpha \approx 0.8$, h is the local water depth.

2.2. Choice of $C_d$ Parameterization Methods

This study investigates the influence of $C_d$ parameterization methods on modeling wind-induced current. The five methods used in this study are described below.

1. Flather (1976)-F76 [15]:

   (7)$\begin{array}{c}\hfill \begin{array}{c}\hfill C_d\times {10}^3=\left\{\begin{array}{c}0.565,0\mathrm{m}/\mathrm{s}\le U_{10}\le 5\mathrm{m}/\mathrm{s}\hfill \\ -0.12+0.137U_{10},5\mathrm{m}/\mathrm{s}\le U_{10}\le 20\mathrm{m}/\mathrm{s}\hfill \\ 2.513,20\mathrm{m}/\mathrm{s}\le U_{10}\hfill \end{array}\right.\end{array}\end{array}$

   Flather (1976) [15] parameterized $C_d$ based on the observed data during a storm surge in the northwest European continental shelf. This parameterization was incorporated into the default TELEMAC2D and TELEMAC3D model and is recommended by the Institute of Oceanographic Sciences (Surrey, UK).

2. Large and Pond (1981)-LP81 [16]:

   (8)$\begin{array}{c}\hfill \begin{array}{c}\hfill C_d\times {10}^3=\left\{\begin{array}{c}1.2,0\mathrm{m}/\mathrm{s}\le U_{10}\le 11\mathrm{m}/\mathrm{s}\hfill \\ 0.49+0.065U_{10},11\mathrm{m}/\mathrm{s}\le U_{10}\hfill \end{array}\right.\end{array}\end{array}$

   The $C_d$ parameterization method developed by Large and Pond (1981) [16] is the default option in the FVCOM model. It showed good performance when measuring ocean circulation and wind-induced current in Okeechobee Lake (Ji 2008) [36] and Michigan Lake (Chen 2004) [37].

3. Large and Yeager (2004)-LY04 [21]:

   (9)$C_d\times {10}^3={\displaystyle \frac{2.7}{U_{10}}}+0.142+0.076U_{10}$

   Large and Yeager (2004) [21] used a nonlinear formula to describe the relationship between the $C_d$ value and wind speed. Later, this equation was chosen by the Community Climate System Model to represent the momentum transfer rate between the atmosphere and the free surface.

4. Andreas (2012)-A12 [22]:

   (10)$C_d\times {10}^3={\left({\displaystyle \frac{0.239+0.0433\{(U_{10}-8.271)+{\left(0.12{(U_{10}-8.271)}^2+0.181\right)}^{0.5}\}}{U_{10}}}\right)}^2$

   Andreas (2012) [22] used the friction velocity coefficient versus the neutral-stability wind speed of 10 m and sea roughness to test the approach given by Foreman and Emeis (1993) [38]. This method has been widely used in the study of ocean circulation and lake modeling (Vieira 2020 [39], Worsnop 2017 [40]).

5. Gao (2020)-G20 [26]:

   (11)$\begin{array}{c}\hfill \begin{array}{c}\hfill C_d\times {10}^3=\left\{\begin{array}{c}{\displaystyle \frac{C_d(U_{10}=5\mathrm{m}/\mathrm{s})}{0.00086}}{(0.0283+{\displaystyle \frac{0.00513}{U_{10}}})}^2,0\mathrm{m}/\mathrm{s}\le U_{10}\le 5\mathrm{m}/\mathrm{s}\hfill \\ -0.65+0.384l_n\left({U_{10}}^{1.537}{F}^{0.177}{d}^{-0.026}\right),5\mathrm{m}/\mathrm{s}\le U_{10}\le 25\mathrm{m}/\mathrm{s}\hfill \end{array}\right.\end{array}\end{array}$

   Gao (2020) [26] parameterized $C_d$ by fitting the experimental data from a wind tunnel and the observed data from shallow lakes considering the influence of wind blowing fetch and water depth. The purpose of this parameterization method is to study wind-induced currents in shallow lakes. However, it has not been used in numerical simulations yet, and therefore analyzing its applicability in this context is necessary.

Figure 1 shows the $C_d$ value as a function of wind speed with different parameterization methods. It was found that at low wind speeds ($U_{10}\le$ 5 m/s), the differences between the $C_d$ values are greater than at wind speeds larger than 5 m/s. The lower the wind speed, the higher the uncertainties. When the wind speed was constant, the $C_d$ value calculated by the G20 method showed a positive correlation with the wind blowing fetch and a negative correlation with the water depth. The specific form of the G20 parameterization method was determined by the wind blowing fetch and water depth of each grid point.

3. Numerical Simulation

3.1. Study Area

Lake Taihu (30°55^′40^″~31°32^′58^″ N, 119°52^′32^′~120°36^′10^″ E), located in the core area of the Yangtze River Delta, is a typical limited-fetch shallow lake. The wind-induced hydrodynamic process plays an important role in its environment and ecosystem. However, the complex topography and shoreline of the lake present great challenges to the simulation of wind-induced current modeling (Hu 2006) [41]. Lake Taihu covers an area of 2338 ${\mathrm{km}}^2$ and supplies the drinking water for the surrounding cities; it also has key roles for fisheries, aquaculture, and tourism. The average length of the lake is 68.5 ${\mathrm{km}}^2$ in the south–north direction and 34.0 ${\mathrm{km}}^2$ in the east–west direction. The lake shore is smooth in the southwest and tortuous in the northeast. The prevailing wind directions of Lake Taihu are southeast in summer and north in winter. A total of 172 tributaries connect to the lake. Most rivers with small discharges have no obvious impact on the flow field. For convenience, only Wangyu River, Tiaoxi River, and Taipu River were considered in this numerical simulation (Li 2011) [3].

3.2. Model Setup

The domain of our model encompasses Taihu Lake, including Wangyu River, Tiaoxi River, and Taipu River. Tributary inflow/outflow was treated in terms of source and sink. Numerical simulations were carried out from 1 to 31 August 2019. As shown in Figure 2, the model results were evaluated against the water elevations measured at the Wulou (WL), Jiapu (JP), Dapukou (DPK), Wangting (WT), Xukou (XK), and Dongtingxishan (DTXS) stations, and the velocity was measured at the V1 and V2 stations. The surface velocity was measured using two Acoustic Doppler Current Profilers (ADCPs). The initial conditions for the numerical model were a constant lake surface height of 3.55 m (the average value for the whole region) and zero velocity fields. The main driving forces for the numerical model were surface wind stress and atmospheric forcing. The wind data used to compute wind stress were acquired from 10 automatic weather stations (Taihu Basin Hydrological Information Service System) and 4 observation stations possessing sensors for detecting wind speed and direction (W1 to W4). The wind field was computed using the inverse distance weighting method with the observed data. The precipitation data were acquired by averaging the data from the Dushan and Dongtingxishan stations. Three sets of unstructured grids were set for the model (800 m × 800 m, 400 m × 400 m, and 100 m × 100 m) to analyze the effect of grid resolution sensitivity on the numerical results for the five $C_d$ parameterization methods. The grids used 20 layers in vertically stretched terrain following the sigma coordinate to more accurately simulate the vertical structure and changes of water bodies, as well as the variations of water flow and sediment at different depths, thereby improving the accuracy and reliability of the simulation. The specific grid parameters are shown in Table 1. The average maximum slopes of water depth for both grids were less than 0.33, avoiding pressure gradient errors from the sigma transformation (Mellor 1994) [42]. The global time step was set to 1 s in order to satisfy the Courant–Friedrichs–Lewy stability criterion (Mao 2008) [4]. The advection and diffusion of the three-dimensional variables were computed using the multidimensional upwind residual distribution scheme (MURD scheme). A calibration with respect to the Manning coefficients $\Delta m$ was carried out against the observed data (as shown in Section 4.1).

The effectiveness of the numerical result is evaluated using the $SS$ coefficient (Murphy 1992) [43] given by

(12)$SS=1-{\displaystyle \frac{{\sum }_{i=1}^N{(mod\left(i\right)-obs\left(i\right))}^2}{{\sum }_{i=1}^N{(obs\left(i\right)-ob{s}_{mean})}^2}}$

where $mod\left(i\right)$ are the numerical results; $obs\left(i\right)$ are the observed values; $ob{s}_{mean}$ are the mean of observed values. The performance can be classified as excellent ($SS$ > 0.65), very good (0.5 < $SS$ < 0.65), good (0.2 < $SS$ < 0.5), and poor ($SS$ < 0.2), respectively.

4. Results

4.1. Model Calibration

For calibration, the numerical model underwent a 4-day sequence (from 5 August to 9 August 2019) of atmospheric forcing and surface wind stress. Winds were predominantly southeasterly. The calibration of the Manning coefficients $\Delta m$ was performed with respect to the $SS$ performance computed from the numerical model’s results and the data observed at the V1 station (Chen and Liu 2015) [44]. The results are shown in Figure 3. The Manning coefficients are assumed to follow a uniform distribution between 0.012 and 0.048, with values sampled randomly using the Monte Carlo method. The calibration analysis shows that when $\Delta m$ is 0.024, the $SS$ coefficient of surface velocity reaches the maximum value of 0.7412; therefore, the Manning coefficient used in the numerical simulation was 0.024 for all calculated regions.

The results of the current structure calibration are shown in Figure 4. The observed data are given by Wu (2018) [45]. The development of wind-driven currents can cause two current structures when transporting water from the southeast portions of the lake to the northwest. Although the model current structures do not agree with the observed directions in some stations, the trend shows good agreement with the observed data. The numerical model can replicate the main feature of wind-induced circulation well.

4.2. Effects of $C_d$ Parameterization Methods

Figure 5 and Figure 6 show the numerical model results of the time series of water levels and velocity magnitudes produced from the five wind drag coefficient parameterization methods compared with the observed values at six water level stations and two velocity stations. The effect of the wind stress parameterization methods on the wind-induced current model is evident. The G20 shows a better performance than other four parameterization methods with the JP, DTXS, WL, and XK stations. The F76 parameterization method shows the worst performance among the six stations. For the simulation of velocity magnitude, the results of the G20 parameterization method are the closest to the observed data of the V1 and V2 velocity magnitude stations. The results of the LP81, LY04, and A12 methods are similar, while the results of the F76 parameterization method are worse at two stations. In the case of low wind speeds (<5 m/s), the velocity calculated by the LY04 parameterization method is obviously higher than the observed data due to its unrealistically high $C_d$ value.

Table 2 and Table 3 show the $SS$ coefficients of different cases between the numerical simulations and observed data. The data of the six water level stations and two velocity stations are averaged. When the wind speed is less than about 5 m/s, the LY04, A12, and G20 $C_d$ values are negatively correlated with wind speed; therefore, the velocity and water level data are divided into two different classes. For the water level simulation, when the wind speed is less than 5 m/s, the $SS$ coefficients for the five $C_d$ parameterization methods are 0.6339, 0.7072, 0.7121, 0.7373, and 0.7908. The closest coefficient to 1 is obtained by the G20 parameterization method, while the farthest is obtained by the F76 method. The results are similar when the wind speed is larger than 5 m/s. For the velocity simulation, the results are similar to those of the water level simulation; the closest $SS$ coefficient is produced by G20 parameterization ($SS$ = 0.7940) and the farthest is produced by F76 parameterization ($SS$ = 0.6894). From the $SS$ coefficient analysis, it can be found that the G20 parameterization method shows the most satisfactory performance in both the water level and velocity simulations while the F76 parameterization method shows worse performance, and the other three parameterization methods produce similar simulation results; the A12 parameterization method shows the second best performance. Taken together, G20 parameterization is more suitable for wind-induced current simulation in shallow lakes compared to the other methods.

Figure 7a shows the wind stress calculated using the F76, LP81, LY04, A12, and G20 parameterization methods. It can be seen from the figure that the surface velocity at the V1 point is positively correlated with local wind speed and wind stress. A higher wind speed leads to higher wind stress and higher surface velocity. As is shown in Figure 7b, the wind stress growth rate calculated by LY04 is much higher than the other formulas with a decrease in wind speed in the low wind speed range. LY04 also shows a much higher surface velocity than that calculated by the other methods. The wind stress calculated by the F76 parameterization method is the lowest under low-wind-speed conditions, and the results of the G20, A12, and LP81 methods are similar. According to the calculation results of the two velocity stations (see Table 3), LY04 parameterization overestimates wind stress in the low wind speed range, while F76 parameterization underestimates it. A12 and G20 can simulate wind stress more effectively, and G20 parameterization has the highest accuracy.

4.3. Effects of Grid Resolution

Several experiments were conducted to analyze the influence of the grid resolution on wind-induced current in shallow lakes. Here, all calculation settings were kept the same except for the grid resolution. The wind stresses were calculated using the F76, LP81, LY04, A12, and G20 parameterization methods. Figure 8 compares the velocity magnitude computed using the 800 m, 400 m, and 100 m grid resolutions. Figure 8a shows the simulated velocity magnitudes under different grid resolutions at the V1 measurement station using the LP81 parameterization method. The grid resolution has a great influence on the velocity magnitude at high wind speeds, but not obviously at low wind speeds. Figure 8b shows the simulated velocity magnitudes under different grid resolutions at the V1 measurement station using the G20 parameterization, With an increase in grid resolution, the velocity magnitude gradually decreases, and the decreasing amplitude at high wind speeds is significantly higher than that at low wind speeds.

Table 4 shows the average $SS$ values of the five $C_d$ parameterization methods. From the $SS$ analysis, it can be found that with an increase in grid resolution, the SS coefficients calculated by F76 and LP81 gradually decrease and the SS coefficients calculated by LY04 and A12 gradually increase. For the above four $C_d$ parameterization methods, the $SS$ coefficients do not evidently change with a change in grid resolution. When the grid resolution increases from 800 m to 100 m, the SS coefficient of G20 rises from 0.5138 to 0.7934. For the coarse grid (800 m), the G20 parameterization method shows an undesired performance ($SS$ = 0.5138), while the SS coefficients calculated by the other four $C_d$ parameterization methods are all above 0.67, with A12 showing the best performance. For the medium grid (400 m), although the accuracy of the G20 is still worse than that other parameterization methods, the difference in the SS coefficients between G20 and the other four $C_d$ parameterization methods is very small. For the fine grid (100 m), the results are the opposite, with G20 parameterization showing the optimum result. It was found that an increase in grid resolution can improve the accuracy of the numerical simulation, and the change in grid resolution has little effect on the accuracy of the other four $C_d$ parameterization methods.

In the numerical model, the velocity magnitude is computed as the mean of the grid boxes. In a grid with high resolution, the box is small, which means that the velocity magnitude is averaged over a small area. In particular, in complex shallow lakes like Lake Taihu, the grid resolution effect can be important because high-resolution grids can better resolve the bathymetry, geography, and topography of the lake. The higher the grid resolution, the closer the bathymetry, geography, and topography obtained from the numerical simulation are to the real values, and the more accurate the flow field results. For the G20 parameterization method, the resolution effect is more important; the $C_d$ obtained from this method is not only related to the wind speed, but also to the wind blowing fetch and water depth at each grid node. Therefore, when the grid resolution is low, the spatial distribution of wind blowing fetch and local water depth is quite different from the real situation, which leads to a significant decrease in the accuracy of the numerical simulation’s results. It was found that the combination of a high-resolution grid and the G20 parameterization method can greatly improve the prediction of wind-induced currents in shallow lakes, while the F76, LP81, LY04, and A12 parameterization methods are not sensitive to changes in grid resolution.

4.4. Computational Efficiency

In addition, Table 5 shows the computing time taken for the different $C_d$ parameterization methods and grid resolutions. It is found that the computing time for the F76, LP81, LY04, and A12 calculations is similar, which is much lower than that of the G20 calculation. Table 5 shows that the higher the grid resolution, the higher the computational accuracy of G20. The calculated water depth and fetch length are closer to the real situation of Lake Taihu in a higher resolution, the value $C_d$ of G20 is thereby more reasonable compared to the values calculated by other four parameterization methods. When G20 parameterization is used, the wind blowing fetch on each grid node must be obtained using the wind direction before each iteration starts. As the local wind stress is calculated by combining water depth and wind speed, the G20 parameterization method takes much longer to use than the other four parameterization methods. According to Table 4 and Figure 8, when the grid is coarse, the A12 calculation time is short and its accuracy is the highest among the five parameterization methods. When the grid resolution is high, the calculation time of the G20 parameterization method is much higher than the other four parameterization methods, but its accuracy is the highest. In terms of computational efficiency, the combination of the A12 formula and a coarse grid can meet the requirements of achieving accurate simulations wind-induced current in fetch-limited shallow lakes. However, when the calculation conditions permit, the combination of the G20 formula and a high-resolution grid can significantly improve the calculation accuracy.

5. Conclusions

Wind stress is the main driving force that determines the performance of wind-induced current modeling in fetch-limited shallow lakes. $C_d$ parameterization, as a physical parameter used to characterize wind stress over a free surface, has a great influence on the results of numerical simulations. This study aims to analyze the influence of wind stress drag coefficient parameterization methods and grid resolutions on wind-induced currents in Lake Taihu.

A series of numerical simulations were conducted for three grid resolutions (800 m × 800 m, 400 m × 400 m, and 100 m × 100 m) using five $C_d$ parameterization methods. The result showed that the choice of different $C_d$ parameterization methods can lead to significant differences in both water level and velocity values. G20 showed a better performance in both water level and velocity simulations, while F76 showed a worse performance than other parameterization methods. In the experiments using different grid resolutions, it was found that when the grid was coarse (800 m × 800 m), the accuracy of the G20 parameterization method was at the bottom of the five $C_d$ methods. The combination of a high-resolution grid and the G20 parameterization method can greatly improve the prediction of wind-induced currents in shallow lakes, while the F76, LP81, LY04, and A12 parameterization methods are not sensitive to changes in grid resolution.

Apart from G20, other four parameterization methods ignore the influence of blowing wind fetch length and water depth on the wind drag coefficient, $C_d$. Distinguished from that in the marine environment, $C_d$ is incompletely developed in shallow lakes. G20 fully considers the impacts of limited-fetch on undeveloped $C_d$ and thus has a higher computation accuracy than other formulas in Lake Taihu. However, we found that in the case of G20 formula, the fetch of each grid point in Lake Taihu needs to be computed before calculating $C_d$ at each time step, which significantly increases the amount of calculation. Therefore, for each grid resolution, the calculation time of G20 was the longest. Considering the accuracy and efficiency of calculation, the combination of the A12 formula and a coarse grid can meet the requirements for accurately simulating wind-induced currents in fetch-limited shallow lakes. We recommend the A12 parameterization method when computing resources are limited and grid resolution is low. However, when the calculation conditions permit, the combination of the G20 formula and a high-resolution grid can significantly improve the calculation accuracy.

Author Contributions

Conceptualization, C.Z. and C.W.; methodology, A.G.; validation, C.Z.; formal analysis, C.Z.; investigation, C.Z.; resources, Y.D. and G.Y.; data curation, Y.D. and J.Y.; writing—original draft preparation, C.W.; writing—review and editing, C.Z. and J.L.; visualization, C.Z.; supervision, Z.Z.; project administration, C.Z.; funding acquisition, C.Z. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

Data are contained within the article.

Conflicts of Interest

Author Jiaqing Liao was employed by the company Hangzhou Dingchuan Information Technology Co., Ltd. The remaining authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. Wind stress drag coefficients [Forumla omitted. See PDF.] as a function of wind speed with different parameterizations. F76 (Equation (7), black), LP81 (Equation (8), red), LY04 (Equation (9), green), A12 (Equation (10), blue) and magenta lines with different forms of symbols represent G20 formula under different wind blowing fetches and water depths (Equation (11)).

Figure 2. The location of the study area, surrounding rivers, observational stations in and around Lake Taihu. Black circles represent the location of water level stations (WL, Wulou station; JP, Jiapu station; DPK, Dapukou station; WTTH, Wangting station; XK, Xukou station; DTXS, Dongtingxishan station). Black triangles represent the location of velocity stations (V1, ADCP1 station; V2, ADCP2 station). Black squares represent the location of wind observed stations (from W1 to W4).

Figure 3. Calibration of manning coefficients based on comparison of the numerical model results and observed data at V1 station. The shaded area represents the envelope of results for manning coefficients from [Forumla omitted. See PDF.] = 0.012 (upper line) to [Forumla omitted. See PDF.] = 0.048 (lower line); best model result ([Forumla omitted. See PDF.] = 0.024, the maximum [Forumla omitted. See PDF.] value of 0.7412) is plotted as blue line, the red dot represents the observed data measured with a one-hour interval.

Figure 4. (a) Surface layer current and (b) bottom layer current structures colored with velocity distribution and observed data (Wu 2018 [45]) in Lake Taihu, where black arrow denotes the modeled current structure, red arrow denotes the observed flow field.

Figure 5. Numerical model results of time series of water levels produced from F76, LP81, LY04, A12, and G20 parameterization at (a) JP, (b) DTXS, (c) WL, (d) XK, (e) WTTH, (f) DPK measurement station.

Figure 6. Numerical model results of time series of surface velocity magnitude produced from F76, LP81, LY04, A12, and G20 parameterization at (a) V1 and (b) V2 different stations.

Figure 7. (a) Wind speed [Forumla omitted. See PDF.] and (b) wind stress [Forumla omitted. See PDF.] produced from F76, LP81, LY04, A12, and G20 parameterization at V1 station during simulated periods from 9–15 August 2019.

Figure 8. Comparisons of velocity magnitude computed using the 800 m, 400 m, and 100 m grid resolution and different [Forumla omitted. See PDF.] parameterizations of (a) LP81 and (b) G20, respectively.

Table 1

Grid parameters with different grid resolutions.

Grid Resolution	Nodes	Elements	Vertical Layers
800 m × 800 m	8222	4367	20
400 m × 400 m	29,472	15,163	20
100 m × 100 m	105,083	53,681	20

Table 2

$SS$ coefficients between numerical simulations and observed data. SS coefficients are obtained by averaging the results from 6 water level measurement stations (XMK, DTXS, WL, XK, WTTH, and JP).

Indicator	F76	LP81	LY04	A12	G20
SS ($U_{10}\le$ 5 m/s)	0.6339	0.7072	0.7121	0.7373	0.7908
SS ($U_{10}\ge$ 5 m/s)	0.4107	0.5692	0.5809	0.5750	0.6520
SS (Total)	0.5592	0.6697	0.6950	0.7050	0.7623

Table 3

The numerical model results of time series of velocity magnitude produced from different wind drag coefficients formulae compared with observed data at 2 velocity stations.

Indicator	F76	LP81	LY04	A12	G20
SS ($U_{10}\le$ 5 m/s)	0.7276	0.7862	0.7539	0.7890	0.8333
SS ($U_{10}\ge$ 5 m/s)	0.6180	0.6498	0.6738	0.6760	0.7897
SS (Total)	0.6894	0.7191	0.7241	0.7468	0.7940

Table 4

$SS$ coefficients between numerical simulations and observations. Five parameterizations and different grid resolutions are compared. SS coefficients are obtained by averaging the results from 2 velocity stations (V1 and V2).

${\mathit{C}}_{\mathit{d}}$ Parameterization	Indicator	100 m × 100 m	400 m × 400 m	800 m × 800 m
F76	SS coefficient	0.6894	0.6997	0.6788
LP81	SS coefficient	0.7191	0.7204	0.7021
LY04	SS coefficient	0.7241	0.7225	0.7001
A12	SS coefficient	0.7468	0.7333	0.7055
G20	SS coefficient	0.7934	0.6824	0.5138

Table 5

Comparisons of computing time using different $C_d$ parameterizations and different grid resolutions.

${\mathit{C}}_{\mathit{d}}$ Parameterization	100 m × 100 m	400 m × 400 m	800 m × 800 m
F76	11.51 h	20.51 h	57.67 h
LP81	11.75 h	20.21 h	57.92 h
LY04	12.54 h	21.26 h	60.35 h
A12	12.95 h	21.84 h	60.33 h
G20	52.87 h	87.50	143.5 h

References

1. Tsanis, I.K. Simulation of wind-induced water currents. J. Hydraul. Eng.; 1989; 115, pp. 1113-1134. [DOI: https://dx.doi.org/10.1061/(ASCE)0733-9429(1989)115:8(1113)]

2. Wu, J. Wind-induced drift currents. J. Fluid Mech.; 2006; 68, pp. 49-70. [DOI: https://dx.doi.org/10.1017/S0022112075000687]

3. Li, Y.; Acharya, K.; Yu, Z. Modeling impacts of Yangtze River water transfer on water ages in Lake Taihu, China. Ecol. Eng.; 2011; 37, pp. 325-334. [DOI: https://dx.doi.org/10.1016/j.ecoleng.2010.11.024]

4. Jingqiao, M.A.O.; Chen, Q.; Chen, Y. Three-dimensional eutrophication model and application to Taihu Lake, China. J. Environ. Sci.; 2008; 20, pp. 278-284.

5. Li, Y.; Jalil, A.; Du, W.; Gao, X.; Wang, J.; Luo, L.; Li, H.; Dai, S.; Hashim, S.; Yu, Z. et al. Wind induced reverse flow and vertical profile characteristics in a semi-enclosed bay of large shallow Lake Taihu, China. Ecol. Eng.; 2017; 102, pp. 224-233. [DOI: https://dx.doi.org/10.1016/j.ecoleng.2017.02.022]

6. Chu, D.; Zhang, J.; Wu, Y.; Jiao, X.; Qian, S. Sensitivities of modelling storm surge to bottom friction, wind drag coefficient, and meteorological product in the East China Sea. Estuar. Coast. Shelf Sci.; 2019; 231, 106460. [DOI: https://dx.doi.org/10.1016/j.ecss.2019.106460]

7. Kranenburg, C. Wind-driven chaotic advection in a shallow model lake. J. Hydraul. Res.; 2010; 30, pp. 29-46. [DOI: https://dx.doi.org/10.1080/00221689209498945]

8. Li, Y.P.; Tang, C.Y.; Yu, Z.B.; Acharya, K. Uncertainty and sensitivity analysis of large shallow lake hydrodynamic models. Adv. Water Sci.; 2012; 33, pp. 271-277.

9. Wróbel-Niedźwiecka, I.; Drozdowska, V.; Piskozub, J. Effect of drag coefficient formula choice on wind stress climatology in the North Atlantic and the European Arctic. Oceanologia; 2019; 61, pp. 291-299. [DOI: https://dx.doi.org/10.1016/j.oceano.2019.02.002]

10. Zhang, C.; Chen, L. A review of wind-driven hydrodynamics in large shallow lakes: Importance, process-based modeling and perspectives. Camb. Prism. Water; 2023; 1, e16. [DOI: https://dx.doi.org/10.1017/wat.2023.14]

11. Morales-Marín, L.A.; French, J.R.; Burningham, H. Implementation of a 3D ocean model to understand upland lake wind-driven circulation. Environ. Fluid Mech.; 2017; 17, pp. 1255-1278. [DOI: https://dx.doi.org/10.1007/s10652-017-9548-6] [PubMed: https://www.ncbi.nlm.nih.gov/pubmed/32009858]

12. Xu, T.; Zhang, M.; Jiang, H.; Jiang, H.; Tang, J.; Zhang, H.; Qiao, H. Numerical investigation of the effects of aquatic plants on wind-induced currents in Taihu Lake in China. J. Hydrodyn.; 2019; 31, pp. 778-787. [DOI: https://dx.doi.org/10.1007/s42241-018-0091-9]

13. Oliveira, H.; Fernandes, E.; Möller, O., Jr.; García-Rodríguez, F. Relationships between Wind Effect, Hydrodynamics and Water Level in the World’s Largest Coastal Lagoonal System. Water; 2019; 11, 2209. [DOI: https://dx.doi.org/10.3390/w11112209]

14. Sheppard, P.A.; Tribble, D.T.; Garratt, J.R. Studies of turbulence in the surface layer over water (Lough Neagh). Part I. Instrumentation, programme, profiles. Q. J. R. Meteorol. Soc.; 1972; 98, pp. 627-641.

15. Flather, R.A. Results from a Storm Surge Prediction Model of the North-West European Continental Shelf for April, November and December 1973; Report No.24 Institute of Oceanographics Sciences: Surrey, UK, 1976.

16. Large, W.G.; Pond, S. Open Ocean Momentum Flux Measurements in moderate to strong winds. J. Phys. Oceanogr.; 1981; 11, pp. 324-336. [DOI: https://dx.doi.org/10.1175/1520-0485(1981)011<0324:OOMFMI>2.0.CO;2]

17. Geernaert, G.L.; Larsen, S.E.; Hansen, F. Measurements of the wind stress, heat flux, and turbulence intensity during storm conditions over the North Sea. J. Geophys. Res. Oceans; 1987; 92, pp. 13127-13139. [DOI: https://dx.doi.org/10.1029/JC092iC12p13127]

18. Ataktürk, S.S.; Katsaros, K.B. Wind Stress and Surface Waves Observed on Lake Washington. J. Phys. Oceanogr.; 1999; 29, pp. 633-650. [DOI: https://dx.doi.org/10.1175/1520-0485(1999)029<0633:WSASWO>2.0.CO;2]

19. Guan, C.; Xie, L. On the Linear Parameterization of Drag Coefficient over Sea Surface. J. Phys. Oceanogr.; 2004; 34, pp. 2847-2851. [DOI: https://dx.doi.org/10.1175/JPO2664.1]

20. Yelland, M.J.; Taylor, P.K. Wind Stress Measurements from the Open Ocean. J. Phys. Oceanogr.; 1996; 26, pp. 541-558. [DOI: https://dx.doi.org/10.1175/1520-0485(1996)026<0541:WSMFTO>2.0.CO;2]

21. Large, W.G.; Yeager, S.G. Diurnal to Decadal Global Forcing For Ocean and Sea-Ice Models: The Data Sets and Flux Climatologies; National Center for Atmospheric Research: Boulder, CO, USA, 2004.

22. Andreas, E.L.; Mahrt, L.; Vickers, D. A New Drag Relation for Aerodynamically Rough Flow over the Ocean. J. Atmos. Sci.; 2012; 69, pp. 2520-2537. [DOI: https://dx.doi.org/10.1175/JAS-D-11-0312.1]

23. Janssen, P.A.; Peter, A.E.M. Quasi-Linear Theory of Wind-Wave Generation Applied to Wave Forecasting. J. Phys. Oceanogr.; 1991; 21, pp. 1631-1642. [DOI: https://dx.doi.org/10.1175/1520-0485(1991)021<1631:QLTOWW>2.0.CO;2]

24. Tolman, H.L.; Chalikov, D. Source Terms in a Third-Generation Wind Wave Model. J. Phys. Oceanogr.; 1996; 26, pp. 2497-2518. [DOI: https://dx.doi.org/10.1175/1520-0485(1996)026<2497:STIATG>2.0.CO;2]

25. Mueller, J.A.; Veron, F. Nonlinear Formulation of the Bulk Surface Stress over Breaking Waves: Feedback Mechanisms from Air-flow Separation. Bound.-Layer Meteorol.; 2008; 130, pp. 117-134. [DOI: https://dx.doi.org/10.1007/s10546-008-9334-6]

26. Gao, A. Study on the wind drag coefficient of shallow lakes considering wind-wave-flow characteristics. Ph.D. Thesis; Nanjing Hydraulic Research Institute: Nanjing, China, 2020.

27. Mastenbroek, C.; Burgers, G.; Janssen, P.A.E.M. The Dynamical Coupling of a Wave Model and a Storm Surge Model through the Atmospheric Boundary Layer. J. Phys. Oceanogr.; 1993; 23, pp. 1856-1866. [DOI: https://dx.doi.org/10.1175/1520-0485(1993)023<1856:TDCOAW>2.0.CO;2]

28. Smith, S.D.; Anderson, R.J.; Oost, W.A.; Kraan, C.; Maat, N.; De Cosmo, J.; Katsaros, K.B.; Davidson, K.L.; Bumke, K.; Hasse, L. et al. Sea surface wind stress and drag coefficients: The hexos results. Bound.-Layer Meteorol.; 1992; 60, pp. 109-142. [DOI: https://dx.doi.org/10.1007/BF00122064]

29. Ebuchi, N.; Toba, Y.; Kawamura, H. On the wave dependence of sea-surface wind stress. J. Oceanogr. Soc. Jpn.; 1990; 46, pp. 177-183. [DOI: https://dx.doi.org/10.1007/BF02125578]

30. Gao, A.; Wu, S.; Chen, L.; Liu, S.; Xu, Z.; Zhao, Y.; Wang, F.; Zhu, S. Experiemental study on the wave height distribution of wind-induced waves in the growth stage under finite water depth. Water Supply; 2019; 19, pp. 2079-2270. [DOI: https://dx.doi.org/10.2166/ws.2019.108]

31. Elfouhaily, T.; Chapron, B.; Katsaros, K.; Vandemark, D. A unified directional spectrum for long and short wind-driven waves. J. Geophys. Res. Ocean.; 1997; 102, pp. 15781-15796. [DOI: https://dx.doi.org/10.1029/97JC00467]

32. Fang, P.; Zhao, B.; Zeng, Z.; Yu, H.; Lei, X.; Tan, J. Effects of Wind Direction on Variations in Friction Velocity with Wind Speed Under Conditions of Strong Onshore Wind. J. Geophys.-Res.; 2018; 123, pp. 7340-7353. [DOI: https://dx.doi.org/10.1029/2017JD028010]

33. Moon, I.J.; Hara, T.; Ginis, I.; Belcher, S.E.; Tolman, H.L. Effect of Surface Waves on Air–Sea Momentum Exchange. Part I: Effect of Mature and Growing Seas. J. Atmos. Sci.; 2004; 61, pp. 2321-2333. [DOI: https://dx.doi.org/10.1175/1520-0469(2004)061<2321:EOSWOA>2.0.CO;2]

34. Costi, J.; Marques, W.C.; de Paula Kirinus, E.; de Freitas Duarte, R.; Arigony-Neto, J. Water level variability of the Mirim—São Gonçalo system, a large, subtropical, semi-enclosed coastal complex. Adv. Water Resour.; 2018; 117, pp. 75-86. [DOI: https://dx.doi.org/10.1016/j.advwatres.2018.05.008]

35. Kopmann, R.; Markofsky, M. Three-dimensional water quality modelling with TELEMAC-3D. Hydrol. Process.; 2000; 14, pp. 2279-2292. [DOI: https://dx.doi.org/10.1002/1099-1085(200009)14:13<2279::AID-HYP28>3.0.CO;2-7]

36. Ji, Z.-G. Hydrodynamics and Water Quality: Modeling Rivers, Lakes, and Estuaries; John Wiley: Hoboken, NJ, USA, 2008; Volume 89, 676.

37. Chen, C.; Wang, L.; Ji, R.; Budd, J.W.; Schwab, D.J.; Beletsky, D.; Fahnenstiel, G.L.; Verploeg, H.; Eadie, B.; Cotner, J. Impacts of suspended sediment on the ecosystem in Lake Michigan: A comparison between the 1998 and 1999 plume events. J. Geophys. Res.; 2004; 109, pp. 1-18. [DOI: https://dx.doi.org/10.1029/2002JC001687]

38. Foreman, R.J.; Emeis, S. Revisiting the Definition of the Drag Coefficient in the Marine Atmospheric Boundary Layer. J. Phys. Oceanogr.; 1993; 40, pp. 2325-2332. [DOI: https://dx.doi.org/10.1175/2010JPO4420.1]

39. Vieira, V.M.; Mateus, M.; Canelas, R.; Leitão, F. The FuGas 2.5 Updated for the Effects of Surface Turbulence on the Transfer Velocity of Gases at the Atmosphere–Ocean Interface. J. Mar. Sci. Eng.; 1989; 8, 435. [DOI: https://dx.doi.org/10.3390/jmse8060435]

40. Worsnop, R.P.; Bryan, G.H.; Lundquist, J.K.; Zhang, J.A. Using Large-Eddy Simulations to Define Spectral and Coherence Characteristics of the Hurricane Boundary Layer for Wind-Energy Applications. Bound.-Layer Meteorol.; 2017; 165, pp. 55-86. [DOI: https://dx.doi.org/10.1007/s10546-017-0266-x]

41. Hu, W.; Jørgensen, S.E.; Zhang, F. A vertical-compressed three-dimensional ecological model in Lake Taihu, China. Ecol. Model.; 2006; 190, pp. 367-398. [DOI: https://dx.doi.org/10.1016/j.ecolmodel.2005.02.024]

42. Mellor, G.L.; Ezer, T.; Oey, L.Y. The Pressure Gradient Conundrum of Sigma Coordinate Ocean Models. J. Atmos. Ocean. Technol.; 1994; 11, pp. 1126-1134. [DOI: https://dx.doi.org/10.1175/1520-0426(1994)011<1126:TPGCOS>2.0.CO;2]

43. Murphy, A.H. Climatology, Persistence, and Their Linear Combination as Standards of Reference in Skill Scores. Weather. Forecast.; 1992; 7, pp. 692-698. [DOI: https://dx.doi.org/10.1175/1520-0434(1992)007<0692:CPATLC>2.0.CO;2]

44. Chen, Y.; Liu, Q. Numerical study of hydrodynamic process in Chaohu Lake. J. Hydrodyn. Ser. B; 2015; 27, pp. 720-729. [DOI: https://dx.doi.org/10.1016/S1001-6058(15)60534-0]

45. Wu, T.; Qin, B.; Ding, W.; Zhu, G.; Zhang, Y.; Gao, G.; Xu, H.; Li, W.; Dong, B.; Luo, L. Field Observation of Different Wind-Induced Basin-Scale Current Field Dynamics in a Large, Polymictic, Eutrophic Lake. J. Geophys. Res. Ocean.; 2018; 123, pp. 6945-6961. [DOI: https://dx.doi.org/10.1029/2018JC014040]