Thermal Analysis of a Concrete Dam Taking into Account Insolation, Shading, Water Level and Spillover

Featured Application: Seasonal temperature distributions in concrete dams affect thermal loads, which can cause high stresses in concrete and signiﬁcantly affect the occurrence and behavior of cracks, as well as displacements (especially in arch-gravity and arch dams). Abstract: This study presents a procedure for modeling the heat transfer process in a concrete dam, taking into account the time-varying boundary conditions on the upstream and downstream sides of the dam (i.e., the water level of the reservoir, spillover, insolation, and shading) which affect the temperature conditions of the dam. The large concrete arch-gravity Moste Dam (in North West Slovenia) was analyzed, where an automated system for the measurement of concrete and water temperatures, and for the monitoring of meteorological effects, was installed. Thermal analyses (1D and 2D) for non-linear and non-stationary heat conduction through solids were performed using a ﬁnite element method (FEM) based program, TeEx, which was complemented by two specially developed programs for determining the effects of convection and insolation, considering also the effect of shading. A 15-day period was analyzed, as well as a period of one year. It was found that the results of the performed analyses ﬁtted in well with the experimentally determined concrete temperature measurements. The results showed that at the insolated side of the dam, the temperature gradient was largest in a very narrow area along the concrete surface, but the temperature did not stabilize shallower than at a depth of about 6 m. As part of the thermal analyses, uncertainty analyses of the results of the calculations were also performed.


Introduction
In recent times, the health monitoring of concrete dams has become a topic of great importance and involves monitoring the static and dynamic behavior of large dams [1]. One of the most important parameters that strongly influence the static behavior of concrete dams, especially arch dams, is temperature. For this reason, it is very important to know the seasonal temperature distributions in concrete dams, as they affect the thermal loads, which can cause high stresses in concrete dam structures and significantly affect the occurrence and behavior of cracks, as well as displacements. The effects of climate change and global warming, which may result in a general temperature increase in these structures, need to be considered [2]. Early fundamental research in the field of thermal analysis related to nonlinear and nonstationary heat conduction through a solid, considering suitable boundary conditions, was carried out by Dilger et al. [3] and Carslaw and Jaeger [4], whereas some years later Leger et al. [5,6] presented a methodology based on the finite element method (FEM), which could be used to determine seasonal temperature and stress distributions in concrete gravity dams. Using this methodology, the thermal behavior of concrete dam structures (from arch to gravity) has been much better analyzed, both in the phase of dam design [7,8] as well as in the later phase of monitoring and surveillance of the dam [9,10]. It has been discussed in an increasing number of articles dealing with: stress analysis of concrete structures, which are affected by variable thermal loads [11]; the determination of the periodic temperature field in a concrete dam [12]; the numerical calculation of temperature in a mass concrete [13]; thermal and stress analysis of the discussed roller-compacted concrete dam [14]; and the effect of environmental impacts on thermal stress analysis of the discussed arch concrete dam [15]. Some other authors have discussed the development of new numerical models, such as the model for the analysis of concrete dams due to environmental thermal effects [16], hydrostatic and temperature time-displacement model for concrete dams [17] and hydrostatic seasonal state model for monitoring data analysis of concrete dams [18]. Some others have reported about the sensitivity study of thermal fields in massive structures [19], and about the thermal protection of concrete dams in northern regions [20].
In recent years a large amount of research has been performed in connection with the thermal performance of concrete dams, in which the effects of the changing surrounding conditions were taken into account. Some researchers have focused on detailed studies of displacements dealing with: analysis of concrete dams based on seasonal hydrostatic loading [21]; seasonal thermal displacements of gravity dams, which are located in cold regions [22]; the method of dam deformation due to thermal actions [23]; and correlation between the air temperature and the daily variation of structural response [24]. Other researchers have discussed cracks in the dam structures [25,26], whereas others have analyzed temperature-induced stresses in concrete dams [27,28]. Others again have studied shading effects [29,30] and the effects of solar radiation [31,32], and some reports have discussed temperature fields in the case of super-high arch dams [33][34][35][36].
Some other researchers have recently compared the field monitoring data with the results of numerical analyses of temperature distributions also for other structures (Xia et al. [37] for long-span suspension bridges, Su et al. [38] for supertall structures) and for specific materials (Abid et al. [39] for concrete-encased steel girders).
A literature review of the research on the dam monitoring data analysis methods, which are very important for monitoring the dam safety and where the influence of the environmental variables (e.g., water level, temperature) has been involved, has recently been presented by Li et al. [40].
However, in the above-mentioned articles there seems to be no exhaustive database about the results of temperature measurements over the several decades since the dam's original construction, or about changing thermal fields in the concrete due to daily and seasonal boundary condition changes.
In this article the authors present a relatively simple procedure which can be used to model the heat transfer process in a concrete dam, taking into account time-varying boundary conditions both on the upstream side of a dam, as the effect of oscillation of the water level of the reservoir, as well as on the downstream, insolated, side, such as the effects of insolation, shading and water spilling over the spillway crest. This procedure was verified at the highest arch-gravity concrete dam in Slovenia.

Thermal Analysis of a Concrete Dam
In the case of a concrete dam with incident solar radiation using the effect of shading, a time-varying water level of the reservoir, and occasional spillover, the heat transfer process is shown in Figure 1.

•
Neumann or imposed heat flux (the defined surface heat flux, that is the absorption of solar radiation considering the effect of shading); • Robin or convective flux (the heat flux is linearly dependent on the temperature difference between the surface and the ambient temperature, e.g., heat exchange by air convection or by water convection); • Radiation from the surface of a solid (the heat flux is nonlinearly dependent on the temperature difference between the surface and the ambient temperature). Heat exchange by convection can be determined by the well-known equation presented by Žvanut et al. [41]: where qc is the convective flux (W/m 2 ); hc is the convection coefficient (W/m 2 K); T is the surface temperature of the dam (K); and Ta is the ambient temperature (K).
The effect of radiation was checked and found to be small, so it was neglected. The same was found regarding the effect of cooling of the concrete due to evaporation of absorbed water.
A description of a relatively simple procedure proposed by Dilger et al. [3] for the quite accurate determination of energy flux density due to the absorption of solar radiation is presented by Žvanut et al. [41], where the amount of solar energy absorbed by the body of a system is given by The corresponding equation which relates nonlinear and nonstationary heat conduction in the case of a two-dimensional space and a homogeneous isotropic solid whose thermal conductivity is independent of the temperature [4] is: where T is the temperature (K); x, y are the Cartesian coordinates (m); ρ is the density (kg/m 3 ); c is the specific heat (J/(kg K)); λ is the thermal conductivity (W/(m K)); and t is the time (s). The principal types of boundary conditions are: • Dirichlet or imposed temperature (the defined surface temperature, for instance of the foundation); • Neumann or imposed heat flux (the defined surface heat flux, that is the absorption of solar radiation considering the effect of shading); • Robin or convective flux (the heat flux is linearly dependent on the temperature difference between the surface and the ambient temperature, e.g., heat exchange by air convection or by water convection); • Radiation from the surface of a solid (the heat flux is nonlinearly dependent on the temperature difference between the surface and the ambient temperature).
Heat exchange by convection can be determined by the well-known equation presented by Žvanut et al. [41]: where q c is the convective flux (W/m 2 ); h c is the convection coefficient (W/m 2 K); T is the surface temperature of the dam (K); and T a is the ambient temperature (K). The effect of radiation was checked and found to be small, so it was neglected. The same was found regarding the effect of cooling of the concrete due to evaporation of absorbed water.
A description of a relatively simple procedure proposed by Dilger et al. [3] for the quite accurate determination of energy flux density due to the absorption of solar radiation is presented by Žvanut et al. [41], where the amount of solar energy absorbed by the body of a system is given by where q s is the energy flux density of solar radiation absorbed by the solid (W/m 2 ); I is the energy flux density of solar radiation at Earth's surface (W/m 2 ); a is the solar absorptivity of the surface; and θ is the angle of incidence of Sun's rays ( • ); where I sc is the solar constant (on average 1350 W/m 2 ); and k T is the transparency factor (depends on weather conditions in the atmosphere and the path length of the radiation through the atmosphere); where k A is factor of altitude H (where H is measured in metres); m is the factor of influence of the relative path length of the radiation (depending on the solar elevation angle β s ); and t U is the air pollution factor (between 1.8 for clean air and 9.0 for very polluted air); The angle of incidence of the Sun's rays can be presented by cosθ = −cosδ · cosΩ · sinα · cosτ · sinΦ + cosδ · sinΩ · sinα · sinτ + sinδ · cosΩ · sinα · cosΦ+ +cosδ · cosα · cosτ · cosΦ + sinδ · cosα · sinΦ and cosθ ≥ 0 where δ is the declination (i.e., the angle between the equatorial plane and the direction of the Sun; positive towards the North Celestial Pole); Ω is the azimuth of the normal to the plane (measured clockwise from north); α is the angle between the plane and the Earth's surface; Φ is the geographical latitude; and τ is the hour angle (positive in the morning), which is denoted by where u is the time of day (in 24-h notation; the hour angle changes at a rate of 15 • per hour).
In the case of the vertical plane (α = 90 • ) Equation (8) can be simplified into cosθ = −cosδ · cosΩ · cosτ · sinΦ + cosδ · sinΩ · sinτ + sinδ · cosΩ · cosΦ (10) whereas in the case of the horizontal plane (α = 0 • ) it can be simplified into where β s is the solar elevation angle. The solar azimuth angle (α s ) defines the direction of the Sun, assuming the northclockwise convention, and it can be defined by the expression The heat conduction Equation (1), with appropriately defined boundary conditions, can be solved numerically by the FEM [42]. Such discretization of the model into finite elements results in a system of nonlinear differential equations of the first order. In the case of the searched for time-dependent node temperatures, the system takes the following form: where K is the thermal conductivity matrix; T is the nodal temperature vector; C is the heat capacity matrix; T,t is the time derivative of temperature; and F is the vector of the external actions. The left-hand side of Equation (13) and the coefficients depend, in general, on the searched for temperatures (T), whereas the right-hand side is usually also an explicit function of time (t).

Fieldwork
The large concrete arch-gravity Moste Dam, built in 1952 on the Sava Dolinka River in the north-western part of Slovenia, was analyzed. It is the highest dam in Slovenia, with a structural height of 59.80 m. The dam crest has a length of 72.00 m, whereas the dam volume is 42,000 m 3 . The spillway crest is at an altitude of 523.50 m. However, the steel gates, which have a height of 1.25 m, increase the effective spillway crest to the altitude of 524.75 m. The dam is located at the north latitude of 46.41 • , whereas the azimuth of the symmetrical axis of the downstream side of the Moste Dam is 186 • , so that it practically faces south.
An automated system for the measurement of concrete and water temperatures, and for the monitoring of meteorological effects using a mobile automatic weather station (MAWS), was installed on the dam in July 2013. A system for taking automatic measurements of the water level of the reservoir and the height of the water spilling over the spillway crest had been established in 2000.
Before determining the locations of the concrete temperature gauges and the water temperature gauges, preliminary calculations and analyses of the results of the conduction of heat through the concrete structure were performed, as well as a preliminary analysis of the available data regarding the oscillation of the water level of the Moste reservoir for the past 30 years.
Three boreholes were drilled in the dam and six concrete temperature gauges were installed. In the upper borehole on the downstream side of the dam, which was made where the dam surface inclination is 30 • , four gauges (TC1 to TC4) were installed at different distances from the dam surface: TC1 at 0. The water temperature measurements were performed using four gauges (TW1 to TW4), placed in a reservoir in a steel pipe, which protected them against floating debris and rushing waters, located at altitudes 519 m, 521 m, 522 m and 523 m. The lowest gauge was installed at a depth that was not affected by the surface zone of changing water temperatures [43], which was confirmed by an additional gauge placed in the reservoir at altitude 511 m.
The measurements of meteorological effects included the following parameters: solar radiation, precipitation, air temperature, wind speed and wind direction. The location of the MAWS, at an altitude 527.65 m, was determined after a detailed visual inspection of the area near the dam, bearing in mind easy accessibility, as few topographical obstacles as possible and low noticeability.
An overview of the locations of the concrete and water temperature gauges as well as the MAWS is given in Figures 2-4. of the area near the dam, bearing in mind easy accessibility, as few topographical obstacles as possible and low noticeability. An overview of the locations of the concrete and water temperature gauges as well as the MAWS is given in Figures 2-4.   of the area near the dam, bearing in mind easy accessibility, as few topographical obstacles as possible and low noticeability.
An overview of the locations of the concrete and water temperature gauges as well as the MAWS is given in Figures 2-4.

General
The temperature conditions of the dam were determined by means of the procedure developed for modeling the heat transfer process in concrete dams, taking into account time-varying boundary conditions on both the upstream side of the dam (the effect of the oscillation of the water level of the reservoir) and on the downstream side (the effects of insolation, shading, and water spilling over the spillway crest). Computer programs for the determination of temperature during time, at each node of the finite element mesh, were developed within the programming environments Mathematica [44] and MATLAB [45]. The thermal properties of the mass concrete used in the determining of the temperature conditions of the dam are given in Table 1.

Effect of Spillover
Spilling of water from the reservoir over the spillway crest has an effect on the determination of the parameters of convection (i.e., the temperature and the convection coefficient) and of the heat flux due to insolation, taking into account the effect of shading, on a selected surface, at a given time. In the case that the selected surface is wetted at a certain time due to spillover, the water temperature and water convection coefficient are considered; however, while there is no heat flux due to the insolation. Otherwise, if at a certain time the selected surface is dry (i.e., there is no spillover), the parameters of air convection (i.e., the air temperature and the air convection coefficient) and the measured heat flux due to insolation are taken into account.
During the analyzed year, 1306 h of water spilling over the spillway crest were recorded (from November to May), which is equivalent to 14.9% of the year. The longest period of continuous spillover was recorded between 14 March 2014 at 23:00 and 25

General
The temperature conditions of the dam were determined by means of the procedure developed for modeling the heat transfer process in concrete dams, taking into account time-varying boundary conditions on both the upstream side of the dam (the effect of the oscillation of the water level of the reservoir) and on the downstream side (the effects of insolation, shading, and water spilling over the spillway crest). Computer programs for the determination of temperature during time, at each node of the finite element mesh, were developed within the programming environments Mathematica [44] and MATLAB [45]. The thermal properties of the mass concrete used in the determining of the temperature conditions of the dam are given in Table 1.

Effect of Spillover
Spilling of water from the reservoir over the spillway crest has an effect on the determination of the parameters of convection (i.e., the temperature and the convection coefficient) and of the heat flux due to insolation, taking into account the effect of shading, on a selected surface, at a given time. In the case that the selected surface is wetted at a certain time due to spillover, the water temperature and water convection coefficient are considered; however, while there is no heat flux due to the insolation. Otherwise, if at a certain time the selected surface is dry (i.e., there is no spillover), the parameters of air convection (i.e., the air temperature and the air convection coefficient) and the measured heat flux due to insolation are taken into account.
During the analyzed year, 1306 h of water spilling over the spillway crest were recorded (from November to May), which is equivalent to 14

Effect of Shading
In order to investigate the effect of shading (i.e., solar exposure or non-exposure of the observation point on the downstream insolated side of the dam), a method developed by Žvanut et al. [41,46]-based on terrestrial laser scanner measurements of the topography of the wider area of the dam, and of the use of two computer programs (developed within the computing environments Mathematica and MATLAB, respectively)-was applied. Using this approach it is possible to determine-for any selected observation point-the elevation angles of the terrain at different azimuths (i.e., the contour of the terrain), the position of the Sun, and insolation over time. The proposed method has been confirmed as highly reliable in the case of the downstream side of the Moste Dam [46]. Figure 6 shows the contour of the terrain as seen from the observation point P1, which is located near the top of the downstream side of the dam ( Figure 3) at altitude 520.35 m, and the virtual path of the Sun in both extreme positions during the year (i.e., at the summer and winter solstices).

Effect of Shading
In order to investigate the effect of shading (i.e., solar exposure or non-exposure of the observation point on the downstream insolated side of the dam), a method developed by Žvanut et al. [41,46]-based on terrestrial laser scanner measurements of the topography of the wider area of the dam, and of the use of two computer programs (developed within the computing environments Mathematica and MATLAB, respectively)-was applied. Using this approach it is possible to determine-for any selected observation pointthe elevation angles of the terrain at different azimuths (i.e., the contour of the terrain), the position of the Sun, and insolation over time. The proposed method has been confirmed as highly reliable in the case of the downstream side of the Moste Dam [46]. Figure 6 shows the contour of the terrain as seen from the observation point P1, which is located near the top of the downstream side of the dam (

Solar Radiation
Back-analyses of the solar radiation were performed for the analyzed period of 15 consecutive clear days in the summer at the location of the MAWS, where the calibrated parameters were obtained [41]. Figure 7 shows the calculated absorbed heat flux due to solar radiation at observation point P1 on the downstream side of the dam in the case of year-round insolation and the measured absorbed heat flux obtained by the measurements of the insolation at the MAWS. The calculated values, which represent the theoretical daily maximums of the heat flux, consistently follow the measured values in the analyzed year, but due to different effects (i.e., cloudiness, partial transparency of the atmosphere) are higher.

General
The thermal analyses (1D and 2D) were carried out by using the FEM-based computer program TeEx (written in MATLAB, [47]) which was complemented by two specially formulated programs for determining the effects of convection and insolation, taking into account the effect of shading.
Since the same initial concrete temperatures were not the proper approximation of actual conditions at the beginning of analyses, the following procedure was implemented. The initial concrete temperatures were obtained in such a way that using preliminary initial temperatures of 20 °C, the temperatures after five days (for the 1D analyses) or after one year (for the 2D analyses) were first calculated, and then four repetitions of this period were performed, whereby the respective temperatures at the end of each period were used as the initial values of the next repetition. In order to eliminate the impact of the preliminarily selected initial temperatures, the initial values of the actual calculation were therefore defined as the temperatures calculated at the end of the final repetition.
The 15-days period (Section 3.2.1) and the one year period (Section 3.2.2) were first

General
The thermal analyses (1D and 2D) were carried out by using the FEM-based computer program TeEx (written in MATLAB, [47]) which was complemented by two specially formulated programs for determining the effects of convection and insolation, taking into account the effect of shading.
Since the same initial concrete temperatures were not the proper approximation of actual conditions at the beginning of analyses, the following procedure was implemented. The initial concrete temperatures were obtained in such a way that using preliminary initial temperatures of 20 • C, the temperatures after five days (for the 1D analyses) or after one year (for the 2D analyses) were first calculated, and then four repetitions of this period were performed, whereby the respective temperatures at the end of each period were used as the initial values of the next repetition. In order to eliminate the impact of the preliminarily selected initial temperatures, the initial values of the actual calculation were therefore defined as the temperatures calculated at the end of the final repetition.
The 15-days period (Section 3.2.1) and the one year period (Section 3.2.2) were first analyzed using a 1D model, represented by a line perpendicular to the concrete surface, for the instrumented locations of the Moste Dam (Figure 2). In both periods the changeable cloudiness with the effect of shading of the dam surface, and the effect of the varying water level of the reservoir with the occasional spilling of water over the spillway crest were taken into account [48]. The heat flux of the absorbed solar radiation was determined from the results of the MAWS measurements. Secondly, the temperature fields of the Moste Dam, during the year, were determined using a 2D model of the symmetrical crosssection of the complete dam, where all the previously mentioned effects were taken into account (Section 3.3). Thirdly, uncertainty analysis of the calculated temperature field of the Moste Dam on 1 August 2014 at 12:00 was performed, using normally distributed random variables (Section 3.4). Figures 8-11 show the measured and calculated concrete temperatures at three downstream gauges (TC1 and TC4 located in the upper borehole, and TC6 located in the lower borehole), with and without the effect of spillover, and at one upstream gauge (TC5), considering the water level of the reservoir, and also the measured air temperatures (or water temperatures, Figure 10) and the height of spillover during the analyzed 15-days period, from 10 May to 25 May 2014.

15-Days Period
It can be seen from Figures 8 and 9 that the calculated concrete temperatures at the downstream gauges TC1 and TC4, when taking into account the effect of spillover, are in very good agreement with the measured concrete temperatures during the whole period (the root mean square error (RMSE) is 1.95 at TC1 and 1.73 at TC4, whereas the mean absolute error (MAE) is 1.39 at TC1 and 1.67 at TC4).
The measured and calculated concrete temperatures at the upstream gauge TC5 (Figure 10), during the period of spillover, or more precisely, during the period when the water level of the reservoir is above the altitude of gauge TC5, are in very good agreement (RMSE is 0.45 and MAE is 0.31).
Appl. Sci. 2021, 11, x FOR PEER REVIEW 11 of 20 (TC5), considering the water level of the reservoir, and also the measured air temperatures (or water temperatures, Figure 10) and the height of spillover during the analyzed 15-days period, from 10 May to 25 May 2014. It can be seen from Figures 8 and 9 that the calculated concrete temperatures at the downstream gauges TC1 and TC4, when taking into account the effect of spillover, are in very good agreement with the measured concrete temperatures during the whole period (the root mean square error (RMSE) is 1.95 at TC1 and 1.73 at TC4, whereas the mean absolute error (MAE) is 1.39 at TC1 and 1.67 at TC4).
The measured and calculated concrete temperatures at the upstream gauge TC5 (Figure 10), during the period of spillover, or more precisely, during the period when the water level of the reservoir is above the altitude of gauge TC5, are in very good agreement (RMSE is 0.45 and MAE is 0.31).
From Figure 11 it can also be seen that there is very good agreement between the measured and calculated concrete temperatures at gauge TC6, taking into account the effect of spillover, assuming that the height of the spillover must be at least 0.1 m, otherwise the exposed surface of the lower downstream borehole is not wetted (RMSE is 1.28 and MAE is 1.07). Figure 3 shows the minor amount of spillover (when the height of spillover is less than 0.1 m), where it can be seen that the exposed surface of the upper borehole, with gauges TC1 to TC4, is wetted irrespective of the height of spillover, whereas the exposed surface of the lower borehole, with gauge TC6, is not wetted at a spillover height of less than 0.1 m, since it is located outside of the central cross-section of the dam and therefore not reached by the wetted area of the downstream side of the dam.
The maximum (or minimum) temperature at the specific point in concrete occurs later than the maximum (or minimum) air temperature. This time difference is termed the time lag and depends on the distance from the concrete surface. At a depth of 0.1 m (gauges TC1 and TC6) the time lag was fairly short (at maximum temperatures 1 to 2 h, and at minimum temperatures 3 to 4 h). At a depth of 0.2 m (gauge TC2) the time lag was slightly longer (at maximum temperatures 2 to 3 h, and at minimum temperatures 4 to 5 h). At a depth of 0.3 m (gauge TC5) the time lag at the maximum temperatures was 5 to 6 h, and 6 to 8 h at the minimum temperatures. At a depth of 0.5 m (gauge TC3) the time lag was already 10 to 12 h, and at a depth of 1.0 m (gauge TC4) the time lag was even several days.             From Figure 11 it can also be seen that there is very good agreement between the measured and calculated concrete temperatures at gauge TC6, taking into account the effect of spillover, assuming that the height of the spillover must be at least 0.1 m, otherwise the exposed surface of the lower downstream borehole is not wetted (RMSE is 1.28 and MAE is 1.07). Figure 3 shows the minor amount of spillover (when the height of spillover is less than 0.1 m), where it can be seen that the exposed surface of the upper borehole, with gauges TC1 to TC4, is wetted irrespective of the height of spillover, whereas the exposed surface of the lower borehole, with gauge TC6, is not wetted at a spillover height of less than 0.1 m, since it is located outside of the central cross-section of the dam and therefore not reached by the wetted area of the downstream side of the dam.
The maximum (or minimum) temperature at the specific point in concrete occurs later than the maximum (or minimum) air temperature. This time difference is termed the time lag and depends on the distance from the concrete surface. At a depth of 0.1 m (gauges TC1 and TC6) the time lag was fairly short (at maximum temperatures 1 to 2 h, and at minimum temperatures 3 to 4 h). At a depth of 0.2 m (gauge TC2) the time lag was slightly longer (at maximum temperatures 2 to 3 h, and at minimum temperatures 4 to 5 h). At a depth of 0.3 m (gauge TC5) the time lag at the maximum temperatures was 5 to 6 h, and 6 to 8 h at the minimum temperatures. At a depth of 0.5 m (gauge TC3) the time lag was already 10 to 12 h, and at a depth of 1.0 m (gauge TC4) the time lag was even several days. Figures 12 and 13 show the measured and calculated concrete temperatures at two of the downstream gauges (TC3 and TC4, located in the upper borehole), with and without the effect of water spilling over the spillway crest, and the height of spillover during the analyzed one year period, from 1 September 2013 to 1 September 2014. It can be seen from these figures that the calculated concrete temperatures at the downstream gauges TC3 and TC4, when considering the effect of spillover, are in very good agreement with the measured concrete temperatures during the whole year (RMSE is 2.50 at TC3 and 2.00 at TC4, whereas MAE is 2.20 at TC3 and 1.84 at TC4).

One Year Period
The annual concrete temperature oscillations for the analyzed year, at gauges TC1 to TC6, as well as the corresponding air (or water) temperature oscillations are shown in Table 2, Table 3, Table 4.      The annual concrete temperature oscillations for the analyzed year, at gauges TC1 to TC6, as well as the corresponding air (or water) temperature oscillations are shown in Tables 2-4.
The annual air temperature oscillations at the location of the MAWS were, for the selected period, between −9.5 and 35.6 °C, so their amplitude was 45.1 °C. At gauge TC1, the analyzed annual oscillations were significantly smaller (measured 35.3 °C and calculated 35.9 °C), whereas at gauge TC2 the oscillations were 30.7 and 29.4 °C, respectively.    The annual air temperature oscillations at the location of the MAWS were, for the selected period, between −9.5 and 35.6 • C, so their amplitude was 45.1 • C. At gauge TC1, the analyzed annual oscillations were significantly smaller (measured 35.3 • C and calculated 35.9 • C), whereas at gauge TC2 the oscillations were 30.7 and 29.4 • C, respectively. At gauge TC3 the measured annual oscillations were still 24.0 • C, and the calculated oscillations up to 20.6 • C. At gauge TC4 the measured oscillations, over the analyzed annual period, were 18.4 • C, and the corresponding calculated oscillations 16.9 • C.
The annual air temperature oscillations at gauge TW2 were between −1.8 • C and 32.1 • C, so that their amplitude was 33.9 • C, whereas the annual water temperature oscillations were from 2.7 to 19.0 • C and the amplitude was 16.3 • C. However, the annual water temperature oscillations at gauge TW1 were from 2.9 to 16.6 • C, and the amplitude was 13.7 • C. At gauge TC5 the measured annual oscillations were 13.8 • C, and the corresponding calculated oscillations 16.8 • C.
The annual air temperature oscillations close to gauge TC6 were, for the selected period, between −5.8 and 30.1 • C, so that their amplitude was 35.9 • C. The annual concrete temperature oscillations at gauge TC6 were significantly smaller (measured 24.7 • C and calculated 27.5 • C).

Temperature Field of the Dam
Based on the results of the 2D thermal analyses of the complete dam for the one year period, the temperature fields of the dam during the year were defined, using time periods of two months, starting from 1 October 2013 at 12:00. The temperature fields were plotted using the computer program DIANA [49], where an identical model of the dam as in MATLAB was constructed (i.e., it had the same arrangement and the same denotation of finite elements and nodes), so that when making the temperature fields in DIANA, the results of the thermal analyses calculations from MATLAB (the nodal temperatures at the selected times) were used. Comparing the calculated temperature fields, it was found that the largest gradients in the concrete temperature were recorded in the area close to the surface of the insolated upper part of the downstream side of the dam (Figure 14). This was in summer, when the highest air temperatures occurred and when not much shade was recorded on the downstream side of the dam. For determining the temperature gradient by depth, the two typical cross-sections, A-A and C-C, were chosen (Figures 2 and 14). The four selected points (nodes) coincided with the locations of the temperature gauges.  The concrete temperatures at the two characteristic cross sections in summer (1 August 2014 at 12:00) and in winter (1 February 2014 at 12:00) are presented in Figure 15. It can be seen from this figure that in cross-section A-A at the downstream insolated side of the dam in the summer the calculated temperature from the concrete surface to a depth of 0.2 m was reduced by 9.8 • C, whereas it was reduced at a depth of 1.0 m by 2.4 • C, and by an additional 2.0 • C at a depth of 2.0 m. This suggests that the temperature gradient is largest in a very narrow area along the surface of the concrete. The calculated concrete temperatures in cross-section A-A also indicate that the temperature stabilizes at a depth of about 6.0 m. In the analyzed case, in winter the temperature gradient of the calculated temperatures close to the surface was much smaller, but the strong effect of temperature can be observed up to a depth of 2.0 m. Both in summer and in winter the measured values at three points (at depths of 0.2, 0.5 and 1.0 m) were somewhat higher (by 1.3 to 2.5 • C) than the calculated values.

Uncertainty Analysis
As part of the thermal analyses, uncertainty analyses of the results of the calculations of the concrete temperatures were also performed, where six normally distributed random variables were used to determine the dispersion of the results. Their expected values and standard deviations are given in Table 1. Thermal analysis was carried out for a typical one-month period in the summer (from 1 July 2014 at 12:00 to 1 August 2014 at 12:00), where appropriate temperatures, obtained when calculating the temperature fields of the dam during a one-year period, were considered as the initial temperatures in the analysis (Section 3.3). The concrete temperatures (at 12 selected points on and inside the dam; Figures 2 and 14, Table 5) were analyzed at the end of a period of one month. Sixty repetitions of the calculation (for 60 different combinations of values of six normally distributed random variables), which resulted in 12 different samples, were performed. The calculated values of seven statistics of selected samples of the concrete temperatures are shown in Table 6, whereas the corresponding box-and-whiskers plots are presented in Figure 16. In cross-section C-C at the upstream side in the summer, when the water level of the reservoir was above the analyzed cross-section, the calculated temperatures from the concrete surface to a depth of 0.3 m were reduced by 2.4 • C, whereas up to a depth of 1.4 m they were reduced by an additional 1.7 • C. This suggests that the temperature gradient is the greatest up to a depth of 0.3 m, and is significantly lower in comparison with the analyzed case in winter (the temperature gradient up to a depth of 0.3 m was 5.8 • C) when the water level of the reservoir was below the analyzed cross-section, so that the concrete surface was exposed to the air. Both in summer and in winter the measured value at a depth of 0.3 m matched the calculated values very well (the differences were 0.1 • C and 0.6 • C, respectively).

Uncertainty Analysis
As part of the thermal analyses, uncertainty analyses of the results of the calculations of the concrete temperatures were also performed, where six normally distributed random variables were used to determine the dispersion of the results. Their expected values and standard deviations are given in Table 1. Thermal analysis was carried out for a typical one-month period in the summer (from 1 July 2014 at 12:00 to 1 August 2014 at 12:00), where appropriate temperatures, obtained when calculating the temperature fields of the dam during a one-year period, were considered as the initial temperatures in the analysis (Section 3.3). The concrete temperatures (at 12 selected points on and inside the dam; Figures 2 and 14, Table 5) were analyzed at the end of a period of one month. Sixty repetitions of the calculation (for 60 different combinations of values of six normally distributed random variables), which resulted in 12 different samples, were performed.
The calculated values of seven statistics of selected samples of the concrete temperatures are shown in Table 6, whereas the corresponding box-and-whiskers plots are presented in Figure 16. The maximum average temperature of the 12 samples of populations was obtained at the point A1 (30.24 °C), and the minimum at the point B5 (10.40 °C). The maximum standard deviations were also obtained at point A1 (1.15 °C), whereas the minimum standard deviations were identified at points B5, C1 and C2 (0.08 °C). It was found that the average temperatures of 12 samples of the populations, as compared with the calculated temperatures at the selected points, without considering the point of view of uncertainty, differed by 0.1 °C at four points (A1, A3, A4, B1), whereas at the other eight points there was no difference. The maximum temperature of 12 samples of the populations was obtained at the point A1 (32.6 °C) and the minimum at the point B5 (10.3 °C). At point A1 the maximum interquartile range (1.83 °C) and the maximum range (4.5 °C) were also registered. The minimum interquartile range was calculated for the point C1 (0.03 °C), whereas the minimum range was recorded at the points B5, C1 and C2 (0.3 °C).   The maximum average temperature of the 12 samples of populations was obtained at the point A1 (30.24 • C), and the minimum at the point B5 (10.40 • C). The maximum standard deviations were also obtained at point A1 (1.15 • C), whereas the minimum standard deviations were identified at points B5, C1 and C2 (0.08 • C). It was found that the average temperatures of 12 samples of the populations, as compared with the calculated temperatures at the selected points, without considering the point of view of uncertainty, differed by 0.1 • C at four points (A1, A3, A4, B1), whereas at the other eight points there was no difference. The maximum temperature of 12 samples of the populations was obtained at the point A1 (32.6 • C) and the minimum at the point B5 (10.3 • C). At point A1 the maximum interquartile range (1.83 • C) and the maximum range (4.5 • C) were also registered. The minimum interquartile range was calculated for the point C1 (0.03 • C), whereas the minimum range was recorded at the points B5, C1 and C2 (0.3 • C).

Conclusions
This article presents a relatively simple procedure for modeling the heat transfer process in a concrete dam, taking into account time-varying boundary conditions on the upstream and downstream sides of the dam (i.e., the water level of the reservoir, spillover, insolation, and shading) which affect the dam's thermal fields.
The validity of the procedure was verified on the example of the large arch-gravity concrete Moste Dam (the highest such dam in Slovenia), where a sophisticated automated system has recently been installed. This system is used for the measurement of concrete temperatures, as well as to monitor the surrounding conditions (i.e., the water temperatures, the water level of the reservoir, the height of spillover, the air temperature and the solar insolation).
Thermal analyses (1D and 2D) for the analysis of non-linear and non-stationary heat conduction through solids were performed using a FEM-based program, which was complemented by two specially developed programs for determining the effects of boundary conditions during the analyzed 15-day period, as well as during a period of one year. Uncertainty analyses were also carried out to estimate inherent random effect.
It was found that the results of the performed thermal analyses fitted in well with the experimentally determined concrete temperature measurements. The results showed that at the insolated side of the dam, the temperature gradient was largest in a very narrow area along the concrete surface, but the temperature did not stabilize at depths shallower than about 6 m.
The obtained results (i.e., seasonal temperature distributions in the dam) affect the thermal loads, which can cause high tensile stresses in concrete and importantly affect the occurrence and behavior of cracks, as well as displacements; especially in arch dams and arch-gravity dams.

Funding:
The presented article is part of research work that was performed within the scope of the first author's doctoral studies, and was partly funded by the European Union through the European Social Found.

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Not applicable.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author. The data are not publicly available due to the requirements of the company Savske elektrarne Ljubljana Ltd.