1. Introduction
Large shallow lakes, especially those located in highly developed areas, provide multifunctional services for industry, agriculture, navigation, and recreation; unfortunately, they have often suffered from severe eutrophication problems [
1,
2,
3] In general, population density around these lakes is high [
2,
4,
5], which leads to an additional high waste water and associated nutrients load [
6]. These issues have triggered increasing attention for the restoration of the water quality and ecological status of large shallow lakes. The shallowness of lakes is usually emphasized, since the total volume is usually small, and shallow lakes are more sensitive to the effects of wind, evaporation, and human interference, compared to the relatively deep lakes [
7,
8,
9].
There are two common approaches for shallow lake restoration nowadays. One rather effective approach is to control the source, and thus, to decrease the total nutrient load [
10,
11], while the other approach is to increase the hydrodynamic circulation [
12]. Especially in densely populated areas, source control almost reaches the limit of present technology, whereas enhancing the hydrodynamic circulation might offer an important contribution to improve the water quality of shallow lakes. However, comparing research on hydrodynamics in oceans, coastal zones, rivers, and deep lakes [
13], only limited attention has been paid to the hydrodynamics of large shallow lakes. In fact, thorough qualitative studies on horizontal circulation patterns of large shallow lakes are rather scarce, and quantitative studies are even more surprisingly rare. Especially due to their shallowness, the dominant hydrodynamic processes and the corresponding ecosystem in shallow lakes differs very much from that in deeper lakes [
14,
15]. In deep lakes and reservoirs, due to the large depth, the temperature profile is determined by thermal stratification and mixing that dominates the hydrodynamics, especially the vertical exchange [
15], while in large shallow lakes, stratification is seldom observed, which results in substantially different processes.
In large shallow lakes, i.e., with a mean depth <3 m [
15], water quality and eutrophication problems are closely related to advection and diffusion processes driven predominantly by wind forcing [
16,
17,
18]. Momentum transferred by the wind via surface shear stresses generates waves, currents, and associated turbulence [
19,
20] While these processes are essentially three-dimensional, the shallowness allows for a depth-averaged, two dimensional representation for some process features [
21,
22,
23,
24]. Currents induced by wind forcing with velocities at the 10 cm/s scale can lead to lake-wide horizontal circulation patterns with the potential of creating vertical circulations, such as Langmuir circulations [
20]. Furthermore, momentum transferred to the bottom will stir up sediment and keep it suspended by turbulence [
18,
25]. During the suspension and resuspension of sediments, pollutions and nutrients attached to the sediment are released into water column and then transported and mixed by the large-scale circulation [
7,
26,
27]. Thus, the spatial and temporal large-scale, shallow lake horizontal circulation is essential for system understanding before we move to water quality and aquatic ecosystem issues.
In this paper, the focus is therefore on the spatiotemporal wind driven circulations in Taihu Lake, an unusual, extremely shallow and geometrically complex lake including bays and islands (
Figure 1). Enhanced anthropogenic emissions in recent years, have had a huge impact on water quality and strongly motivated the eutrophication [
28]. Quite some studies have been carried out to seek solutions for the water quality issue of Taihu Lake. One of the most famous engineering interventions is the water transfer project, which diverts water from the Yangtze River of better water quality but more suspended sediment though to the lake, to dilute the excess nutrients and pollutions in the lake water. However, whether the water transfer project has succeeded in improving the general water quality in Taihu Lake remains unclear, since a positive influence could only be observed in some parts of the lake [
29,
30]. These facts indicate that a better understanding of the hydrodynamics of Taihu Lake is urgently required for future water quality management and restoration of the lake ecosystem. In this paper, numerical models are used to study the hydrodynamics and water quality of Taihu Lake under steady and unsteady wind conditions. Even though, over 20 studies have been carried out using two or three-dimensional numerical models before this study, their focus is either on the resulting water quality index or on the ecological status, and none of them are dedicated to a thorough, quantitative description of (wind induced) large-scale hydrodynamic circulation itself, nor to the implication of hydrodynamic circulation to environment or ecology in this lake [
25,
31,
32].
In this research, hydrodynamic circulation in shallow lakes is defined as the spatially heterogeneous large-scale movement of water. Velocity vectors and particle tracers are used to indicate the hydrodynamic circulation patterns. Timescales are usually from days, weeks, to seasons, and spatial scales can be a few kilometers. Thus, barotropic seiches (~1 day), wind-driven short waves (~seconds), and other processes of smaller timescale are not included in this study.
The overall goal of this study is to gain a better understanding of the wind-induced hydrodynamics, and thereby to provide essential knowledge for the design and implementation of future lake restoration measurements, using state-of-the-art numerical models as a quantitative assessment tool.
Thus, our objectives in this work are: 1. To investigate the rich structure of spatial and temporal varying hydrodynamic circulation (i.e., direction, intensity and position) in a large shallow lake with complex geometry and irregular shape; 2. To quantify wind induced changes in hydrodynamic circulations (volume exchange between subbasins and vertical variations) on spatial scales; 3. To discuss implications of anthropological effects, such as large scale water transfer, on hydrodynamic circulations.
3. Model Description
3.1. Numerical Model Description
The open source three-dimensional shallow water numerical model Delft3D (Delft, The Netherlands), developed by Deltares, is used in this study. Delft3D consists of various modules, covering different physical processes, ranging from flow, sediment transport, morphology to water quality, aquatic ecology, and particle tracking, etc. The model has been extensively applied worldwide in the fields of hydrodynamics, sediment transport, morphology, and water quality in fluvial, lacustrine, estuarine, and coastal environments. In this study, the model is used to simulate the hydrodynamics of the lake for the purpose of illustrating the spatial and temporal large scale hydrodynamic circulation of Taihu Lake induced by the wind shear stress, and the discharge from the tributary rivers given the complex geometry and shallow bathymetry of the lake.
3.2. Model Setup
A Cartesian rectangular computational grid is used with grid resolution of 1000 m and total grid number of 9660 (
Figure 1). Five vertical sigma layers are defined uniformly for three-dimensional scenarios. For Taihu Lake, the deepest point of the lake appears in the lake center area at 2.66 m, while the shallow points are located close to shorelines, having a typical depth of 0.8 m.
The reference scenario simulates the hydrodynamics of the lake over the entire year of 2008. The simulation time step is based on the Courant–Friedrichs–Lewy (CFL) number and the grid size. To ensure model stability and accuracy, the time step is set to be 5 min. The output time step for observation points is 60 min.
3.3. Tracer Redistribution Calculation
Conservative tracers are used in this numerical study to simulate and illustrate the temporal and spatial scale of the hydrodynamic circulation throughout Taihu Lake induced by wind. Conservative tracers do not decay with time and space, and are only passively transported and mixed with hydrodynamics. In this study, when initially released, they are uniformly distributed in the water column. In Delft3D, the principle of simulating the tracer movement is to describe both the advection and the diffusion processes. The lateral term could be realized by giving an additional random movement other than the advective movement in each time step. Size and direction of the movement is proportional to the horizontal and vertical diffusion. Tracers were released in each subbasin in the model area at a constant rate of 5 kg/m
3 for duration of 3 months. The tracer model results will be further discussed in
Section 5.
3.4. Model Calibration
The numerical model calibration uses meteorological data and boundary river discharge data for the entire year of 2008. For bottom roughness, we use the Chezy formula, with the roughness coefficient assumed constant at 65 m
1/2/s for the entire lake. Horizontal viscosity is set to be 0.002 m
2/s based on literature [
35]. For vertical turbulent transport process, a standard k-ε turbulence model is used.
For model calibration, in situ measurements of water levels from five monitoring stations across the lakes (for locations see
Figure 1; data source Taihu Basin Authority) are compared to the simulated results in
Figure 3. These five monitoring stations are, namely, Wangting station, Dapukou station, Jiapu station, Xiaomeikou station, and Xishan station. Water levels were observed once a day at each station with occasional interruptions.
Figure 3 show that although the five monitoring stations are located in different parts of the lake, the water levels show very similar trends. Considering the water level difference due to wind setup over the largest fetch length in Taihu Lake with constant 5 m/s wind at ~5 cm scale, and given the comparison of tributary discharge and average water level trends (
Figure 2) with the mass balance check result (
Table 1), water levels in Taihu Lake are predominately modulated by tributary discharge, precipitation, and evaporation. The water levels had a slight increase until mid-February, and reached an annual peak value at the beginning of July, then gradually fell to 0 at the end of the year. The highest water level of the whole lake reaches around 1 m, while the lowest water level is around 0 m at the beginning of April. A summary of quantitative model performance indicators are listed in
Table 2.
Model performance indicators indicate that the simulation results and in situ measurement for each observation station differ only by a few centimeters, which demonstrates that the model well reflects the water level trend due to the influence of various sources, including tributary discharge, precipitation, and evaporation. The largest variance between modelled and measured water level is in Wangtingtai Station, potentially because of the irregular lake margin near the station. Further sensitivity analysis shows variation of other parameters, like bottoms roughness contributes little to the water level model results.
5. Discussion
Wind influence on the hydrodynamic circulations has been discussed above. The Delft3D simulation model is shown to be able to reproduce the spatial and temporal hydrodynamic circulation patterns in Taihu Lake. However, given the uncertainty in numerical and process input parameters and the unique geomorphology of Taihu Lake, these model results should be carefully scrutinized. To evaluate input errors, model sensitivity analyses is provided in
Section 5.1 and
Section 5.2, concerning grid size and bed roughness, respectively. Velocity vorticity, as the key indicator of hydrodynamic circulation, predominantly modulated by wind, depth, bathymetry gradient etc., is discussed through both a theoretical analysis and flat bottom model tests in
Section 5.3. Furthermore, Lagrangian-based tracer tests are used to evaluate emergency pollution/leakage effects in
Section 5.4, and water transfer effects in
Section 5.5.
5.1. Grid Size Effects
Before model calibration, sensitivity tests have been carried out to test the influence of numerical parameters. One of the major issues is grid size. In this study, the grid size is chosen to be 1000 m. A grid size of 500 m is utilized in most earlier numerical studies of Taihu Lake. To test the influence of grid size on the model accuracy and efficiency, a 500 m grid size model was set up and model behavior was analyzed (
Table 3). In both cases, the chosen time step of 10 min fulfils the courant number requirement.
Model performance for the water level of five monitoring stations with the finer grid of 500 m is listed in
Table 3. Results of the sensitivity analysis show the water level variations and hydrodynamic circulation patterns hardly change with grid sizes. Model performances of both models are equally good, while the coarse grid has an advantage in model efficiency. Thus, the 1000 m grid size was chosen.
5.2. Bed Roughness
The bottom roughness coefficient is another common tuning parameter in hydrodynamic models. To be noticed, in Delft3D-FLOW, bottom roughness coefficient with other formulas is converted to Chezy coefficient before simulation. Sensitivity analysis on Chezy coefficient ranging from 40 to 80 has been carried out. Results show little change on the hydrodynamic circulation patterns, water level variations, and velocity profiles. Thus, we conclude sensitivity of bottom roughness coefficient is not significant in this case. A default value of Chezy coefficient 65 m1/2/s is chosen. However, if the model is further applied to ecology or sediment dynamics studies, attention should be paid on the bottom roughness since the correlated bed shear stress have a critical impact on the sediment resuspension process.
5.3. Velocity Vorticity: Key Indicator of Hydrodynamic Circulation
Vorticity can be used to describe the spatially varying rotational character of a flow field [
38]. In this section, we focus on the vorticity of the horizontal velocity field in shallow water as the indicator of hydrodynamic circulation. The appearance of hydrodynamic circulation in the form of vorticity in the Lagrangian horizontal velocity field requires driving forces, which can be the wind shear stress gradients due to inhomogeneity of the wind velocity field, existence of submerged and emergent plants, the Coriolis effects, and the bathymetry variations. The influence of bathymetry variations is noteworthy compared to the other factors mentioned above [
39]. In Lake Ontario, Csandy (1973) schematized the study area into a long and narrow lake with parallel depth contours, and proposed the idea of “topographical gyres”, where the depth-averaged current direction is identical to wind direction in the shallow area, and opposite in the deeper area. This phenomenon was also observed in other shallow lakes around the world [
19,
40,
41]. However, these studies mainly focus on long and narrow lakes. For shallow lakes with other shapes and more rugged topography, where lateral differences are more significant, studies are lacking. With different bathymetry and bathymetry gradient, behavior of hydrodynamic circulation in the lake center area and the littoral zone varies (
Figure 4).
In previous studies, an analytical solution for the influence of bathymetry on the vorticity has been derived [
19,
42]. Hereby we present a brief review. Assuming a steady state with constant wind, i.e., the time derivative of horizontal velocity is 0, and the depth-averaged shallow water equations are
where
is the horizontal velocity vector with x component U and y component V,
is the Coriolis factor,
is the gravitational acceleration,
is the surface elevation,
is the total depth,
is the surface wind shear stress,
is the bottom shear stress,
is the density of water, and
is the horizontal eddy viscosity.
is the Del operator and
is the Laplace operator.
Introducing the vorticity of the horizontal velocity,
, and taking the y-derivative of Equation (5) and x-derivative of Equation (6), together with Equation (4), the governing equation for the vorticity becomes
The first term on the right-hand side represents the vortex stretching, the second term represents the vorticity input from wind, the third term is the vorticity sink due to bottom shear stress, and the last term is eddy-viscosity related. Compared to wind shear stress term, the latter two terms are one or more magnitudes smaller, thus, they are negligible.
With constant water density, the wind vorticity input term becomes
The first term in the bracket is related to wind stress curl, with constant wind, this term is 0, and a so-called barotropic topographic gyre is generated; the second term represents the velocity gradient influence on the vorticity.
Results from the numerical model also confirm the analysis above. With constant wind from northwest, the depth-averaged vorticity is shown in
Figure 13. In the littoral area, where the total water depth gradient is perpendicular to the wind speed, wind generates a large vorticity value. Along the southwest shore in the Taihu Lake, the current direction is the same as the wind, thus, it will produce a counter clockwise horizontal circulation. Even in the same gyre, the vorticity is larger, with the depth becoming smaller. Similar phenomena could be observed in the northeast boundaries of the Lake, the vorticity is negative, which means a clockwise gyre and the shallower part has the larger absolute value of vorticity. For locations where the depth gradient is small, for example, Dongtaihu Bay at the southeast part of Taihu Lake, the vorticity is almost zero.
These conclusions can be further illustrated by a numerical test with a flat bottom. Changing the bottom depth to a uniform value of 2 m and maintaining all other parameters yields model results, as shown in
Figure 14.
With constant depth, except for the grids near the boundary, depth gradients at all locations are zero. Thus, according to Equation (5), the velocity vorticity at each layer should be 0 as well, as illustrated in
Figure 14. Along most of the boundaries, for example, the southwest shoreline, the velocities of surface layer and bottom have the same magnitude but opposite directions; thus, the vorticity is positive and negative, respectively. Since the depth-averaged vorticity is the superposition of each layer’s vorticity, the values near the boundaries are approaching 0 (see
Figure 14c).
It may be concluded, both from the theoretical analysis and numerical tests results above, that the vorticity in the velocity field is related to directions of wind, current, and bathymetry gradient, as well as the value of water depth and wind shear stress. When the current velocity is in the same direction as wind, the major bathymetry influence is the depth gradient and depth itself. Wind has the largest influence on the vorticity when its direction is perpendicular to that of the bathymetry gradient. Wind related vorticity change is larger with larger wind speed and smaller water depth.
5.4. Transport Due to the Horizontal Circulation
Circulation and redistribution of lake water throughout the lake is further illustrated by the model results of the tracer experiment for the first half-year of 2008, with steady wind scenarios. At the beginning of the simulation, the tracer concentration all over the lake is 0. After that, tracers are released persistently in the first three months of the simulation, at a rate of 5 kg/s. Tracers are released at a chosen subbasin for each scenario. Here in this study, depth-averaged tracer concentration is presented.
Tracer movements reflect well the characteristics of the hydrodynamic circulation under steady wind conditions, discussed in
Section 4.1, e.g., when taking the tracer releasing experiment in Eastepigeal as a representative example (
Figure 15 and
Figure 16).
With constant southeast wind, initially, tracers move northwards via the shallow east margin of Epigeal Zone, where the depth-averaged flow velocity in the wind direction is the largest. Subsequently, the tracers penetrate into Gonghu Bay, with the counterclockwise circulation near the entrance of Gonghu Bay (
Figure 4). After that, tracers move westerly with the same circulation pattern, and enter the semi-closed subbasins of Meiliang Bay and Zhushan Bay, with smaller subbasin scale circulations. Finally, the tracers move southwards along the deeper center, and spread all over the lake. The entire time for the advection and diffusion process to fully mix the tracers in the entire lake is around 2 months’ time. Similarly, with constant northwest wind, tracers’ movements follow the hydrodynamic circulation pattern.
However, as shown in
Figure 16 on 30 April, tracer concentration is still high at the southern boundary of the lake, which means that for different locations in the lake, the duration of high tracer concentration due to redistribution of lake water by wind effects is different. For lake water quality management, it is crucial to assess the temporal influence at a critical spot of a random pollution emission. Thus, altogether, six important spots are chosen around the lake, including one drinking water intake point [
28] and five tourist places, as shown in
Figure 17, to analyze the influence of pollution release in the numerical experiment.
Taking spot A as an example, time series of tracer concentration with different release spots are shown in
Figure 18. Significant variance in the tracer concentration trend is also the direct result of hydrodynamic circulation, which determines moving trajectory of tracers between the interested location and the releasing spots. Under the constant southeast wind condition, water from Gonghu Bay, Eastepigeal Zone, and Center Zone has the largest influence. Location A lays near the joint point of Meiliang Bay and Gonghu Bay, where a subbasin scale circulation occurs (
Figure 4). With southeast wind, this circulation gyre will first drive water from inside the Gonghu Bay to point A. Then, several days later, water from Eastepigeal Zone and Center Zone will join this circulation gyre, and move towards location A, which could also explain why the peak of concertation curve of the latter two releasing is several days after the release stop time on 31st of March.
With constant northwest wind scenarios, which are represented as blue lines in
Figure 18, the situation is similar. Tracers from Meiliang Bay reach location A first, however, the distance travelled in this situation is longer than tracer released in Gonghu Bay with southeast wind. Due to mixing and diffusion processes, the peak concentration at location A with northwest wind is smaller than that with southeast wind. Water from Zhushan Bay and Center Zone arrives later, and introduces delayed concertation peaks. Specially to be noted, tracers released from Dongtaihu Bay show less influence due to longer travel distance, and weak hydrodynamics inside Dongtaihu Bay.
For lake authorities, the analysis of hydrodynamic circulation with tracer experiments under constant wind conditions may provide a quick tool to assess emergency pollution impact on critical area before carrying out a more complicated hydrodynamic and water quality model. Also, the model could help to design temporary observation locations based on the environmental conditions. For example, given the emission location at Eastepigeal Zone and southeast wind, an observation points located to the north of Xishan Island is more useful than to the south.
5.5. Is Large-Scale Water Transfer Effective?
Diluting and flushing of Taihu Lake with water pumped from the Yangtze River has been the largest engineering intervention to improve water quality in Taihu Lake after 2007. Previous studies with field data collection have been carried out to identify the variation in water quality and biological indices before and after the water transfer. Results show that the effectiveness of water transfer project remains of concern. Positive effects occurred only in parts of Taihu Lake, while in some other areas, water quality even worsened [
29,
30,
43].
To investigate the hydrodynamic circulation condition change due to the project, numerical experiments are carried out using the original discharge data of 2008, and the water transfer rate from literature (Table 1 in [
30]). The water transfer in 2008 lasted from 23 January to 9 June, with a total transferred water volume of 870.2 million m
3, which is around 6.7 million m
3/d.
Model results suggest that, except for Dongtaihu Bay and Gonghu Bay, no significant change in hydrodynamic circulation condition occurred with and without water transfer. Total discharges of Gonghu Bay and Dongtaihu Bay are shown in
Figure 19.
For other open cross-sections, such as southwest and northwest, volume exchanges slightly fluctuate, but for the entrance cross-sections of semi-closed subbasins like Meiliang Bay, the total discharge remains the same.
Model results of total volume exchanges across the cross-sections show that transferred water neither significantly changed the water amount within each subbasin nor enhanced the volume exchange between subbasins. However, it indeed stimulated the hydrodynamic circulation in Gonghu Bay and Dongtaihu Bay.
Although the water transfer project contributes little to changing the total volume in most subbasins, hydrodynamic circulations still transfer and mix water persistently throughout the lake. Considering that the water transferred from the Yangtze River contains less nutrient or pollution than the average of Taihu Lake, with the amount of transferred water that stayed in the subbasins, the nutrient concentration would consequently change; a tracer experiment is conducted to further elucidate the redistribution of transferred water. This scenario lasts for the entire year of 2008, and tracers are continuously released on Wangyu River boundary at a concentration of 5 kg/m3 between 23 January and 9 June 2008.
Model results are shown in
Figure 20. As can be observed from the top right subplot, after 129 days on 9 June, a large amount of transferred water has already reached Gonghu Bay, Meiliang Bay, and Center Zone. Later, the tracers move southeastwards before they leave the lake from Taipu River. To be especially noticed, almost no tracers pass through the entrance of Zhushan Bay, making tracer concentration virtually zero inside this subbasin. Volume exchange between Zhushan Bay and other parts of Taihu Lake is more severely blocked, partially due to a weak hydrodynamic circulation near the entrance of the subbasin.
6. Conclusions
Wind-induced hydrodynamic circulations and associated transport and mixing processes in large shallow lakes play a significant role in their environmental and ecological processes. Knowledge of these physical mechanisms helps understanding, with the hope to solve eutrophication problems like algae blooms, which happen increasingly frequently in large shallow lakes, like Taihu Lake. Previous studies emphasized mostly the environmental, biological, and ecological aspects of eutrophication problems. However, the driving forces of nutrient transport, algae scums, pollutant mixing, and the underlying physical processes, such as wind-driven or anthropological horizontal circulation, are usually overlooked, but they might be key factors leading to algae blooms, and hitherto ambiguous. In this study, hydrodynamic circulation in shallow lakes is defined as the large-scale movement of water in the lake basin. A three-dimensional, numerical Delft3D model of Taihu Lake, driven by steady and/or unsteady wind, river discharge, rainfall, and evaporation, is used to quantitatively illustrate the complex hydrodynamic circulation and their effects in transporting and mixing within the lake.
A relative stable circulation pattern is found to be formed after 2 days, on average, with steady wind, where the overall hydrodynamic circulation structure, i.e., direction, intensity, and position, is determined by wind direction, wind speed, and initial water level. Vertical variations of horizontal velocity are found to be related to the relative shallowness of water depth. In the shallow marginal area, flow at the bottom and surface layers has the same direction, and the surface flow has a larger velocity, while in the deeper area, the bottom flow reverses, opposite to the wind direction and with a larger velocity. Volume exchange between subbasins, influenced by wind speed and initial water level, differs due to the complex topography and irregular shape. With unsteady wind, these findings are still valid to a high degree. Vertical variations in hydrodynamic circulation are found to be very important in explaining the surface accumulation of algae scums in Meiliang Bay in summer. Vorticity of current velocity, as the key indicator of hydrodynamic circulation, is determined by wind direction, bathymetry gradient, and water depths, while the maximum change of velocity vorticity happens when wind direction and bathymetry gradient are perpendicular to each other. Furthermore, we use Lagrangian-based tracer tests to estimate emergency pollution/leakage effects and to evaluate water transfer effects. One emergency pollution leakage point is added to the model to demonstrate the effect on five tourism hotspots and one drinking water intake point, suggesting that the model application may serve as an operational management tool. The water transfer project shows that even a large-scale water transfer (about 1/5 volume of total lake volume in 138 days from Yangtze) does not alter the hydrodynamic circulation and volume exchanges between subbasins significantly, but it succeeds to transport and mix the imported Yangtze River water to the majority of Taihu Lake area.
This study may well be extended further to provide insight into the spatial and temporal biological process corresponding to wind induced hydrodynamic circulation in Taihu Lake, and similar large shallow lakes, to support ecologically sound design and implementation for lake restoration.