Thermal functioning of a tropical reservoir assessed through three-dimensional modelling and high-frequency monitoring

Urban lakes and reservoirs provide important ecosystem services. However, their water quality is being affected by anthropogenic pressures. The thermal regime is a strong driver of the vertical transport of nutrients, phytoplankton and oxygen. Thermal stratification can modify biogeochemical processes. In this paper, a three-dimensional hydrodynamic model was implemented and validated with high-frequency measurement of water temperature. The simulation results were in agreement with the measurements. For all simulation period, the model performance was evaluated based on hourly values, presenting a maximum RMSE of 0.65 oC and Relative Error of 2.08%. The results show that high-frequency measurement associated with a three-dimensional model could help to understand and identify the reasons for the changes in the thermal condition of a shallow urban lake. The impact of the stream inflow on the temperature was highlighted, showing that during higher discharge events, when the river temperature is colder than the lake water, it flows into the lake deeper layers. The inflow water sank to the deeper layers where the lake morphology changes. The model showed an impact along the entire lake, showing the importance of monitoring the inflow water temperature. This modelling tool could be further used to study specific patterns of reservoir hydrodynamics.


INTRODUCTION
In urban areas, lakes and reservoirs provide many ecosystem services. Environmental changes at a local scale, within their watersheds, and at a global scale are affecting the water quality and ecological functioning of urban lakes (Birch & McCaskie, 1999;Friese et al., 2010;Li-Kun et al., 2017), resulting in an increase in phytoplankton biomass and the occurrence of potentially toxic cyanobacteria blooms (Davidson et al., 2016;Ho & Michalak, 2015;Jöhnk et al., 2008;Whyte et al., 2014).
Most urban lakes are medium-size and shallow (Birch & McCaskie, 1999;Gong et al., 2016). They may react strongly and rapidly to external hydrological and meteorological forcing, presenting many alternation periods of thermal stratification and mixing. Thus, high-frequency measurements are necessary to capture transient changes and to help understanding the relationship between factors that may impact ecological status. However, high-frequency monitoring is rarely used, especially in urban lakes in tropical regions.
Due to the complex interplay between ecological and hydrodynamic processes in lakes, with high heterogeneity in space and time, three-dimensional models help in better understanding the processes involved (Gong et al., 2016;Soulignac et al., 2017Soulignac et al., , 2018Zbiciński & Ziemińska-Stolarska, 2017). Usually in hydrodynamic three-dimensional models, hydrodynamic processes are represented in a sub-model that provides water temperature and velocity in time and space to an ecological sub-model. Three-dimensional models are increasingly implemented to better understand the lake ecosystem response to external pressures, like climate forcing, nutrient loading, and pollutant input (Aparicio Medrano et al., 2013;Chanudet et al., 2012;Soulignac et al., 2018;Vinçon-Leite & Casenave, 2019;Yi et al., 2016). However, the use of three-dimensional modelling in urban lakes is much less reported (Soulignac et al., 2017), especially in tropical water bodies which are subjected to distinct conditions of external forcing and may present different responses compared to temperate water bodies. Furthermore, the use of high-frequency measurements in three-dimensional modelling studies in urban lakes is rarely found. Chanudet et al. (2012) used a three-dimensional model with high-frequency water temperature measurements in the shallow Nam Theun 2 Reservoir (Laos). However, the authors obtained high errors between simulated and measurement values in the reservoir deepest layer, especially during strong thermal stratification events. Soulignac et al. (2017) used three-dimensional model with high-frequency water temperature measurements in the shallow Lake Créteil (France). The authors could well represent the water temperature dynamics over all depths. Jin et al. (2000) used three-dimensional model with high-frequency water temperature measurements in the shallow Lake Okeechobee (USA). The authors obtained a considerable root mean square error of water temperature of 9.18%. Jin et al. (2000) highlighted the wind-induced wave and sediment deposition/resuspension. However, none of them highlighted the impact of the inflow water temperature. Belico (2017) used a unidimensional model in the lake Pampulha (study site of the present article). The one-dimensional model could not rightly simulate the temperature in the deeper layers of the lake. Other modelling studies were previously applied to Lake Pampulha (Silva, 2014;Silva et al., 2016Silva et al., , 2019 however, none of them used three-dimensional modelling with high-frequency measurement in different depths.
This paper presents a rare application of three-dimensional modelling and water temperature high-frequency measurements in a shallow and medium-size tropical urban lake using hourly meteorological and inflow data. The main objective of this study is to associate high frequency measurements and three-dimensional hydrodynamic modelling in order to investigate the influence, during thermal stratification periods, of inflow water discharge and temperature and lake morphology on the lake thermal behaviour.

Study site
Lake Pampulha is a reservoir located in Belo Horizonte, a city with 2.5 million inhabitants (Instituto Brasileiro de Geografia e Estatística, 2014) in the southeast of Brazil. Lake Pampulha and its modern urban ensemble, was designed by the architect Oscar Niemeyer and built in the 1940s. The climate characteristics of the lake region are of tropical altitude type (Companhia de Pesquisa de Recursos Minerais, 2001) and summarized in Table 1. Dry season occurs from April to September. Annual mean inflow in Lake Pampulha is 1.30 m 3 /s and 2-years return period inflow is 106.7 m 3 /s (Campos, 2004,   Plec et al.

3/25
Located at an altitude of 801m, Lake Pampulha has a surface area of approximately 2.0 km 2 and a volume of about 10 million m 3 at 801 m (spillway crest) with a residence time of 89 days (considering the annual mean inflow). The lake is a shallow reservoir with a maximum depth of 16 m and an average depth of 5 m. The lake has a polymictic thermal regime, in which the thermal stratification of the water column during a few days alternates with mixing episodes. Lake Pampulha catchment, composed of 8 main tributaries, has a surface area of about 98 km 2 (Companhia de Pesquisa de Recursos Minerais, 2001), with more than 70% urbanised (Matos et al., 2017) and approximately 500,000 inhabitants. Lake Pampulha was built for many purposes, including drinking water supply. However, due to an intense urbanisation process in the lake catchment from the 1970s on, the lake water quality was severely compromised resulting in a hypereutrophic reservoir, leading to the interruption of drinking water provision from the 1980s, with frequent cyanobacteria blooms (Silva, 2014) and presenting a hypoxia condition at the bottom of the lake (Friese et al., 2010;von Sperling & Campos, 2011).
In 2016, the lake has become a UNESCO World Heritage site, which boosted public sector to seek for water quality improvement actions. The municipality of Belo Horizonte (PBH) committed to invest about R$ 30 million to recover the water quality of Lake Pampulha over 12 months, starting from March 2016 (Prefeitura de Belo Horizonte, 2016). Phoslock® and Enzilimp® were applied in the lake in order to remove phosphorus from the water column and to reduce coliforms and the Biochemical Oxygen Demand -BOD (Barçante et al., 2020). The treatment was renewed for another year at the cost of about R$ 16 million. Efforts to improve sanitation infrastructure in the Pampulha catchment were also carried out aiming at collecting 95% of sewage at an investment of R$ 875 million from 2004 to 2016 (Prefeitura de Belo Horizonte, 2017).

Available data
Concerning the lake water inflow, since 2011 monitoring stations integrate a Hydrological Monitoring System (Figure 1) operated by PBH (Siqueira et al., 2019). Monitoring stations are present only on the two main tributaries, Sarandi stream (station SAR18F, 2.8 km upstream the lake), located at 19º 52' 01" south and 44º 0' 33" west and Ressaca stream (station RES17F, 2,0 km upstream the lake), located at 19º 52' 11" south and 43º 59' 57" west. These streams represent 70% of the reservoir inflow (Companhia de Pesquisa de Recursos Minerais, 2001). These stations measure water level. Nogueira (2015) determined the rating curves, i.e. the water level x discharge relationship, for both streams at the monitoring stations. In this study, the rating curves were used to obtain the inflow data from water level measurements during the simulation period. The other tributaries are not monitored. Thus, time series of the inflow discharges were estimated by multiplying Thermal functioning of a tropical reservoir assessed through three-dimensional modelling and high-frequency monitoring the inflow time series of the reference stations by the ratio of each tributary catchment area divided by the catchment area of the reference stations.
Inflow water temperature was monitored for a short period of two months (April and May 2013) downstream the confluence of Ressaca and Sarandi streams (Silva et al., 2016). In this research we used hourly water inflow temperature. Two methods to estimate hourly water inflow temperature were evaluated. The first one was a linear regression with air temperature, used in previous studies (Chen & Fang, 2015;Hébert et al., 2015;Stefan & Preud'homme, 1993). The second method was the Modified Sine and Sinusoidal Wave Functions Model (MSSWF), proposed by Chen & Fang (2016).
Hourly meteorological data (air temperature, air humidity, wind intensity and direction, and cumulative hourly solar radiation and precipitation) were provided by the National Institute of Meteorology of Brazil (INMET) and used as uniform values in the computation domain. The nearest automatic weather station (a521) is located inside the campus of the Universidade Federal de Minas Gerais at 19º 53' 02" south and 43º 58' 10" west, about 3 km from Lake Pampulha ( Figure 1). The cloud coverage data are recorded at the meteorological station named Belo Horizonte (83587), located at 19º 56' 04" south and 43º 59' 43" west, 9.5 km from Lake Pampulha (Figure 1). This station provides nebulosity data three times a day (0h, 12h, and 18h) and these values were linearly interpolated to hourly values for simulation.
The most recent Lake Pampulha bathymetry was measured by PBH in November 2014. Based on the bathymetry data (Figure 2), the reservoir can be divided into four regions: (i) the upstream channel of Ressaca and Sarandi streams, (ii) an upstream very shallow region (about 1.5 m depth), (iii) a transition region, where the depth varies from 1.5 m to almost 13 m, and (iv) a deeper region further downstream. The inflow channel of the Ressaca and Sarandi streams was represented as thin dams in the model. In Delft3D-Flow, thin dams are infinitely thin and defined on the edges of elements. They prevent flow exchange between the two contiguous computational cells keeping the total wet surface and the volume of the model (Deltares, 2014). At station P 1 (10 m depth, Figure 1), located at 19º 51'5" south and 43º 58' 35" west, high-frequency water temperature monitoring (Table 3) was performed at the surface (0.5 m) at a time step of 30 minutes and at 2.5 m, 5.5 m and 9.5 m depth at a

Hydrodynamic modelling
The three-dimensional hydrodynamic model Delft3D-Flow was used. It solves in three dimensions the Reynolds Averaged Navier-Stokes (RANS) equations for an incompressible fluid, under the shallow water and the Boussinesq assumptions. In the vertical, a hydrostatic pressure equation is considered, and the velocities are computed from the continuity equation. The k−ε model, a second-order turbulence closure model was chosen. The model computes the horizontal and vertical turbulent coefficient by transport equation for the turbulent kinetic energy and kinetic energy dissipation based in Kolmogorov and Prandtl concept. A full description can be found in Deltares (2014) section 9.4.2.
Hourly data of water temperature were used to calibrate the parameters of the heat exchange module included in the three-dimensional model. Thus, the heat exchange main equations (Deltares, 2014), will be presented. In the model, the heat exchange module is named as "Ocean". The heat exchange at the free surface is considered as resulting from the separate effects of solar (short wave) and atmospheric (longwave) radiation, and heat loss, represented in terms of back radiation, evaporation and convection. The total heat flux ( The net incident solar radiation flux ( sn Q ) that will reach water surface is calculated as a function of short wave radiation for a clear sky condition, in which a part is reflected by the water surface (surface Albedo) and the magnitude is reduced by cloud cover (Equation 2).
Where a is the albedo reflection coefficient (dimensionless), sc Q is the solar radiation for a clear sky condition (J/m 2 .s) and C F is the sky clouds that cover the sky (%). The albedo coefficient used was 0,06 (default value) and the time series of solar radiation and cloud cover are input data.
Concerning the long wave radiations, the difference between the back radiation ( ) where S T is the water surface temperature (ºC) and a e is the vapour pressure (mbar).
Not all the radiation is absorbed at the water surface. A portion is transmitted to deeper water. The absorption of radiation in the water column is calculated by an exponential decay function (Lambert-Beer law). The evaporative heat flux ( ) ev Q is composed by a forced convection term (Qev,forced), caused by wind, and by a free convection term (Qev ,free ), due to buoyant forces from density differences between air and water (Equation 4). , , ev ev forced ev free The forced convection of latent heat is calculated according Equation 5.
where L v is the latent heat of vaporisation (J/kg), ρ a is the air density (kg/m 3 ), q s (dimensionless) and q a (dimensionless) are the specific humidity of respectively saturated air and remote air (10 m above water level), s T is the water surface temperature (ºC), a T is the air temperature (ºC) and ( ) where e c is the Dalton number (dimensionless) and 10 U is the wind velocity (m/s) measured at 10 m height.
The free convection of latent heat is calculated according Equation 7.
The convective heat flux (sensible heat), is also composed by a forced (Equation 8) and free term (Equation 9).
Where a ρ is the air density (kg/m 3 ), p c is the specific heat of air (0.24 cal/(g ºC)), H c is the Stanton number (dimensionless) and 10 U is the wind velocity (m/s) measured at 10 m height.
where g is the acceleration of gravity (m/s 2 ), air ϑ is the air viscosity (m 2 /s), 0 a ρ is the saturated air density (kg/m 3 ), 10 a ρ is the remote air density (10 m above water level) (kg/m 3 ), a ρ is the average air density (kg/m 3 ), p c is the specific heat of air (J/(kgC)). Wind intensity and direction measurements, used as uniform in the computation domain (spatial variations of wind speed were not considered), are carried out 10 m above the ground. Therefore, for the wind at lake surface ( U Lake surface ) an evaluation of a correction factor (k) of wind intensity was applied (Equation 10). Wind factor values were tested in a range of 40% higher and 40% lower than measured values, varying by 20% for each simulation. According Jöhnk et al. (2008), due to lake geometry and terrain roughness 6/25 Thermal functioning of a tropical reservoir assessed through three-dimensional modelling and high-frequency monitoring between the lake and the meteorologic station the wind intensity required to be calibrated. Dalton coefficient, Stanton coefficient, horizontal viscosity and diffusivity were also estimated during the model calibration.

Model setup
The model requires a computation domain and the specification of boundary and initial conditions. Lake Pampulha was discretized in ten vertical layers and 4019 horizontal grid elements (40190 elements in total). The mesh size varies from 50x20 to 10x10 m ( Figure 2). Near the dam and Ressaca and Sarandi inflow the mesh was refined to account for likely higher water current velocities due to inflow and outflow regions.
Bathymetry of Lake Pampulha ( Figure 2) greatly varies from upstream to downstream. During calibration of the hydrodynamic model the σ-grid and Z-grid options were used to compare the model performance and the influence of the grid type. The Z-grid option does not change the thickness of each vertical layers, which may help to represent the thermocline depth and water stratification. The σ-grid option has the advantage to preserve the number of active layers along the whole grid. Thus, in the shallower part of the lake the vertical resolution was higher (0.15 cm). At monitoring point P1, the Z-grid and σ-grid were very close, with 1-m thickness vertical layers (Figure 3). The σ-grid option allows to analyse the hydrodynamics at the lake bottom in a same layer in continuity. Otherwise, using the Z-grid option the bottom of the lake would be represented in different layers.
Outflow was concentrated in the free spillway located at the dam ( Figure 2). The open downstream boundary was set as the spillway rating curve (Felisberto et al., 2015). Thus, the outflow discharge is defined according to the water level in the spillway. The simulated periods were performed during the dry season, with low inflow (maximum inflow was smaller than 2-year return period flow), representing a water level of about 10 cm higher than the spillway crest. Thus, the initial water level of the lake was considered as the same level of the spillway for all simulated periods.
The initial conditions for the reservoir water temperature at the beginning of each simulation period were obtained from the lake monitoring station. The simulations were initialized with vertically uniform temperature. For velocities, the water was supposed to be at rest as initial condition. For all simulations, the following constants were used: acceleration of gravity at 9.81 m 2 /s, water density at 1000 kg/m 3 , air density at 1 kg/m 3 and water salinity at 0.15 ppt (calculated with the mean of field conductivity data from 2015 to 2016 in surface water using the PSS-78 (Practical Salinity Scale of 1978 or Sal78 function) formulated and adopted by UNESCO (1981). For bottom roughness, the Chezy formulation was chosen with uniform values of 65 m 1/2 .s -1 with free slip condition in wall roughness.
According to Deltares (2014), turbulent viscosity and turbulent diffusion are calculated differently in the horizontal and vertical directions. In turbulent flows, the diffusion mechanism is not only performed by molecular motion, but also intensified by eddy motion. Thus, the diffusion coefficient is computed by adding turbulent viscosity to molecular. Vertical terms are very small compared to horizontal terms.
The Delft3D-Flow model uses background horizontal eddy viscosity and diffusivity coefficients to represent unresolved hydrodynamic processes not calculated by the turbulent model.
Measured Secchi depth presented almost no variation (between 0.20 and 0.35 m, with a mean of 0.30 m). Secchi depth was set as 0.3 m.

Simulation periods
The simulation periods were selected according to the following criteria: (a) vertical uniform temperature of the water  The calibration period lasted 18 days, from May 16 th to June 3 rd 2016, corresponding to 440 hourly data of water temperature for each of the four depths (0.5, 2.5, 5.5 and 9.5 m). The water temperature time series measured during this period at four depths at station P1, which has a maximum depth of 10 m, is plotted in Figure 4. The period started with a mixed thermal condition of about 22.9 ºC. This temperature value was set as the initial temperature condition for the lake simulation for all grid elements.
The validation period lasted 16 days, from May 29 th to June 14 th 2016, corresponding to 388 hourly values of water temperature ( Figure 4). The period started with a mixed thermal condition of about 22.2 °C. Tributary inflows during the calibration and validation period are presented in Figure 5. A third period, from May 15 th to August 10 th 2015 (90 days, 2107 hourly values), was simulated using the calibrated parameters in order to evaluate the model performance on a longer period. For this period, only subsurface temperature data is available (Figure 6). Tributary inflows during the period are presented in Figure 7. The meteorological variables, respectively in 2016 (calibration and validation) and 2015 (long validation period) are presented in Figure 8 and Figure 9. Figure 10 presents the box plot with minimum, first quartile, median, third quartile, and maximum values for each meteorological variable. During calibration and validation periods no precipitation were measured at the meteorological station near the lake (station A521, Figure 1). During the period of 2015 no significant precipitation was measured (maximum of 14.2 mm/h on 27 th July). The effect of precipitation on the water temperature is not taken into account by the model. According to Figure 11, during the simulated periods, east wind direction is dominant on the reservoir.  Thermal functioning of a tropical reservoir assessed through three-dimensional modelling and high-frequency monitoring

Performance indicators
The hydrodynamic model performance during calibration and validation periods was assessed through comparison between measured and simulated water temperature. Four mathematical indicators were used to compare simulated and measured water temperature: the mean absolute error (MAE), the root mean square error (RMSE), the coefficient of determination (R 2 ) and the relative error (RE). These indicators were computed using hourly model outputs and hourly-averaged measurements.
Regarding the inflow water temperature estimation, the Nash-Sutcliffe efficiency coefficient, NSE (Equation 11) was used for assessing the goodness of fit. The NSE shows how well a model represents observed data variance compared to the mean value of observed data. NSE varies from -∞ to 1, the closer to the unity, the better the model performance. A NSE below zero indicates model estimation worse than the mean value of the observed data (Bennett et al., 2013).
where n is the number of observations and the number model results, y i is the ith observation, y i is the mean value of observations and y' i is the ith model result.
To analyse the impact of the stream inflow in the lake thermal stratification, Schmidt stability index and Lake number were computed using Lake Analyzer (Read et al. (2011).

Inflow water temperature
The hourly water inflow temperature was estimated from air temperature using linear regression and MSSWF method. Root mean square error (RMSE), the coefficient of determination (R 2 ) and the Nash-Sutcliffe efficiency (NSE) were used to compare the results of linear regression and of the MSSWF method (Chen & Fang, 2016). Compared to MSSWF, linear regression presented better values for all indicators (Table 5). The hourly inflow temperature was computed using a linear regression with air temperature (Equation 12).
where WT is the water temperature (ºC) and AirT is the air temperature (ºC).

Model calibration
In the calibration step, wind factor values were tested for a range from 40% higher to 40% lower than 1, varying by 20%. For all depths, the better adjustment between measured and simulated temperature values were obtained for the lower wind factor (0.6), however with no significant differences ( Table 6). At the water surface, the simulated values presented a higher variation compared to the measured values for all wind factors. For the 2.5, 5.5 and 9.5 m depths, a systematic error occurs after May 25 th 2016. The simulated water temperature is colder than the measured values.
Setting wind intensity equal to the measured values (wind intensity factor equal to 1), the Dalton and Stanton coefficients   The Z-grid option was also used to compare the influence of the grid type on the model performance. Using the Z-grid option, a non-homogenous vertical distribution was necessary to provide more active layers in the lake upstream region, avoiding numerical instabilities due to inflows. The results for the Z-grid and σ-grid were interpolated to the measuring depth.
The mathematical indicators are very similar (Table 7). The indicators using σ-grid are slightly better in surface (0.5 m) and 2.5 m depth for RE and MAE. Thus, the σ-grid was chosen.
The final parameter values are shown in Table 8. The calibration period started on May 16 th 2016 with a uniform water temperature of 22.9 °C. Over the following ten days, the lake water column stratified with warming in the surface layer to around 25.3°C. Mixing events (temperature difference between surface and bottom less than 0.25 °C) occurred on May 25 th and 26 th 2016, with temperature decreasing to 22.6 °C. After May 26 th , the temperature at 9.5 m depth decreased by almost 1º C and the surface temperature displayed daily warming. A third stratification period lasted three days, ending with a mixing event on May 29 th at 22.1 °C. During the whole period, the temperature difference between the surface and 9.5 m depth, at maximum stratification, was 2.9 °C (May 18 th at 15:00).
Observed and simulated water temperatures at the four depths during the calibration period are presented in Figure 12. Figure 13 shows the difference between surface and bottom layer for measured and simulated water temperature. Figure 14 shows the temperature for the four depths in which measurements were carried out. Table 9 shows the number of days with a difference between surface (0.5 m) and bottom (9.5 m) for measurement and simulated temperatures of 0.25 °C to 2.0 °C. According to both Figure 13 and Figure 14 and to Table 9, the model could accurately represent the main observed patterns: (i) the daily time scale of surface warming; (ii) the maximum stratification time, (iii) the cooling at the lake bottom on May 26 th ; (iv) the timing

Validation Period in 2016
The model was run with the calibrated parameters values (Table 8). For all depths, the RMSE and RE between the measured and simulated values were good (Table 10 and Figure 15), with maximum values of 0.45 °C and 1.51%, respectively. The R 2 ranged from 0.41 to 0.81, also considered as a good result. The MAE presented a very good performance, with a maximum value of 0.34 ºC at 0.5 m. According to the mathematical indicators, the model can be considered validated.
The validation period in 2016 started on May 29 th with a uniform water temperature of 22.2 °C. Over the first five days a daily stratification and mixing alternation was observed. According to Figure 16, on May 30 th the difference between surface and bottom temperature was 0.35 °C for measured values and 0.07 for simulated values. On May 31 st , this difference was 0.53 and 0.31 ºC, on 1 st June it was 0.77 and 0.44 ºC and 2 nd June 0.11 and 0.46 °C, respectively for measured and simulated values. Until 5 th June, comparing with the measured values, the simulated water surface presented a good adjustment, however, after this date the simulated values presented higher amplitude with lower low values (Figure 15). This discrepancy resulted in higher differences between simulated surface and bottom water temperature ( Figure 16).
On 3 rd June, the water temperature column was uniform around 22.0 °C. Then, four hours later, at 09:00 a strong cooling at 5.5 and 9.5 m depth occurred. At 9.5 m depth, the water temperature decreased by 1.8 °C over 4 hours. At 5.5 m depth, the decrease was 0.8 ºC over 12 hours. After this cooling of the bottom layers, the water at 5.5 m depth remained at a relatively constant temperature of 21.5 °C for the whole period. In contrast,     warming was observed at 9.5 m depth at a constant rate of 0.08 °C per day. At 2.5 m depth, the temperature remained constant until 7 th June when it began to increase. This temperature dynamics could also be represented by the simulation (Figure 17).

Validation Period in 2015
A longer period of 88 days (from May 15 th to August 10 th 2015) was run with the calibrated parameters. For this period, only surface temperature data are available. The simulation results and data measurements during the period are presented in Figure 18. Even with a slightly warmer simulated water temperature, the MAE between the measured and simulated values was very good (0.45 °C), the RMSE (0.65 °C), the RE (2.08%) and the R 2 (0.55) also presented good values.

Impact of the river inflow
On June 3 rd 2016, during the validation period, the water temperature strongly decreased by 1.8 °C over 4 hours at 9.5 m and by 0.8 °C over 12 hours at 5.5 m (Figure 15 and Figure 17) at station P1. This strong temperature variation observed only in the deepest layers was not expected. The model was used as a tool to investigate what could be the reason for this cooling observed only below 5.5 m depth and more intense at 9.5 m depth. According to the computed water inflow temperature (Figure 19), the tributary water, just before the water cooling, presented a lower temperature (16.1 °C) and, therefore, a higher density, compared to the lake temperature (uniform at 22 ºC). This colder inflow entered the reservoir during 15 hours ( Figure 5), with a high discharge of 50 m 3 /s of the streams Ressaca and Sarandi. On 1 st June at 17 h and 2 nd June at 10 h other inflow high discharge events were observed, however with lower discharge (25.0 m 3 /s of Ressaca and Sarandi streams) and smaller duration (one to two hours). The inflow water temperature, respectively for June 1 st and 2 nd was 22.0 ºC and 16.2 °C. Even though these events were secondary, compared to the event of June 3 rd , they occurred in different condition of the water column stratification. On June 1 st , the lake was stratified, with water surface at 23.9 ºC and bottom at 21.9 °C. On June 2 nd the difference of water surface and bottom was less than 0.50 ºC, with water surface at 22.3 ºC and bottom at 21.85 °C, therefore, the lake was warmer that the water inflow.
During the peaks of water inflow on June 1 st and 2 nd , no decrease of temperature was measured and simulated at station P1, probably due to a lower inflow discharge and duration that was dissipated before reaching point P1. At the upstream very shallow region, on June 1 st , no water temperature decrease was observed. However, the water inflow on 2 nd June resulted in a water temperature decrease of 1.7 ºC at 2.0 m depth (reaching 17.6 ºC), 1.6 ºC at 1.1 m depth (reaching 19.0 ºC) and 0.6 ºC at 0.5 m depth (reaching 20.7 ºC). Figure 20 shows the water temperature profiles at a point in the shallow region indicated in Figure 21 as UP point. Even being lower, the event on June 2 nd impacted the bottom temperature that remained colder until a higher inflow arrived on June 3 rd , in which water temperature reached 16 ºC at this point. This water temperature dynamics, during June 1 st and 2 nd , shows that the first water high inflow was dissipated fast and the second one impacted a wider region of the shallow part of the lake. The colder water remained in the shallow region until the high discharge arrived on June 3 rd . In the supplementary material videos show the lake temperature at surface (Video 01) and bottom (Video 02) layers of the lake during the validation period. Videos also show the impact of the Tijuco and Mergulhão tributaries (Figure 2) that were not so intense due to a lower inflow rate ( Figure 5).
Even with higher inflow, the simulated water level variations at station P1 were less than 12 cm during the passage of the higher inflow values (Figure 22). This small water variation is in agreement with the downstream condition, in which a water elevation of 0.5 m would spill 121 m 3 /s, much higher that 50 m 3 /s, the maximum inflow discharge during the simulated period.
On June 3 rd 2016 at 13:00 hours, when the water velocity module in the bottom layer at P1 station is maximum, a  (Figure 21) of the simulated water temperature shows colder temperature in the lake deeper layers (Figure 23). The model suggests that the propagation of the inflow water on 3 rd June occurred through the bottom of the lake (Figure 23 and Figure 24). The inflow water sank to the deepest layers of the lake. The plunging point is located 1500 m far from the inflow channel, where the lake shallow region ends. In the supplementary material a video shows the lake temperature evolution along the longitudinal transect (Video 03). Unfortunately, the stream inflow water temperature was not measured in this period to confirm this hypothesis. Therefore, to evaluate the influence of the inflow water temperature and to assess the impact of different intensity of inflow water conditions, two simulations considering the water inflow at a constant temperature of 22 ºC were performed for two different scenarios (Figure 19). In the first scenario (D1), the water inflow temperature was considered constant (22º) on June 3 rd . In the second scenario (D2), the water inflow temperature was considered constant (22º) on June 2 nd to 3 rd . In these simulations the main difference in the temperature profile was on the bottom of the lake ( Figure 24) and more evident in the scenario considering the water inflow at a constant temperature during two days (D2). In the simulations with the computed water inflow temperature (validation period) the lake bottom temperature presented a decrease of 1.8 ºC at station P1.
In scenario D1, with the inflow water temperature of 22.0 ºC, the decrease was 0.8 ºC. Changing the inflow water temperature to 20.0 ºC, during 3 rd June, the decrease was 1.4 ºC, and, changing it to 18.0 ºC the decrease was the same of the validated simulation (1.8 ºC). Therefore, the variation of the inflow water temperature  directly impacts the water stratification. For the scenario D2, with the inflow water temperature of 22.0 ºC, the decrease of the lake bottom temperature was only 0.2 ºC. The difference of water temperature decreases between scenario D1 and D2, with the inflow water temperature of 22.0 ºC, shows that the event of 2 nd June also influences the water temperature at station P1, even if indirectly. This probably happened because the shallow region of the lake remained at a lower temperature until the inflow peak arrives on 3 rd June. Concerning the time evolution of the water velocity module, as expected, the higher values occurred on 3 rd June, and due to the water sank to the deeper layers of the lake the higher values were at the bottom. At the moment of higher decrease of the water temperature, the bottom velocity module was 7 cm/s and the surface velocity was 1.5 cm/s. In the scenario D1, the bottom velocity module was 5.8 cm/s and the surface velocity 1 cm/s. Changing the inflow water temperature to 20 ºC and 18 ºC the water velocity module presented a very slight increase in almost all depths with a higher increase at 9.5 m depth (Table 11). In the scenario D2, the bottom velocity was 4.5 cm/s and the surface velocity 0.5 cm/s. Table 11 presents the water velocity module in different depths for the presented scenarios.
Close to the spillway, according to the model results, the lake presented a water temperature decrease with a higher intensity at the bottom ( Figure 25). This water cooling happened 6 hours later that the cooling observed at station P1 and with a lower intensity (decrease of 1.1 ºC). Analysing the profile of the water velocity module, the highest value was at bottom (6.1 cm/s). The simulated surface velocity module was 2.5 cm/s. This point (DOWN in Figure 21) is located about 3.0 km far from the Ressaca and Sarandi inflow channel. Figure 26 and Figure 27 present the profiles of the water velocity module for the UP and DOWN points ( Figure 21). According to the simulation values, the dynamics of the upstream region of the lake is different of the downstream region. At the upstream region the lake presents higher velocity module values at surface. At the downstream region the bottom layers present higher velocity values. In the supplementary material videos show the lake velocity module at surface (Video 04) and bottom (Video 05).
Schmidt Stability represents the resistance to mechanical mixing (Idso, 1973). It corresponds to the energy required to completely mix a stratified lake. For the validation period, Schmidt stability shows a pattern of day-night oscillations ( Figure 28), induced by the heating and cooling alternation period. It corresponds to the same pattern of the water temperature difference between surface and bottom depths ( Figure 16). However, just before of the cooling of the bottom depth on June 3 rd (Figure 4 and Figure 17) a mixing of the water column occurs that strongly reduces the value of Schmidt Stability index. This mixing event reduced the daily variation on June 3 rd , which shows a smaller amplitude. Before the    mixing event, the daily variation was around 14 J/m 2 , on 3 rd June it reduced to 4 J/m 2 . Figure 28 also shows that the model could represent this Schmidt stability evolution. However, after June 5 th , the agreement between measurements and simulated results is weaker (Figure 15), with larger differences between surface and bottom depth (Figure 16), and directly impacting the Schmidt stability index computed from simulated results. The Lake Number, defined by Imberger & Patterson (1990), is the ratio of the strength of stratification (Schmidt Stability) to the effect of the wind stress. Lake number below 1 means that the lake stratification is weak with respect to wind stress and the lake will probably mix. According to Figure 28 and Figure 29, the Lake Number peaks were coincident with the Schmidt Stability peaks for measured and simulated values. The high inflow discharge on 3 rd and 4 th June resulted in the cooling of the bottom of the lake. On 5 th June, after the cooling of the bottom depth the water difference between surface and bottom depth increases (Figure 16), which increases the Lake Number to values higher than 600 J/m 2 .   Thus, the stream inflow impacted the lake stratification strength and the three-dimensional model could represent this pattern.

DISCUSSION
A literature review of papers published from 2000 to 2015, carried out by Meinson et al. (2016), showed that more than two-thirds of the studies that used high-frequency monitoring in lakes were carried out in northern temperate zone, either in North America or in Europe. The authors, also showed that water temperature was the most common measured parameter, because of the sensor reasonable price and low maintenance and the most important, because water temperature is a controlling factor of biological, ecological and chemical processes (Meinson et al., 2016).
In our paper, a three-dimensional hydrodynamic model, Delft3D-FLOW, was implemented in Lake Pampulha, a shallow urban tropical reservoir. The calibration period lasted 18 days, from May 16 th to May 29 th 2016, corresponding to 440 hourly data at 0.5 m, 2.5 m, 5.5 m and 9.5 m depth. The model could accurately reproduce the alternation of stratification and mixing conditions along the period. The validation period lasted 16 days, from May 29 th to June 14 th 2016, corresponding to 388 hourly data at the same depths of calibration period. Performance indicators showed that the model performed well. The model was also validated for a longer period from May 15 th to August 10 th 2015, corresponding to 2107 hourly data at 0.5 m depth.
Mean absolute error (MAE) during all simulation periods for water temperature varied between 0.15 and 0.45 ºC and was similar to other three-dimensional model studies. Soulignac et al. (2017) obtained values in a range between 0.25 and 2.34 ºC for Lake Créteil (surface area of 0.4 km 2 and mean depth of 4.5 m,) using hourly temperature values at five depths (0.5, 1.5, 2.5, 3.5 and 4.5 m) over 8 months (around 5500 values). In Lake Yilong (surface area of 28.4 km 2 and mean depth of 3.9 m), Zhao et al. (2013) obtained a R 2 between 0.65 and 0.74 calculated from measured data of 15 campaigns throughout one year, in the same range of our values (0.55 to 0.93). In Lake Minnetonka (maximum depth of 25.0 m and surface area of 8.01 km 2 ), Missaghi & Hondzo Thermal functioning of a tropical reservoir assessed through three-dimensional modelling and high-frequency monitoring (2010) obtained a R 2 between 0.91 and 0.98 using 6 months of bi-weekly measured profiles of temperature for calibration and validation period. Belico (2017) investigated the thermal regime of Lake Pampulha during rain events from 2011 to 2015 with a onedimensional model using hourly measurements of water temperature at surface and monthly vertical profiles. The unidimensional modelling provided a Root Mean Square Error (RMSE) of 1.08 ºC at surface. The authors highlighted that the model could not rightly simulated the lake bottom layers. Other studies were also applied to Lake Pampulha (Silva, 2014;Silva et al., 2016Silva et al., , 2019 however, none of them used three-dimensional modelling with high-frequency measurement in different depths.
In the present paper, high-frequency measurement could detect sudden changes of water temperature with different amplitudes depending on the depth.
The water temperature measurements in the centre of the lake were performed at surface (0.5 m), 2.5 m, 5.5 m and bottom (9.5 m) and were compared to three-dimensional hydrodynamic model results, allowing a thorough analysis of the lake hydrodynamics. In our paper, the RMSE varied between 0.15 and 0.65 ºC and highfrequency measurement coupled with three-dimensional model could appropriately represent high vertical temperature gradient.
The good simulation results indicate that the meteorological data, collected at a weather station located 3 km from the lake (and from 9.5 km for the cloud cover data) reproduce well the meteorological conditions over the lake most of the time. However, the performance of the model may have been affected due to high sensibility of the lake hydrodynamics to external forcing such as nebulosity and wind intensity and direction. Another source of uncertainties was the inflow water temperature that was estimated using air temperature.
Using high-frequency time series, a good model performance was obtained using short simulation period, which is suitable in terms of computation time demand. However, the simulations were not performed for all seasons. We focused on the dry season because it is the period with higher episodes of cyanobacteria blooms. The model must be validated over longer periods representing other hydro-meteorological conditions. Despite using data from only one monitoring station to calibrate the model and despite the limitation in water inflow temperatures measurements, the model was useful to understand and evaluate the impact of forcing scenarios. During high discharge events, when the river temperature is colder than the lake water, it flows through the deepest layers of the lake. In the simulated episode, the inflow water sank to the bottom of the lake, 1500 m far from the inflow channel. The model showed that the colder water underflow in the deepest part of the lake is still visible almost 3.5 km far from the river entrance. The 3D hydrodynamic model could be further used for investigating in other hydrometeorological conditions the complex spatio-temporal dynamics of Lake Pampulha or similar reservoirs.
Lake Pampulha water quality is severely impaired, with recurring hypoxia in the bottom layers (Friese et al., 2010;von Sperling and Campos, 2011). The water current caused by the inflow, through the deepest layers of the lake, may help to renew and oxygenate the water. It also may impact the lake ecological status by affecting the sediment resuspension and nutrient release. Thus, providing a modelling tool which can help to investigate the physical conditions of hypolimnion aeration is an important outcome of this study.
Due to the importance of water stratification and hydrodynamics on phytoplankton growth, our results may also help to understanding algal blooms. The good results of the hydrodynamic model allow, for the next steps, to study phytoplankton and nutrient dynamics by coupling the water quality module.

CONCLUSION
In this paper, a three-dimensional hydrodynamic model, Delft3D-FLOW, was implemented in Lake Pampulha, a shallow urban tropical reservoir. For the three simulated periods, the agreement between the modelled and hourly measured temperatures were considered very satisfactory.
This paper highlights the complex hydrodynamics of a shallow reservoir and the importance in monitoring the inflow water temperature and discharge. The model outcomes made it possible to investigate the cause of sudden changes in the water temperature. Future studies may focus on extending the simulation period to include rainy season and obtaining data in other monitoring stations in the lake to provide a further assessment of the 3D hydrodynamic model performance.
Investigating the inflow impact on the lake hydrodynamics may contribute to assess the efficiency of restoration techniques, supporting the decision-making for water management of urban lakes.