A Detailed Analysis on Hydrodynamic Response of a Highly Stratified Lake to Spatio-Temporally Varying Wind Field

Abstract

Wind is generally considered an important factor driving the transport and mixing processes in stratified enclosed systems such as lakes and reservoirs. Lake Abashiri is one of the instances of such a system. For these systems, typically, the temporally unsteady but spatially uniform nature of wind has been assumed for simplicity. However, the spatial non-uniformity of wind could significantly alter compound hydrodynamic responses. In this study, such responses were investigated under the continuous imposition of different inhomogeneous wind conditions by applying numerical models and integrated analysis. The resultant tracer transport in both uniform and non-uniform wind cases was insignificant for the total study period of 9 days. However, under the short interval of $T_i$, where $T_i$ is the internal fundamental period, different behaviors of both surface particle transport and the internal wave field were identified. Particularly, surface mass transport responses to higher spatial wind variance were obviously different from those in the uniform case. In addition, internal wave spectra under strong wind magnitude, which has low spatial variances, became identical to that of uniform wind; however, there were some discrepancies in the non-uniform case in the wave spectra under the influence of weak-to-moderate wind of high spatial variances. The results could provide an in-depth understanding of the lake’s hydrodynamic response to inhomogeneous wind which could improve water management in lakes and reservoirs.

1. Introduction

Water quality in lakes and reservoirs has implications for its beneficial applications, such as fisheries’ management, water supply quality, phytoplankton population, and control of the manganese level [1,2,3]. Such conditions are controlled by both nutrient loadings and physical factors such as mixing, and transport processes driven by winds [4]. As a major source of energy in lakes, coastal waters, and some estuaries, wind forcing plays an important role in the formulation of lake mixing, circulation, and internal seiche excitation [5,6,7]. Wind blowing over a lake rapidly distributes its momentum to water and consequently induces surface current and transport of floating objects such as phytoplankton along the wind direction. In addition, turbulent mixing induced by wind can affect the structure of the stratification, leading to changes in hydrodynamics and, in turn, species composition, proliferation of toxic algae, and dwindling of drinking water supplies [8,9].

Information on wind-induced transport, mixing, and internal waves is typically explained by numerical approaches, whereby the source of wind input is frequently retrieved from a point-based stationary wind field to spatially represent the phenomenon [5,7,10,11]. However, such an assumption is often oversimplified, although its dependency on the surrounding orography is crucial, resulting in a lack of spatial representation of the inhomogeneity of actual meteorological conditions [12]. Moreover, spatial wind fields often vary across large lakes due to pressure gradients over the water surface [13], and across small lakes due to sheltering effects [14,15]. Consequently, wind-induced lake hydrodynamics are insufficiently explained with the inaccurate surface boundary condition. Therefore, the coupling effects of the horizontally non-uniform wind and the real bathymetry of the concerned study sites are needed to realistically simulate the phenomena [16,17]. This could be achieved through a stand-alone uncoupled hydrodynamic model interacting with spatially interpolated wind [18,19], climate products [20], and numerical downscaling methods [12,21].

Generating the spatial variability of wind either by observation or by numerical simulation illustrated improvements in hydrodynamic simulations that were underpinned in a number of studies. For example, the simulation in Lake Kinneret in Israel improved the wave amplitudes and reduced the overall phase error after the spatially inhomogeneous wind stress was appropriately estimated [22]. The flow structure in Lake Belau changed drastically under the induction of spatial variation in the wind [23]. In Lake Ledro, TKE and its dissipation rate were reduced because of the representation of spatially inhomogeneous wind, as they are sensitive to the wind intensity [14]. The assumption of wind uniformity may often fail to reproduce the natural circulation condition and its characteristics in Spermberg Reservoir [24]. The spatial variations in the wind field separated the natural internal seiche excitation into different modes, as well as acting as a controller for higher horizontal modes in Lake Iseo [15]. Furthermore, the spatial variations in wind-generated wind curl, formulating seasonal regimes of Ekman layers resulting from the accumulation of particulate organic matter (POM) at the Tokyo Bay head, are driven by the pumping velocity [25]. Normally, comparisons of the response of lake hydrodynamics to spatial uniform and non-uniform wind fields are done with visualization, without being quantified. Quantification of these phenomena based on wind inhomogeneity levels is therefore necessary to elucidate their in-depth features.

To understand the complicated response of hydrodynamics to spatial and temporal variations in wind, meteorological wind-induced lake hydrodynamics were numerically simulated in Lake Abashiri from 21 November to 29 November 2011. To the best of our knowledge, this study is the first to quantify the interactions between spatial variations in wind forcing and lake hydrodynamics in Lake Abashiri. Therefore, under different spatial and temporal wind conditions, there are three main objectives in this study: (1) to quantify and understand the effects of wind forcing on surface mass transport and water current; (2) to quantify the effects of wind forcing on tracer movement; and (3) to estimate the potential mixing and to quantify the effects of wind forcing on internal seiches. The overall aim is to elucidate and quantify the differences in lake hydrodynamics between spatially uniform and spatially non-uniform wind.

2. Materials and Methods

2.1. Site Background

Lake Abashiri in Hokkaido (Figure 1a) is a small lake with a surface area of 32.87 km², a shoreline of 42 km, a maximum depth of 16.1 m, a mean depth of 7.2 m, and an altitude of 0.4 m. Lake Abashiri connects Abashiri River at the southernmost part of the catchment area and estuary at the northeast end (Figure 1b). Seawater frequently flows backward from the sea to the Abashiri River, so that the lake is stratified mainly with salinity. Lake Abashiri is an example of a highly stratified system linked by upstream freshwater intrusion and downstream seawater intrusion. This complex system forms a strong halocline, distinguishing the lake’s layers. The lower layer of saline water is deoxidized up to <0.1%Sat due to a long residence time [26]. The lake provides important recreational activities, such as fishing and sports, throughout the season. However, eutrophication has recently become a concern due to the nutrient loading [27]. Lake Abashiri can emit 11% of the total methane in Lake Hokkaido and produce dissolved methane (DM) concentrations five times higher than the estimated one-dimensional methane model (DM concentrations were estimated by surface area) [27]. Such an issue could be problematic if the lower-layer gas was released from the lake. As an important factor of such a system, continuous wind exceeding 10 m/s could upwell the lake, releasing CH4 to the environment [7].

2.2. Observation Data

Hourly values of wind speed and direction taken from Tokoro station (Figure 1b) from 21 November 2007, to 30 November 2007 were obtained from the Automated Meteorological Data Acquisition System (AMeDAS) which was then used to validate the WRF model. AMeDAS is a collection of Automatic Weather Stations (AWS) and is monitored by the Japan Meteorological Agency (JMA) [28].

2.3. Weather Research and Forecasting Model (WRF)

This study employed WRF for numerically downscaling and simulating wind inhomogeneities [29] with four nested domains (Figure 1a), with 200 × 200, 190 × 190, 178 × 178, and 166 × 166 cells and horizontal resolutions of 10 km, 3.3 km, 1.1 km, and 333 m, respectively. The temporal resolution was set for every 20 min.

Table 1. Weather research and forecasting model configurations.

Domains	Domain 1	Domain 2	Domain 3	Domain 4
Resolution	10 km × 10 km	3.3 km × 3.3 km	1.1 km × 1.1 km	333 m × 333 m
No. of grid points	200 × 200	190 × 190	178 × 178	166 × 166
Microphysics	WSM3	WSM3	WSM3	WSM3
Longwave radiation	RRTM	RRTM	RRTM	RRTM
Shortwave radiation	Dudhia	Dudhia	Dudhia	Dudhia
Cumulus parameterization	Kain–Fritsch	Kain–Fritsch	Kain–Fritsch	None
Land-surface option	Noah-MP	Noah-MP	Noah-MP	Noah-MP
Surface layer option	Monin–Obukhov	Monin–Obukhov	Monin–Obukhov	Monin–Obukhov
Planetary boundary layer	Mellor–Yamada–Janjic	Mellor–Yamada–Janjic	Mellor–Yamada–Janjic	Mellor–Yamada–Janjic

illustrates WRF configurations and the domain setup for the simulation. In the resolution from the outermost to the innermost domain, the former covers almost the entire country of Japan to capture overall atmospheric phenomena, while the latter includes Lake Abashiri surrounded by a relatively flat area. In this context, hourly ERA5 [30] reanalysis data were dynamically downscaled from 0.25° spatial resolution. ERA5 is the most recent global reanalysis product that has benefited from a decade of advancing model forecasting and assimilation methods [31]. In addition, studies have reported that ERA5 is good at improving simulation in complex terrains and outperforms several other global reanalysis products in terms of surface wind [32,33,34]. We used the rapid radiative transfer model (RRTM) scheme for longwave radiation and the Dudhia scheme for shortwave radiation. Furthermore, we applied the Monin–Obukhov scheme [35], the Noah-MP land-surface model, and the Mellor–Yamada–Janjic TKE scheme for the surface layer, land surface, and planetary boundary layer, respectively. In addition, substitution of higher land-use resolution taken from Global Map Japan [36] and observation-based dynamical nudging to the initial conditions were implemented to further improve model performance. Each run began at 00:00 UTC and generated outputs every 20 min.

2.4. Fantom-Refined Hydrodynamic Model

The Fantom-Refined [37,38] model was developed following equations of continuity and 3-D Navier–Stokes with incompressible and Boussinesq approximation together with thermal flux (details of the model can be found in Appendix A). The model is further improved with a local mesh refinement (LMR) technique based on a structured grid which uses rectangular grids similar to a tree-based grid. This method uses two different perspectives of a horizontal grid, macro, and micro views. The macro view determines the refining unit (hereafter the “container”), and each container is discretized in a structured manner with arbitrary horizontal resolution (hereafter the “micro view”). This method introduced the concept of a container that defines local groups with constant grid resolution in each direction. With this container, we may specify any computational cell in the domain with a simple, modifiable LMR and straightforward indexing. The study, however, applied only the macro view for simplification and for computational efficiency.

The horizontal resolution for the hydrodynamic model was set at 200 m in each direction (Figure 2a). In terms of vertical resolution, finer grids (0.1 m) were set at the middle of the lake’s depth to represent the halocline fluctuations visibly and accurately, while coarser grids were used at the remaining depths (Figure 2b). Vertical temperatures were distributed evenly along the depth, and salinity was initialized with measured data, in which the halocline existed at a depth of approximately 5 m below the water surface.

Complexity exists when a clear understanding of lake hydrodynamics is challenged by both spatially varying drivers (wind and bathymetry), whereby lake hydrodynamics potentially react differently under the spatial distribution of the former, which is controlled by the surrounding orography. Meanwhile, variations in depth and geometry of the latter are of great concern in further alternating the mechanism. A deep comprehension of how a specific main driver affects hydrodynamics is the key to gaining a step-by-step understanding of the complex physical behavior of lake hydrodynamics. Therefore, flat bathymetry of 16 m was assumed to separate the complex effects of the real bathymetry.

To investigate the effect of spatial and temporal variations in wind on hydrodynamics in the lake, a comparative study was conducted between unsteady but spatially uniform wind, and unsteady and spatially non-uniform wind (

Table 2. Fantom simulation durations setup. For every hour of time increment, the simulations of surface particle transports and internal seiches were reset every 8 h ($T_i)$ for the entire 9 days of the simulation. Tracer transport and another case of internal seiche were simulated with no reset for the entire 9 days of the simulation.

	Duration of Simulation		Time Steps
	Spatial Uniformity	Spatial Variations
Surface particle transports	$T_i$ (8 h)	$T_i$ (8 h)	10 s
Tracer transport	9 days (tracer was placed at the 21st)	9 days (tracer was placed at the 27th)	10 s
Internal seiche	9 days and $T_i$ (8 h)	9 days and $T_i$(8 h)	10 s

). By a simple interpolation technique, the horizontal and temporal resolutions of the wind field extracted from the fourth domain of the WRF model were refined to 200 m and 10 s, respectively, to match those of the hydrodynamic model. The uniform wind field was created by taking wind data at the center grid of the lake. Wind shear stress applied over the water surface was computed explicitly with the function of (U10) [38].

$${\tau }_s={\rho }_a{C}_D\left|U_w\right|U_w$$

(1)

where ${\rho }_a$ is the air density, $C_D$ is the drag coefficient, and $U_w$ is the wind magnitude at 10 m above the water surface, i.e., U10. The model applied a bulk formula [39] to estimate sensible and latent heat fluxes, and short and long wave radiations from meteorological weather conditions.

2.5. Indicators for Spatial Variations in Wind Heterogeneity

Quantification of wind inhomogeneities is challenging due to the complex distinction between wind direction and wind speed. Recently, several studies have attempted to quantify the spatial variations in wind. For example, smaller root mean squared errors of vectors u and v indicated uniformly spatial variations and vice versa [14]. This technique, however, using the vectors u and v, cannot quantify the invariance of wind speed or wind direction, which could be directly addressed using 2D mean filter analysis [40].

The spatial consistency of the wind was determined based on the integrated coefficient of variance and a 2D mean filter [40].

$$f{\left(x,y\right)}_n=\frac{1}{{\left(2n+1\right)}^2}\sum _{\left(s,t\right)\in S_{xy}}g{\left(s,z\right)}_{n-1}$$

(2)

where $f{\left(x,y\right)}_n$ is the filtered image, $n$ is an integer and can be selective with incremental size (e.g., 1, 2, 3, …, m), m is arbitrarily chosen for grid size, $g{\left(s,z\right)}_{n-1}$ is the value in each pixel of the grid of the previous filtered image; and s and z are the corresponding coordinates of the grid itself. As noise reduction comes at the expense of blurring the image (a 2D visual representation of any pixel), if the resulting image ($f{\left(x,y\right)}_n$) has features in common with all the previous filtered image ($f{\left(x,y\right)}_0, \dots , f{\left(x,y\right)}_{n-1}$), the wind field will be considered to be spatially homogeneous. Details of the computation procedures are described in Appendix B.

2.6. Indicators for Temporal Variations in Wind Heterogeneity

The metric Pearson correlation was applied to measure the linear relationship between two sets of data. It is the ratio of two variables’ covariance to the sum of their standard deviations.

$$r_{xy}=\frac{{\sum }_{i=1}^o\left(x_i-\overline{x}\right)\left(y_i-\overline{y}\right)}{\sqrt{{\sum }_{i=1}^o{\left(x_i-\overline{x}\right)}^2}\sqrt{{\sum }_{i=1}^o{\left(y_i-\overline{y}\right)}^2}}$$

(3)

where overbar $x$ and $y$ are the average value of the sample population, and $o$ is the sample size. To analyze the homogeneity of temporal wind field, a correlation map constructed between the center grid (the 2D time series data at the center of the chosen dataset) and the surrounding ones was calculated during the simulation period.

2.7. Wedderburn Number

The Wedderburn number, a non-dimensional and systematic relationship between the wind-induced shear stress and the baroclinic pressure gradient, was computed to identify the mixing condition in the reservoirs and lakes [41].

$$W=\frac{g^{\prime }{h}_e^2}{u_{*}^{2}L}$$

(4)

in which $h_e$ is the epilimnion depth, $L$ is the wind fetch length, $u_{*}$ is the friction of wind velocity, and $g^{\prime }$ is the reduced gravity.

2.8. Surface Particle Transport, Tracer Simulation, and Internal Seiche

In principle, wind triggers two main important physical phenomena, lake circulation and surface and internal waves, both of which control the mass transport in lakes. To investigate the movement of the former, virtual particles are often placed and numerically simulated together with the atmospheric forcing interaction. This approach can be visually quantified by the Lagrangian tracking method, whereby the simulated particles often representing phytoplankton, fish eggs or larvae, are released, and the dynamics of the particle trajectories are recorded in time and space over the lake [12,19,21]. In this study, surface particles were placed uniformly on the water body to elucidate the effect of spatial and temporal variations of the wind field. In addition to the Lagrangian tracking method, tracer tracking is frequently adopted for its representation of either two- or three-dimensional transport of fluid containing pollutant [10,20]. To this aim, a normalized tracer concentration of 1.0 was set at the northern part of the lake.

Investigation on the development of the two-layer model of the temperature profile [42] was the key to further intensive empirical studies in an effort to classify the mixed-layer dynamics in lakes [43] and to readily determine the first vertical internal wave mode V1H1 [44]. To get a better insight into the internal waves that are generated by both homogeneous and inhomogeneous wind, a temporal variation of density interface for each simulation period was extracted. After that, spectral analysis was applied to evaluate and highlight the periodic evolution of internal seiches in both uniform and non-uniform cases. This technique can easily reveal hidden oscillatory motions in a time series structure or can detect the dominant frequencies of a periodic time series in terms of both spatial and temporal resolution [45]. Welch’s method was used to estimate power spectral density (PSD), where the Hanning window function was adopted to weight the time series dataset.

Standing internal waves in confined basins can be categorized by nodal points on the vertical (V) and horizontal (H) components, V_aH_b, where $a$ and $b$ are the fundamental vertical and horizontal mode components, respectively [46]. Finally, we computed analytical internal waves based on the calculation described in [6] to theoretically identify the period of each mode. We also compared the performance of the numerical model forced by different simulation periods with the analytical model in a variety of wind conditions. Based on the procedure described in [6], the periods of internal seiches were estimated by a 1D multi-layer hydrostatic linear model as

$$T_{a,b}=\frac{2L}{b\sqrt{g{\beta }_a}}$$

(5)

where $a$ and $b$ are the fundamental vertical and horizontal modes, respectively, $g$ is the gravitational acceleration, and ${\beta }_a$ is the eigenvalues assuming a two-layer model system ($a$ = 1). The fundamental period can also be calculated following the Merian formula:

$$T_i=\frac{2L}{\sqrt{g^{\prime }\frac{h_1{h}_2}{h_1+h_2}}}$$

(6)

where $g^{\prime }=g\frac{{\rho }_2-{\rho }_1}{{\rho }_2}$ is the reduced gravity of the two-layer system, ${\rho }_1$ and ${\rho }_2$ are the densities of the two layers, respectively, $h_1$ is the thickness of the upper and $h_2$ is the thickness of the lower layer. Additionally, the fundamental period $T_i$ can also be estimated as $T_{a=1,b=1}$.

3. Results

3.1. WRF Model Validation

Wind data taken from AMeDAS (Figure 3a) illustrated a strong north-easterly wind with a clear diurnal-cyclical alteration during the simulation period. Strong winds evolved during the early morning and weakened from the beginning of the afternoon until nocturnality, with a high frequency of maximum and minimum intensities of approximately 10 m/s and around 2 m/s, respectively. These periods with strong magnitudes and clear periodicity could potentially influence the synoptic-scale occurrences related to the lake’s hydrodynamics. From the 24th to the 26th, low to moderate magnitudes were noted, with a tendency toward mild fluctuations. Validation between the model and observed data showed the general overestimation of model data, with RMSE and MBE values of 6.46 m/s and 0.17 m/s, respectively. There was also an overestimation of wind magnitudes on the first half of the 21st and the night before the 22nd. WRF showed minor phase differences from the 23rd, from the early morning to the late afternoon. A relatively consistent trend between the two datasets was seen from the 24th onward. In addition to wind magnitude, wind direction showed an overall consistent distribution between model and observation data (Figure 3b,c). Fortunately, the disagreement between numerical and observed data did not impose significant effects on the study, since the aim of this study was to distinguish uniform and non-uniform wind-induced lake hydrodynamics instead of perfectly achieved modelled data.

3.2. Spatio-Temporal Wind Quantification

The spatial inhomogeneity levels of wind speed and direction were calculated with nine levels of the 2D mean filter (m = 9). The overall results are presented in (Figure 4a). The study period indicated general homogeneous wind direction in accordance with strong wind speed. Moreover, strong, and weak winds corresponded to low and high averaged coefficient of variance (Appendix B), ${\overline{C}}_v$ scores, revealing an inverse relationship between the wind magnitude and the quantification score. According to a rough wind inhomogeneity estimation [40], values of ${\overline{C}}_{v\_dir}$ less than 0.01 and ${\overline{C}}_{v\_spd}$ less than 0.66 are considered spatially homogeneous, while values of ${\overline{C}}_{v\_dir}$ higher than 0.01 and ${\overline{C}}_{v\_spd}$ higher than 0.76 demonstrate spatial inhomogeneity. In addition to the inverse relationship mentioned above, moderate wind (3–5 m/s) illustrated an ambiguous pattern in spatial distribution. As a result, the complexity of spatial distribution could potentially signify the difference in the lake’s hydrodynamics. Based on the wind magnitude and the spatial inhomogeneity indicator, Figure 4b–d represent the spatial variations or spatial uniformity in wind. On the last day of the simulation (29 November 2007 23:50), very weak spatial wind (~1 m/s) produced insufficient momentum, disturbing the direction of the wind over the surface water body. Additionally, light air (~2 m/s) combined with moderate breeze (~4 m/s) (24 November 2007 15:00) formed a wind curl as a result of the shear stress differences in space. Figure 4c,d illustrated light (~4 m/s) to fresh breeze (~10 m/s), respectively. As we can see, near-spatial uniformity was noted in the wind over the lake’s surface. Therefore, we can assume a spatially uniform wind field during the strong wind periods, while the moderate wind-induced periods can alter the lake’s hydrodynamics due to the spatial variation of winds.

Pearson correlation coefficients were estimated in the 2D time series wind data between the central grid and its surrounding grids over the simulation period (Figure 5). Correlation coefficients were calculated using Equation (3). For the wind direction, temporal inhomogeneity exhibited an overall positively moderate correlation from the middle point to its nearby grids, while weaker relationships were located farther from the center grid (Figure 5a). Therefore, insignificant spatial wind similarity covered the entirety of Lake Abashiri. This could be explained by the spatially disturbed wind direction during weak-to-moderate wind periods. On the other hand, surprisingly, wind speed illustrated relatively spatially homogeneous correlation coefficients all over the area, with moderate correlation coefficients. In conclusion, wind speed represented acceptable similarity (~0.5) but wind direction illustrated higher spatial disturbance over the simulation period.

3.3. Surface Particle Transport and Tracer Transport

Figure 6 shows the numerical simulations of wind-induced surface particle transport. Theoretically, particles traveled unidirectionally following the wind direction under spatially uniform wind. This is consistent under strong and moderate wind speeds, which inversely correlate to the low ${\overline{C}}_v$. Figure 6c,f display a constant movement of particles following the water currents driven by the strong and moderate winds, respectively. Therefore, low ${\overline{C}}_v$ generated translational motions in surface particles along the wind directions. In contrast to the strong wind, weak breezes generated the distinguishable motions of the water currents during a short time scale of the fundamental period, $T_i$ (8 h), because of the large ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$, showing the tendency of dispersing particles (Figure 6k). Compared to the spatial uniform wind, spatially distorted weak wind (~1 m/s) initiates a short distance of movement after the $T_i$ period. However, particles started to travel omnidirectionally after providing sufficient drag force, causing the particles to scatter in different directions (Figure 6k). This means that after a long period (>$T_i$) of spatial variations in winds (large ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$), the movement of surface particles would be far more different from that of spatially uniform wind.

Tracer mass transport was initially placed at the small estuary where the substances from the seaside were first provided in this lake. Figure 7 shows the concentration and distribution of the tracer under different wind conditions. Although the effects were relatively similar in both cases due to the prevailing wind direction, there were some discrepancies noted in the first few days of the simulation. For the first two days of strong prevailing NE wind, the tracer dispersion was mainly transferred to the estuary mount. As a result, the mass transport was blocked by the right-hand side boundary. The tracer was stirred up by the vertical circulation, which resulted in the tracer spreading over half of the water surface. Under uniform wind, the stable and homogeneous wind increased the concentration of tracer further south, while a lower concentration level was observed in the case of non-uniform wind (Figure 7b). The results of the inhomogeneous wind, which generated spatially disturbing water currents compared to that of spatially uniform wind, illustrated changes with less concentration in the southernmost part of the lake. By the end of the simulation, tracer mass transport was seen to be completely mixed up in the two wind cases, leaving no traits of a dominantly concentrated area.

3.4. Internal Seiche

In general, the uniform and non-uniform wind shared similarities in inducing internal waves under the strongly stratified lake system (Figure 8). The depth of halocline computed under both wind cases did not present any complete upwelling under strong wind imposition. This was further confirmed by accessing the relationship between wind and the response of stratification (Figure 8c). There were roughly 96% of the total periods that are classified as “regime three” [43]. This indicated that the most prominent feature is internal seiching. Under this regime, strong wind durations are less than the entrainment time, illustrating the rare event of complete vertical mixing [43]. Regime four accounted for the remaining periods, showing the overall domination of buoyancy force. However, the complete dominance of regime three during the simulation period emphasized the importance of internal waves in this lake system.

A Wedderburn (W) (Figure 8a) number of 4 to 10 was often observed during the first period of strong prevailing NE wind force which gradually transferred energy and induced internal waves for approximately 5 days (Figure 8b). Strong oscillation was shown after the 26th and continued to highly fluctuate for roughly two days. This pattern was similar from the 28th onward. During this period, there was some time when W < 1 indicating potential mixing. This could explain why the resuspension events eventually happened under high wind speed [47]. However, the short periods of strong wind induced could not entirely mix up the strongly stratified lake system. After that, the imposition of weak-to-moderate wind did not generate enough energy, indicating the domination of buoyancy and damping effects.

The depth of the halocline in both wind cases under strong wind illustrated noteworthy features. On the first few days of the simulation period, both spatially uniform and non-uniform wind exhibited their similar features because of the strong wind speed and low ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$ (Figure 4a). The minor discrepancies started to exhibit on the next few days under the effects of weak-to-moderate wind speeds. Between the 26th and 28th, strong wind imposed on lake displayed similarities in the high fluctuations following the notable Wedderburn numbers. According to Figure 8b,c, after the 28th, the buoyancy and damping effects were much higher than those produced by low wind speeds. It is also noted that more mixings were seen under the imposition of spatially uniform wind (Figure 8b). Therefore, the effects of low-to-moderate wind speed accompanied by high fluctuation in wind direction must not be neglected.

4. Discussion

4.1. The Differences in Surface Mass Transport and Water Current between Spatially Uniform Wind and Spatially Non-Uniform Wind

The weak and spatially inhomogeneous wind forcing generated the uneven water current and could explain the significant differences in the particle distributions (Figure 6) among three cases during the short time scale of $T_i$, as shown in Figure 9.

Low ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$ in the case of moderate to strong wind induced the water current to flow unidirectionally over the water surface. In addition to the high wind magnitudes, the wind direction did not change rapidly as it did in the weak wind condition (Figure 9a,b), which induced enough energy to push almost all the surface particles unidirectionally. This is consistent with the numerical simulation in [12]; the study illustrated that the horizontal transport of particles was stronger in winter Föhn due to the high and persistent magnitude of wind speed. As mentioned earlier, low values of ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$ wind fields can indicate uniform wind fields. This can further reaffirm why wind speed is important and highly correlated with the wind direction as well as the water current.

On the other hand, in spatially non-uniform wind (Figure 9d), wind speeds were mostly below 2 m/s, and the low wind speed was often accompanied by omnidirectional wind direction [48], resulting in a disturbance of the movement of water current omnidirectionally across the lake. In contrast, spatially uniform wind (Figure 9c) produced a unidirectional water current in most parts of the lake, revealing the major differences in water current directions in both wind cases. Therefore, scattering effects would occur in the former wind case.

The assumptions of the flat bathymetry and the enclosed system, additionally, affirm that water-current curl (Figure 9d) can be formed under the imposition of spatial variations in wind due to the shear stress gradients where sheltering effects can be neglected in this lake [20]. A two-gyre circulation developed in Figure 9d, in which the clockwise gyre in the southern part coexisted with the counterclockwise gyre; meanwhile there was no such circulation under the strong wind condition. Depending on the lake size, this phenomenon is often observed when applying wind inhomogeneities as a governing factor [24]. Additionally, because Lake Abashiri is located in a flat and simple orography, the sheltering effects are minor, thus confirming the major effects of wind inhomogeneity. Based on the spatial distribution of water current (Figure 9), although the geometry of Lake Abashiri is simple, littoral zones in the south and at the middle generated clockwise flow, thus increasing the horizontal circulation (Figure 9d); meanwhile, this cannot be seen in spatially uniform wind conditions under the period of $T_i$ (Figure 9c). As mentioned earlier, based on the calculated ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$ scores, the high values when moderate winds were observed would create wind curls due to the gradient. The uneven distribution of wind shear stress can be sufficiently important for understanding the movements of water currents as well as particulate organic matter (POM).

4.2. Correlation between Spatio-Temporal Indicators and Tracer Mass Transport

A strong and prevailing north-easterly wind (NE) was the main factor controlling the similar movement of the tracer down to the south of the lake on the first few days; the tracer was finally mixed up entirely at the end of simulation period (Figure 7). As the tracer was completely mixed up after nine days of simulation, the complex alteration in hydrodynamics might be difficult to observe. To overcome that issue, a tracer placed on the 27th of November revealed the noteworthiness of spatially varied wind imposition (Figure 7c). After one day, the tracer started reducing its concentration at the north boundary, while a higher concentration accumulated under non-uniform conditions and remained at the north of the lake. At the end of the simulation, concentrations remained higher in the case of spatial variations in wind, while the tracer concentration was fully mixed up with uniform wind.

Interestingly, the result obtained from the spatial variations in wind (Figure 7c) indicated that the movement of POM would transport it to the center of the lake, which is similar to the transportation of POM gathering at the center of the Tokyo Bay [25]. High moisture POM tends to accumulate at the center of the bay head [49,50]. Although the effects of wind curls and the assumptions of flat bathymetry in this study can be the factor producing the high concentrations at the lake center, real bathymetry can also increase the effect because it increases the water current velocity at the littoral zones. Further study is needed to understand the formation of positive and negative wind curl, and its effects on upwelling and downwelling under spatial variations in wind in combination with the real bathymetry.

4.3. Effects of Temporal and Spatial Frequencies of Winds to Internal Seiche

Overall, Figure 8 illustrates that the differences in internal waves of both wind cases were minor after nine days of the simulation. Noticeably, under strong wind during the first days of simulation (Figure 8), the height of the halocline responded comparably under both wind scenarios, wherein the differences as well as the major component causing the difference cannot be distinguished. To differentiate, the weak to strong wind-induced internal wave should be tested over a brief period. Simulations of internal waves within $T_i$ were conducted to understand the differences under both wind conditions.

Figure 10 shows that the uneven distribution of wind direction was highlighted, with inverse correlation existing at southern part of the lake and a gradient from the center to the lake boundary in magnitude. Consequently, upwelling points were scattered, reducing the similar features of wind-induced internal wave. Considering the differences in spatial and temporal wind features, the results of spectral analysis generated a huge difference in PSD as well as frequency. While uniform wind coincides with V1H3 mode, non-uniform wind appeared to match to V1H2 mode. These alterations were obviously the main factors changing the behaviors of internal seiche in this system.

An eight-hour simulation of moderate wind (25 November 2007 00:00) indicated the overall stability of the lake’s system. The temporal correlation of wind speed and direction showed a reasonable similarity of both spatial uniformity and non-uniformity (Figure 10). The moderate wind featured North-East (NE) as the general prevailing wind, resulting in most of the upwelling locations being located at the south end of the lake. This feature indicated a similar internal seiche in V1H2 mode, whereas both wind-induced internal wave cases disagreed in V1H3 mode, with a lower PSD in non-uniform wind compared to that of uniform wind. In addition, lower wind speed correlated to the PSD generated by the internal waves.

Finally, for the strong wind (~10 m/s), winds were predominantly easterly (E) during the period of ($T_i$). This could be seen from the upwelling locations accumulated at the western end of the lake boundary. Spectral analysis revealed that the prevailing wind (E) of spatially uniform and non-uniform wind excited internal waves with relatively similar power spectral densities at modes V1H1, V1H2, and V1H3. However, the numerical model did not match the mode V1H4. Moreover, both wind cases excited a similar PSD of internal waves, which reaffirms the similarity in both spatial uniformity and nonuniformity under strong wind fields.

5. Conclusions

The study aimed to investigate the responses of the hydrodynamics of a lake under both spatially uniform and spatially non-uniform wind fields. Lake Abashiri was selected as the field, and the flat bottom of the lake was assumed for the analysis in order to focus on the effect of wind nonuniformity on the lake response. The realistic nonuniform weather was generated by the weather research and forecasting model with four level nesting grids over the area. To quantify the level of inhomogeneity of wind speed and wind direction, the wind products from WRF were then analyzed based on a 2D mean filter. Three of the lake’s physical responses, including surface particle transport, tracer tracking, and internal waves, were examined under the different wind conditions. The analysis highlighted the key factors that effect changes when applying different wind scenarios.

Firstly, strong wind corresponding to a low level of variances both for wind speed and direction, ${\overline{C}}_{v\_spd}$ and ${\overline{C}}_{v\_dir}$, provided sufficient momentum to unidirectionally push the surface particles horizontally to the shoreline within the fundamental period of this stratified lake, $T_i$. On the other hand, wind curl was formed under the imposition of a weak-to-moderate wind field of high values of variance both in speed and direction. The surface particles moved only short distances but not unidirectionally, so they would be transported differently for the succeeding phase in the case of the uniform wind.

Secondly, when appropriately placing the tracer mass transport under wind of high variances, the possible effects of wind non-uniformity could be observed in the transport process in the POM concentration at the bay head of Lake Abashiri. Particularly, in the present study, a much higher concentration was seen in the more realistic wind condition, while a smaller concentration was seen in the case of spatially uniform wind.

Finally, the internal wave spectral analysis in both wind cases showed comparably similar patterns when considering the whole simulated period. Furthermore, the potential mixing, as defined by the Wedderburn number, showed that strong wind provided enough energy to fluctuate and potentially upwell the downwind along the fetch length. However, the system could not be mixed; this was probably due to the short duration of the strong wind and the assumption of uniform bathymetry. The analysis under the shorter time scale ($T_i$), revealed that the weaker the wind, the larger difference in power spectral density. This further elucidated the importance of weak-to-moderate wind effects on lake hydrodynamics.

This study emphasized the significant differences between spatially uniform and spatially non-uniform winds when comparing their interactions with surface mass transport, tracer tracking, and internal waves, excluding the effects of the bathymetry. Although it is important to consider the effects of bathymetry, they can, together with spatially non-uniform wind, alter and contribute to the complexity of lake hydrodynamics. As a result, future studies will include bathymetry effects to further analyze the sensitivity of lake hydrodynamics, coupled with spatial and temporal variations in wind.