Earth rotation effects on the internal wave field in a stratified small lake

Earth rotation effects on the internal wave field in a stratified small lake: Numerical simulations The Princeton Ocean Model was applied to the Sau Reservoir, a medium-sized Mediterranean reservoir (5.7 km2) located in Catalonia, Spain, during the summer season, when the water column is continuously stratified and forced by a breeze regime with velocities of up to 3-4 m/s. Based on our simulations, the internal wave field has been analysed and the numerical results compared with the field data previously analysed by other authors. The model adequately reproduces all the significant modes observed. The simulations show the importance of rotational modes on the internal wave field. The Burger number S for all rotating internal waves is on the order of 10−1, and the internal Rossby radius R is on the order of 102 m, that is, smaller than the width of the lacustrine area of the reservoir (103 m). Specifically, two rotating third vertical modes were found during the analysed period: the first, a 24-hour period, was forced by the wind and the phase rotated clockwise; and the second, a 12-hour period, can be interpreted as a second azimuthal horizontal mode of a Poincaré wave. A second vertical mode of 8 hours was observed to rotate counter clockwise, although in this case the Earth’s rotation appears not to have been of importance because the same results could be obtained without taking into account the Coriolis force. Finally, two first vertical modes, one of 6 hours and the other of 5 hours, were observed, although with no rotational behaviour. These modes correspond, respectively, to the first and second horizontal stationary modes. Further analysis of the simulated velocities shows the existence of a net clockwise flow along the shore as a consequence of the dominant mode of 24 hours, with a mean velocity of 0.5 m/s that reverses at a distance of about 300 m from the coast. Preliminary results based on the simulated velocity field predict horizontal trajectories of passive particles of up to 1 km per day and vertical displacements of up to 5 m across the entire water column.


INTRODUCTION
Recent advances in the numerical modelling of hydrophysical processes in aquatic systems show increasing agreement between simulations and observations (Lemmin et al., 2005;Beletsky & Schwab, 2001) and allow users to be more confident when using numerical modelling to study the dynamics of various water basins, especially when field data are limited or unavailable. Analyses of modelling results combined with measurements (if available) help us to better understand the dynamics of complex aquatic systems and clarify some basic hydrodynamic problems that have yet to be solved, including the influence of wind stress curl and baroclinicity on large-scale circulation in lakes  and the contribution of internal waves to the energy cascade in hypolimnetic waters (Imberger, 1998). In enclosed basins, stationary waves are very important as ubiquitous and even dominating features of general circulation (Mortimer, 1974). In continuously stratified systems, such as Mediterranean reservoirs during the end of the summer season, high-order vertical modes can be observed to dominate the internal wave field (Pérez-Losada et al., 2003), even when these oscillations are not forced by the wind (Roget et al., 1997). The motion of the water body prior to a specific wind force has also been shown to be important in determining the internal wave response, which depends on morphometry, stratification and forcing characteristics (Antenucci et al., 2000). This is the context in which numerical models are useful tools for predicting the complex internal wave dynamics of lakes and for studying their effect on mass transport and particle resuspension (Stocker & Imberger, 2003;Wüest & Lorke, 2003).
Models (Mortimer, 1979;Salvadé et al., 1988) considering one-and two-dimensional continuously stratified layers (Munnich, 1996) have been widely used to reproduce stationary internal modes in lakes, but more recently it has been demonstrated that three-dimensional numerical models can reproduce the internal wave field in large lakes as well (Casulli & Cheng, 1992;Hodges et al., 2000;Rueda & Schladow, 2003). In terms of morphometry, Wetzel (2001) considers a lake to be large when its surface area is larger than 50 km 2 . In this study, we run a 3D hydrodynamic model to simulate internal waves in a medium-sized reservoir (maximum surface area of 5.7 km 2 ), which can be considered a small lake according to Wetzel's criteria.
In small lakes, the velocities of dominant flow patterns are often rather low (∼10 −2 m/s), so the Rossby radius is also relatively small and on the same order as the characteristic horizontal scale of the system, allowing the possible influence of the Earth's rotation to develop. Examples of flows in small and large narrow lakes, which are affected by rotation, have been given by Hamblin et al. (1978) and Roget & Colomer (1996). In this study, we used the Princeton Ocean Model (POM), which accounts for rotation effects and has shown good results not only in oceanographic applications but also in estuaries and large lakes. For instance, the POM has been used to study the topographic effect of wind-induced circulation and its influence on vertical mixing in Lake Mask in Ireland (Bowyer, 2001). There are also wellknown POM simulations of the Great Lakes. In Lake Erie (Kelley et al., 1998) and Lake Michigan , the POM has been used to study general circulation, the annual stratification cycle and the inter-annual variability of that cycle. In interdisciplinary research, water quality models have been coupled with the POM to examine the effect of interior heating on dissolved oxygen (Blumberg & Di Toro, 1990) and sediment resuspension in lakes (Lou et al., 2000;Chen et al., 2004). In Japan, the POM has been used together with a biochemical model to study water quality in Lake Hamana (Taguchi & Nakata, 1998).
Here, we present the results of POM simulations carried out at the end of the summer season in the Sau Reservoir in Catalonia, Spain. This reservoir, with a capacity of 170 hm 3 , is a system of reference for Spanish limnologists and has been studied continuously since 1965 (Armengol et al., 1986;Armengol et al., 1999;Marcé et al., 2007). This study is based on data presented by Vidal et al. (2005) [hereafter VCCS]. Following Pérez-Losada et al. (2003), who studied the Boadella Reservoir, VCCS have shown that the Sau Reservoir can react as a four-layer system dominated by an internal oscillation of 24 hours. Using a two-dimensional model (one vertical and one horizontal direction), which did not consider the Coriolis effect, VCCS reproduced a third vertical mode of an internal seiche of 25.3 hours and interpreted the 24-hour oscillation as a stationary baroclinic mode. Two additional oscillations with periods close to 5 and 6 hours were measured in the lake, although only the 5-hour oscillation was interpreted as a first vertical mode corresponding to the first horizontal mode. With a 2D model, VCCS found that the stationary first vertical and horizontal mode was 5.48 hours. The 2D model was not able to simulate the oscillation period of 12 hours found in the data recorded in the reservoir, and the authors suggested that the 12-hour oscillation probably corresponds to a sub-daily wind speed pattern. Another internal oscillation of approximately 8 hours was identified in the recorded data by VCCS but was not interpreted. In this paper, we will discuss the characteristics and the driving forces of all the modes described by VCCS.
The effect of the Earth's rotation on the internal wave fields of small lakes is usually thought to be negligible (Antenucci & Imberger, 2003). The main objective of this paper is to show that in small lakes, there are internal waves with phases that propagate slowly enough (low frequency internal waves) to feel the effects of the Earth's rotation and so describe a circle of inertia on the order of the characteristic scale of the lake. In this case, crests begin to propagate around the boundary of the lake (Poincaré waves) rather than across the lake. More generally, the comparison of numerical results with field recorded data presented in this paper leads us to consider using the POM, a freely distributed model, to simulate the internal wave field in small lakes under synoptic meteorological conditions. This is reviewed in the discussion session, where we also discuss the possible influence of internal waves on processes of ecological importance, such as particle transport and resuspension (Eadie et al., 2002).

SITE AND METHODS
The Sau Reservoir, at 41 • 58 4.28 N and 2 • 24 46.32 E, has a total length of 18 km, although the lacustrine zone is only 3600 m long with a maximum width of 1300 m. Its maximum height is 426 meters above sea level (m asl) and the deepest bottom is located at 365 m asl. During the analysed period the surface height was 410 m asl and the maximum depth was 45 m.
The bathymetry of the Sau Reservoir used for the simulation is presented in figure 1a. The figure has been rotated (see the arrow pointing north) and the farthest part of the river tail has been eliminated. This is justified because the surface level of the reservoir during the period under study was not as high as it could have been and because the removed area was not a stagnant section that could be associated with the reservoir but was a very narrow river flowing into it.
When the characteristic basin dimension L is close to or larger than the internal Rossby radius, defined as R = c/f , i.e., the ratio of the non-rotating phase speed of the internal wave c to the inertial frequency f , the Earth's rotation can affect the dynamics of internal gravity waves. Additionally, the Burger number, which is defined as S = c/L f , indicates the relevance of the Earth's rotation to internal wave dynamics. When S is less than or equal to unity, the Earth's rotation can be considered important. The phase speed of the internal wave can be estimated as c = λ/T, where T is the period of the internal wave and λ ∼ 2L and λ ∼ L for first and second horizontal modes, respectively, where L is the characteristic length of the lacustrine area of the reservoir (here 2000 m).
The characteristic transversal horizontal scale of the Sau Reservoir is relatively small (∼10 3 m), but the phase speed of the internal waves can also be low (∼ 10 −2 m/s for low frequency modes). Therefore, the Burger number is on the order of 10 −1 , and the internal Rossby radius is on the order of 10 2 m, which is one order of magnitude less than that of the width of the lacustrine area of the reservoir. Accordingly, we can expect low frequency internal waves to be affected by the Figure 1. a) Bathymetric map where contour intervals, in meters, refer to depth below the free surface. Predominant wind directions during day and night are shown. Locations of the measurement stations (Station 1 and Station 2) and some points (Pn; n = 1.6) that are considered throughout the text to discuss the simulated results are indicated. Letter T stands for a transect perpendicular to the shore line from where data are shown in Figure 11. The location where the meteorological station was placed is indicated by the letters MS, and the location of the dam is also indicated; b) Wind velocity at 10 m above the ground recorded at station MS from September 7 to 12, 2003; c) Wind direction recorded at station MS from September 7 to 12, 2003. Time zero corresponds to the beginning of September 9. a) Mapa batimétrico dónde los intervalos de contorno, en metros, son referidos a la profundidad por debajo de la superficie libre. Se muestran las direcciones del viento predominantes durante el día y la noche. Se indican las localizaciones de las estaciones de medición (estación 1 y estación 2), y de algunos puntos (Pn; n = 1.6) los cuáles son considerados a lo largo del texto para discutir los resultados de la simulación. La letra T presenta un transecto perpendicular a la línea de la costa dónde los datos son mostrados en la figura 11. La localización de la estación meteorológica está indicada por las letras MS y la localización de la presa también está indicada; b) Velocidad del viento a 10 m por encima del suelo registrada en la estación MS del 7 al 12 de septiembre del 2003; c) Dirección del viento registrada en la estación MS del 7 al 12 de septiembre del 2003. El tiempo cero corresponde al inicio del día 9 de septiembre.
Earth's rotation, as will be discussed more completely in the following section.
The temperature and velocity data used in this study were continuously measured at two stations in the reservoir from September 9 to 12, 2003. More precisely, at Station 1 (Fig. 1a), an electromagnetic current meter was deployed close to the surface, as was an acoustic Doppler current meter, positioned face down to receive data from 20 depth bins, each 2 m in height. Temperatures at this station were recorded with an 11-sensor thermistor chain. At Station 2 (Fig. 1a), another thermistor chain with four regularly spaced sensors was deployed. Several specific locations (P), marked in figure 1a by the numbers n = 1. . .7, were to be used for the analyses of the simulated results. The location of the dam is indicated on the right hand side of figure 1a by a white arrow.
Due to the prevalent summer anticyclonic conditions, the wind blew toward the dam during the day (Fig 1a), but there was a slight breeze in the opposite direction during the night with a maximum wind speed of about 3 m/s, as shown in figures 1b and 1c. Time equals zero indicates the beginning of September 9 th when the field campaign started. The predominant wind directions corresponding to the lake are shown on the map of the system presented in figure 1a. Because the wind direction varied during the night, we plotted two arrows in figure 1a to indicate the range of this variation. The meteorological station was placed close to the shore at the location indicated by the letters MS in figure 1a (upper right hand corner).
The mean temperature and Brunt-Väisälä frequency (N) profiles during the first day of the studied period at Station 1 are shown in figure 2. A homogeneous surface layer 5 m deep can be observed, followed by a continuously varying stratified layer ranging from 5 to 15 m in depth, where N reaches a maximum value of 0.045 s −1 . From this point down to the lake's floor at 37 m in depth, the stratification continuously diminishes.
The POM is a model based on hydrostatic three-dimensional primitive equations with the Mellor and Yamada closure scheme (level 2.5) and terrain-following (sigma) coordinates. Because the model is widely documented (http://www. aos.princeton.edu/WWWPUBLIC/htdocs.pom/), in-depth physical and numerical details are not reviewed here, and only the parameters most relevant to our simulations are given. The reservoir was presented as a grid of 70 × 70 m rectilinear horizontal cells. The time step, which satisfied the constraints based on the Courant-Friedrich-Levy criterion, was reduced to 0.05 s for the external mode and to 3 s for the internal mode to gain accuracy close to the boundaries. Twelve vertical cells were used for the simulations discussed, although other simulations with a larger number of sigma levels were run showing the same main features. Zero normal velocities were used as the lateral boundary conditions. At the surface, wind stress was calculated as τ = ρ a C D V 2 10 with the drag coefficient C D defined by Hasselmann (1988). No mass input or output was considered in the simulations because during the period of the campaign the river inflow and water withdrawal was sufficiently negligible. The initial temperature profile for all simulations was set to the mean profile during the first day of measurements (Fig. 2). The model was run over 10 days, forced by measured wind data (September 6 to 15), but the results are discussed only for the four day-period (September 9 to 12) when experimental data were available.
The study of the horizontal structure of internal waves was based on filtered deviations from  the equilibrium of the surface level so that the progress of internal motions could be followed. The method resulting from the hypothesis that internal oscillations are accompanied by oscillations at the surface level, typically 1/1000 smaller in amplitude, was initially proposed by Mortimer (1963) and was later used by different authors, e.g., Caloi et al., 1986, Sirkes, 1987and Lemin et al., 2005.

RESULTS
The most important result of our study is that the POM successfully reproduces all of the significant modes in the observed internal velocity field spectra. This can be seen in figure 3a where the spectra of the measured and simulated velocities in a 124 • direction (following the main axis of the reservoir) at Station 1 (see Fig. 1b) are shown. In figure 1a, both experimental and numerical velocity data are taken at a depth of 14 m, where all the modes contribute to the horizontal velocity. Both spectra show modes of about 24, 12, 8, 6, and 5 hours, with the 24-hour mode being the most energetic. The model also reproduces a four-layer structure as was also deduced from experimental field data obtained by VCCS (see introduction). In figure 3b, an instantaneous profile of the simulated horizontal velocity along the 124 • direction at Station 1 is presented, and a 4-layer structure can be clearly observed. However, the water column structure varies in space and time as will be discussed. The persistence of a background four-layer structure is shown in figure 3c. Here, the cross-correlation functions of the 24-hour component of the horizontal velocity across the main axis of the lake at four different depths at Station 1 are presented. More precisely, for a time lag of 0 hours (± an integer times 24 hours), the velocities at 1 and 15 m and at 15 and 26 m are shown to be completely out of phase. Furthermore, the velocities at 15 and 34 m are in phase, so that the flow at depths of 26 and 34 m is out of phase, confirming that horizontal flow velocity reverses three times at a depth with a periodic behaviour of 24 hours. Accordingly, based on POM numerical simulations, the mode of 24 hours can be identified as a third vertical mode, as was also shown by VCCS.
A comparison of numerical and experimental data can also be made from figure 4, where the contour plots of the a) computed and b) measured horizontal velocity components at Station 1 are presented. Note that experimental data were not available above 3 m or below 30 m (indicated in Fig. 4a by dashed lines). As observed, above a depth of 16 m, both contour plots present similar patterns, and there is a main periodicity of about 24 hours. At greater depths, however, the simulated data show that shorter periods dominate, although this is not very clearly observed from the experimental data. A more detailed discussion of the five modes identified in figure 3a that lead to the complex structures observed follows.

The 24-hour mode
As discussed in the Site and Methods section, the 24-hour mode can be analysed by looking at the 24-hour component of the surface elevation, plotted for different times within one period in figure 5. The brighter end of the scale on Figure 5. a-f) Band pass-filtered (24-hour) surface elevation for September 10 at different times when the model was forced with the measured wind. a-f) Elevación superficial filtrada con paso de 24 horas el 10 de septiembre a diferentes momentos dónde el modelo ha sido forzado con el viento medido. the plots indicates the maximum elevation above the horizontal level, while the darker end denotes the maximum elevation below it. A comparison of figures 5a, b, c, d, e, and f shows that the higher vertical displacements are found to propagate along the shoreline in a clockwise direction with a periodicity of 24 hours. Furthermore, within the lacustrine zone, only one horizontal nodal line is observed, indicating that the horizontal structure of this wave corresponds to a first horizontal mode, which is a third vertical mode as was discussed with reference to figure 3c.
The drag effect of the wind is responsible for the 24-hour gyre shown in figure 5, which follows the daily breeze regime. To discuss the forced character of this oscillation, a simulation was run with a theoretical periodic wind of 30 hours and a maximum amplitude of 3 m/s, as is approximately the case with real wind data. The results show that the system also reacts with a rotating oscillation of 30 hours (not shown), which confirms the forced character of the lowest frequency oscillation. In both cases -wind forcing of 24 hours and 30 hours -the corresponding forced modes do not show any rotating behaviour when the Coriolis parameter in the model is switched off ( f = 0). During the study period, the phase velocity for this mode is 0.046 m/s, and the corresponding internal Rossby radius is 470 m, that is, smaller than the amplitude of the lacustrine area of the reservoir (10 3 m), which confirms the importance of the Coriolis effect; the Burger number is 0.47.
In figures 6a, b, c, d, e, and f, the daily flow around the lake can also be identified from the surface elevation along the shore at Stations P1, P2, P3, P5, P6, and P7 (Fig. 1a). Note that the maximum elevation moves progressively between consecutive stations along the shore. Although the velocity between stations is not exactly constant, the mean velocity corresponds to a gyre every 24 hours. Alternatively, in figure 6g, a progressive vector diagram of the horizontal velocity at the centre of the lake also shows clockwise rotation over 3 consecutive days. Note that, apart from a daily gyre, there is a net displacement along the main axis of the lake in the direction of higher wind speed (dashed arrows). From another viewpoint, the existence of a net flow can be deduced from figure 6h, where the daily vertical integrated horizontal velocity along the shore is presented for a section perpendicular to Figure 6. a-f) Surface elevation at consecutive points (P1, P2, P3, P5, P6, and P7) along the shore (see Fig. 1A); g) Progressive vector diagram of the velocity at a depth of 1 m at the centre of the lake from September 9 to 11, 2003; h) Time and depth averaged values of the horizontal alongshore velocity component during one day of the simulation. a-f) Elevación superficial en puntos consecutivos (P1, P2, P3, P5, P6, y P7) a lo largo de la costa (ver Fig. 1A); g) diagrama del vector progresivo de la velocidad a 1 m de profundidad en el centro del lago del 9 al 11 de septiembre de 2003; h) Componente de la velocidad paralela a la costa promediada en el tiempo y en la profundidad durante un día de la simulación. the water edge. This particular section is located near the east shore of the lake and is indicated in figure 1a with the letter T. As shown in figure 6h, a mean velocity of about 0.5 cm/s rotating clockwise is found close to the shore. This velocity decreases towards the lake's interior and reverses at a distance of about 180 m.

The 12-hour mode
The 12 h mode could not be reproduced as a stationary mode with the 2D model used by VCCS and was suggested to be related with specifics of the wind pattern.
The simulated results show that this mode rotates clockwise. Rotation patterns at the surface related to this mode can be visualised from the simulated results based on the measured wind data, but they can be visualised more effectively with the forcing of an idealised sinusoidal wind with a period of 24 hours. In figure 7, we present a sequence of the 12-hour component of the results obtained with the idealised forcing; the same sequence recurs every 12 hours. Figure 7. Band pass-filtered (12-hour) a-d) surface elevation at different times when the model was forced with an idealised sinusoidal wind with a period of 24 hours and maximum velocity of 3 m s −1 and e) simulated horizontal velocity component along a 124 • direction at Station 1 over September 9 and 10. a-d) Elevación superficial filtrada con paso de 12 horas en diferentes momentos donde el modelo ha sido forzado con un viento sinusoidal idealizado con un periodo de 24 horas y con una velocidad máxima de 3 m s −1 y e) Componente de la velocidad horizontal filtrada con paso (12 h) simulada a lo largo de la dirección 124 • en la estación 1 del 9 al 10 de septiembre.
A clockwise rotation pattern can be observed, which presents two nodal lines within the lacustrine zone. This oscillation corresponds to a second azimuthal horizontal mode and not to a first, as it did for the 24-hour rotating oscillation. Both modes are third vertical modes, as shown in figure 7e, where the 12-hour band pass-filtered horizontal velocity component along the main axis of the lake is presented during a period of 48 hours. The flow reverses four times within the water column. Note that in the simulations, the idealised wind forcing did not have any 12-hour pattern that could directly force this mode. Because the period of this oscillation is now shorter than the inertial period (T i = 17.9 h), which was not the case for the oscillation of 24 hours, the 12-hour rotating wave can be interpreted as a second azimuthal mode of a Poincaré wave. As in the 24-hour case, when the model is set at f = 0, no rotating behaviour can be observed, demonstrating the importance of the Earth's rotation for the development of this mode. It should also be remembered that this mode could not be reproduced as a stationary wave using a 2D model. During the study period, the internal Rossby radius for this mode is 470 m, which is smaller than the amplitude of the lacustrine area of the reservoir, and the Burger number is 0.47, which confirms the importance of the Coriolis effect. The phase velocity is of 0.046 m/s.

The 8-, 6-and 5-hour modes
The POM was also able to reproduce the 8hour mode measured in the lake, although VCCS could not reproduce it as a stationary mode using a 2D model. Figures 8 a, b, c, d, and e present the plots of the 8-hour component of the surface elevation for this mode at different times within one period. As was the case for the 24-hour oscillation, higher amplitudes can be seen moving along the shore, now in a counterclockwise direction, Figure 8. Band pass-filtered (8-hour) of a-e) surface elevation at different times when the model was forced with the measured wind on September 11 and f) simulated horizontal velocity component along a 124 • direction at Station 1 over September 11 and 12. a-e) Elevación filtrada con paso (8 h) en diferentes momentos donde el modelo ha sido forzado por el viento medido el 11 de septiembre y f) Componente de la velocidad horizontal filtrada con paso (8 h) simulada a lo largo de la dirección 124 • en la estación 1 del 11 al 12 de septiembre. Figure 9. Band pass-filtered (6-hour) of a-c) surface elevation at different times when the model was forced with the measured wind on September 9 and d) simulated horizontal velocity component along a 124 • direction at Station 1 over September 9. a-c) Elevación filtrada con paso (6 h) en diferentes momentos donde el modelo ha sido forzado por el viento medido el 9 de septiembre y d) Componente de la velocidad horizontal filtrada con paso (6 h) simulada a lo largo de la dirección 124 • en la estación 1 durante el 9 de septiembre.
with a period of 8 hours and consistent with the structure of a first horizontal mode as observed in figures 8a-e, where a single nodal line can be observed. In figure 8f, the 8-hour component of the simulated horizontal velocity at Station 1 is plotted for different depths along the main axis of the lake. The horizontal lines in the plot show a three-layer structure. It can be concluded that the velocity field presents an 8-hour oscillation, rotating counterclockwise, which corresponds to a second vertical mode and a first horizontal one.
Similar results were also simulated by the POM when the Coriolis parameter was switched off, indicating that the Coriolis effect is irrelevant in terms of the rotational character of this mode. Because the period of this flow around the lake is about 8 hours, the daily breeze may energise it every third round and it can be observed throughout the entire measured period. Although the origin of this internal wave is not clear, it is interesting to see that the POM reproduces a rotating oscillation of the same period as the one observed, which could not be identified with a 2D model and which, according to our results, would develop even without the Earth's rotation.The phase velocity of this mode is 0.14 m/s, the Burger number 1.44, and the internal Rossby radius 1440 m.
The periods of 6 and 5 hours observed in both the experimental and simulated spectra in figure 3 are found to be first vertical modes corresponding to the first and second horizontal stationary modes of the lacustrine zone along its main axis. The horizontal structure of surface elevation for both modes is given in figures 9 and 10, together with the corresponding band pass-filtered simulated horizontal velocity component along the 124 • direction at Station 1 for one day. From both figures it can be seen that in this mode the water column behaves like a two-layer system, as it does for the first vertical modes. The estimated phase velocity for the mode of 6 hours is 0.19 m/s, and the internal Burger and Rossby radius are S = 1.95 and R = 1.950 m. For the mode of 5 hours, the estimated phase velocity Figure 10. Band pass-filtered (5-hour) of a-c) surface elevation at different times when the model was forced with the measured wind on September 9 and d) simulated horizontal velocity component along a 124 • direction at Station 1 over September 9. a-c) Elevación filtrada con paso (5 h) en diferentes momentos donde el modelo ha sido forzado por el viento medido el 9 de septiembre y d) Componente de la velocidad horizontal filtrada con paso (5 h) simulada a lo largo de la dirección 124 • en la estación 1 durante el 9 de septiembre. is 0.11 m/s and S = 1.13 and R = 1130 m, which are above the critical values. This explains why, as observed in figures 9 and 10, these modes do not rotate. VCCS identified the mode of 5 hours as the first horizontal mode but did not discuss the presence of the 6-hour mode (Vidal et al., 2005;Fig. 9).

DISCUSSION AND OVERVIEW
Our results show that in small lakes the Earth's rotation can affect the internal wave field, which can present rotating and non-rotating stationary and forced modes. Several high vertical modes can be excited simultaneously in the continuously stratified water column of Mediterranean reservoirs, and the use of a 3D model is needed to study their dynamics. In this work, we have shown that the POM, a widely tested public oceanographic circulation model, can be a useful tool to simulate internal wave fields in small lakes for synoptic meteorological variations. Furthermore, the computational time required is quite reasonable. In our case, it takes about 15 minutes per day using a fairly standard computer by today's standards, with a 4GHz Quad Core and a LINUX UBUNTU 64-bit server.
Although the POM is a useful tool with which to characterise the basic structure of the velocity field for synoptic forcing conditions according to our results, we have not proved that it could be used in medium-sized reservoirs for intermediate time predictions because its mixing closure (Mellor-Yamada, 1982) has not been validated in cases of closed water systems with a small fetch. However, the latest version of the POM (we used pom97) contains a wet and dry scheme that allows the modelling of water surface area variations due to changes in water level, which commonly happens in reservoirs.
The influence of internal oscillations on the ecology of lakes and reservoirs has been illustrated from different perspectives (O'Sullivan & Reynolds, 2004). For example, the significance of the internal wave field in the Sau Reservoir for the ecology of the lake can be examined in terms of particle transport, which can be discussed on the basis of the simulations of velocity fields presented in previous sections. More precisely, it can be considered that the local velocity V determines the displacement of the particle, k, which, depending on the time t, is located at different positions, X i (i = 1, 2, 3 and refers to the Cartesian component), so that X i (k, t + 1) = X i (k, t) + V i (x, y, z, t)Δt, where Δt is the time step of the simulated background velocity field. In our case, as stated in the numerical set-up section, the time step was 3 s. Furthermore, the vertical Cartesian component of the velocity, V 3 (x, y, z, t) ≡ w (x, y, z, t), was obtained from the component of the velocity perpendicular to the sigma coordinate simulated by the POM, ω, using the equation w = ω + σ[u(∂H/∂x) + v(∂H/∂x)], where σ is the sigma coordinate, H the mean depth at that particular point, and u and v the two horizontal velocity components of the flow.
Due to the chaotic dynamics of transport phenomena in lakes (Boffetta et al., 2003), an ac- Figure 11. a and b) Horizontal velocity field when the particles (indicated by stars in the central panels) were released at 8 and 30 m of depth, respectively; c, e, g and i) Horizontal trajectories of the particles within one day; d, f, h and j) Corresponding vertical trajectories of some of the particles. a y b) Campo de la velocidad horizontal donde las partículas (indicadas con estrellas en los paneles centrales) son dejadas a 8 y 30 m de profundidad, respectivamente; c, e, g y i) Trayectorias horizontales de las partículas durante un día; d, f, h y j) trayectorias verticales de algunas partículas. curate forecast of the transport must account for forcing variability and would require a robust statistical study, which is not the aim of this work. However, a few examples of the trajectories of passive particles in the Sau Reservoir are given in figure 11. Several particles were released simultaneously at depths of 3, 8, 13, and 30 meters at the locations indicated with stars in figures 11 c, e, g and i. The instantaneous velocity fields at depths of 8 and 30 meters are presented in Panels a and b. The approximate mean vertical density structure was presented in figure 2. After release, the particles were tracked for one day and their horizontal trajectories reached up to about 1 km (Figs. 11 c,e,g and i). Although statistical analysis is needed to reach a conclusion, at the top layers where the wind drag effect is larger, the horizontal particle displacements would be expected to be longer than in the bottom layers. As expected in a non-stationary velocity field, particle tracks depart from the streamlines at the time of release. The vertical trajectories of the particles released along the main axis of the lake are shown in Panels d, f, h, and j. It can be observed that vertical particle displacements are important in the top and bottom layers due to the existence of internal waves. Experimental data presented by Marcé et al., 2007, and recorded in similar wind conditions in Sau reservoir show internal wave amplitudes of 5 m.
Due to internal waves, horizontal velocities appear above the bottom that can generate resuspension. Depending on the characteristics of the bed (the water content and biological activity of the sediment, if freshly settled materials), the critical flow velocities for resuspension vary considerably. In fact, although the majority of studies of resuspension refer to flow velocities greater than 0.10 m/s, critical velocities for resuspension as low as 0.016 m/s are also reported in the literature (Van Raaphorst et al., 1998). Our simulations do not show velocities above the bottom greater than 1.5 cm/s, so under similar conditions resuspension should not be expected, although it could occur under stronger winds or different stratifications conditions. Finally, internal waves contain a large amount of the entire mechanical energy of the system, which can eventually be partially transformed into turbulence, thereby inducing mixing and enhancement of the mass fluxes in the interior of the lake (Wüest & Lorke, 2003).