Modelling Salt Intrusion and Estuarine Plumes under Climate Change Scenarios in Two Transitional Ecosystems from the NW Atlantic Coast

Abstract

The Minho and Lima estuaries, separated by approximately 20 km, are two important transitional systems from the NW Atlantic Portuguese coast. Considering their contrasting characteristics, the present work aims to investigate climate change impacts in salt intrusion and estuarine plumes interaction, through an integrated hydrological characterization. The three-dimensional numerical model Delft3D was implemented to explore the hydrological response of these two systems taking into account their interaction. The model was applied to investigate two different climate change scenarios (RCP 4.5 and RCP 8.5), where physical variables such as sea level and river outflow were varied following summer and winter conditions. The numerical results revealed typically estuarine hydrological patterns, although some different features were found for each estuary. The head and mouth of the estuaries may be the zones most affected by climate change. In the first zone, the salinity may increase in the future in response to freshwater reduction from the Minho and Lima rivers. Considering the cumulative impact of the sea level rise, the saline intrusion may increase in the future: in the Minho estuary it reaches 14 km for the Present scenario during the dry season, penetrating further under climate change conditions (15 km for the RCP 4.5 scenario, and 17 km for the RCP 8.5 scenario); in the Lima estuary, the salt propagates further inland, reaching 17 km for the Present scenario, 18 km for the RCP 4.5 scenario, and a maximum of 20 km for the RCP 8.5 scenario. In all scenarios, during winter conditions, the estuarine plumes are larger and wider, being predicted that in the future they will reduce their dimensions during the winter period, and remain unchanged during summertime. Concisely, this modelling study provided more expertise about the physical patterns of these two estuaries and their interaction, as well as to anticipate the climate change impacts in their dynamics.

1. Introduction

Identified changes (applying statistical tests) in the climate state mean and/or in its property variability, persisting for long periods (from decades to centuries), is designated as climate change [1,2,3,4,5,6]. The internal dynamics of the climate system associated with external forcing agents are the factors described as responsible for the evolution of the climate system, with impacts on land, water, and atmosphere [3,7], which, in turn, responds through various feedback mechanisms [2]. In the marine environment, the sea surface warming, the sea-ice melting, the changes in circulation and mixing patterns, the ocean acidification, the sea level rise, the salt intrusion, the oxygen depletion, and the frequency and intensity of storms induced by the actual climatic modifications are changing the marine ecosystem’s structure and function [8,9,10].

Estuaries and lagoons are complex transitional environments in coastal areas where interconnected land and ocean hydrodynamic processes occur, exchanging matter, energy, and information between terrestrial and marine ecosystems [11,12,13,14]. These processes are regulated by intense physical, biological, and chemical interactions forced by waves, tides, river discharge and sediment supply, wind stress, and exchanges with the atmosphere [15,16].

In the period from 1901 to 2010, the global mean sea level rose by 0.19 ± 0.02 m. Compared with the period 1986–2005, future projections predict an increase between 0.26 and 0.82 m by the end of the century (2081–2100), so an average increase of 0.50 m is considered more likely pp. 1137–1216 in [2], pp. 321–445 in [5] and [17]. Today, the sea level rise results from an increase in the volume of the water already present in the ocean, since water warms and expands, as well as from an increase in the mass of the water in the ocean due to the melting glaciers [2]. This phenome exacerbates extreme events and coastal hazards, and has numerous adverse impacts on marine coastal ecosystems pp. 321–445 in [5], leading to environmental and human life quality risks associated with floods, shoreline erosion, and other coastal problems [18]. This situation induced by global warming has increased the frequency and intensity of coastal inundation, with significant impacts on water level, current, and wave propagation in estuaries [19], leading to, besides others, alterations in saltwater intrusion and estuarine plumes spatial and temporal gradients variability [20,21].

Seawater intrusion resulting from mean sea level rise is also a key concern for coastal and estuarine environments around the world [21,22,23,24,25,26], influencing water residence times, and inducing important biological modifications in species abundance, composition and distribution, and in the provision of water for coastal communities activities, such as industrial and irrigation uses and soils fertility and cultivation [18,26,27,28,29,30,31].

The estuarine plume is the riverine outflow that can be directly observed near the estuary’s mouth by the different colours between estuarine and adjacent ocean waters [32]. These waters also have lower density and carry with them large quantities of dissolved and suspended river materials [33,34]. Through spreading and advection, these plumes also contribute to the horizontal redistribution over long distances of nutrients and/or pollutants [35]. The dispersion and impact of estuarine material along the coast and over-shelf depends on the interaction of outflow inertia, buoyancy, Earth’s rotation, wind forcing, and alongshore coastal current [15,33,34,36].

The two estuaries investigated in this study are located on the Portuguese northwest coast, adjacent to each other, in the hydrographic region of the Minho and Lima rivers. Their ecological importance is very high, as they provide significant protection statutes related to their high primary production, fish’s shelters and nurseries, recognized by Natura 2000 Network [37,38], Biotope CORINE [39,40], bird’s Special Protection Area [41], and Important Bird Area (IBA) [42,43]. The development of human populations and consequent anthropogenic stress near the estuarine borders turn these regions into vulnerable systems requiring frequent surveillance [44,45]. These water bodies have similar latitudes (separated by only 20 km), being subject to the same oceanic, atmospheric, and environmental conditions, however have different characteristics such as size, shape, depth, and fluvial inputs [32,46].

Due to the importance of these two estuaries, the application of models capable to analyzes climatic impacts is important to understand the changes of estuarine hydrodynamics under future scenarios. Numerical models applications able to reproduce physical, chemical, and biological interactions in estuarine systems emerge to overcome the absence of long-term continuous observations and to deliver forecasts of the natural processes in the coastal systems [47]. This kind of information can be very useful for a handful of applications, such as scientific studies on marine ecosystems, monitoring the quality of the seawater and freshwater, and support of marine and coastal management decisions [48,49].

In this context, the Delft3D model arises to study hydrodynamics, sediment transport, morphology, and water quality for fluvial, estuarine, and coastal environments. The application of this model to investigate, for instance, the hydrodynamic of the Arabian Gulf [50] and the hydrodynamic and sediment transport in the Yangtze estuary (China) [51], the re-distribution of nourished sediments in Terschelling (Netherlands) influenced by tides, wind and waves [52], the concentration of suspended sediments in Germany [53], and the propagation patterns of the Tagus and Sado estuarine plumes in Portugal [54], has proven its good capacities to study different processes in a large number of coastal systems.

Several hydrodynamic and morphodynamic studies, analyzing the dependence of the hydrography on tidal variability, the concentrations of the heavy metal and metalloids in sediments, and the development of some spatial planning instrument to climate change adaptation in Minho and Lima estuaries may be found [38,55,56,57,58,59,60,61,62]. However, none of them discusses the 3D distribution of water temperature and salinity within the estuaries, nor the climatic impact on salt intrusion and plumes dynamics.

Therefore, the main goal of this paper is to use a single three-dimensional baroclinic model application to study the salt intrusion and estuarine plumes dynamics in two adjacent estuaries, researching the modifications induced by climate change. This assessment will provide a baseline against which the extent of current activities can be measured, and a possible solution can be planned under the EU directives.

2. Materials and Methods

2.1. Study Area

The Minho and Lima estuaries are located on the northwest coast of Portugal (Figure 1) and are under the influence of a maritime climate (mesothermal and highly wet weather) [57,63].

The Minho river is an international freshwater body with a special aspect, which forms a natural border between the north of Portugal and the Galicia region (Spain). The river has 300 km of total length and a hydrographic basin area of 17080 km², the upper part (≈230 km) and 95% of the watershed stay in Spain, and the lower part (≈70 km) and 5% of the remaining basin in Portugal [64]. In the lower part of the river is placed the Minho estuary with 35 km long and with a total area of 23 km² (Figure 1). The estuary depth has its maximum near the mouth (approx. 6 m), with an average depth of 2.6 m [65]. The astronomical tide range, with a semi-diurnal and mesotidal regime, varies between 2 m and 4 m (neap and spring tides, respectively) [66]. The Minho estuary is considered a partially mixed system, becoming a salt wedge estuary during high flood periods [64]. The annual average river discharge is 305 m³s⁻¹, varying seasonally between 121 m³s⁻¹ in the summer season and 692 m³s⁻¹ in the winter season [46].

Going to the south, the Lima river is located, which rises in Orense (Spain) and flows to Viana do Castelo (Portugal). This Iberian river has a total length of 108 km (41 km in Spain and 67 km in Portugal), and a basin area of 2460 km² [63]. The Lima estuary is 22 km long and has a total area of 5.5 km² (Figure 1). The mean depth of the estuary is 3 m, although due to the constant dredging activities, the maximum depth is 10 m in the main channel [63]. The estuary’s mouth is defined by a 2 km long jetty, that directs the river outflow to the south. The tide is characterized by a semi-diurnal and mesotidal regime like Minho estuary, with low and high tidal ranges of 1 m and 3.4 m, respectively [62]. Lima is considered a partially mixed estuary, and the annual average river discharge (≈70 m³s⁻¹) is very variable, as it depends on the human control of two dams (Touvedo and Alto Lindoso) located upstream [67,68].

The geographical position and the consequent proximity to the Atlantic Ocean, the influence of the Azores Anticyclone, as well as the shape and arrangement of the main mountain groups that surround it, make this area one of the rainiest in Portugal [69]. In the Minho river basin, the average annual rainfall is around 1600 mm, while near the coast the annual averaged evaporation is 700 mm and inland is 860 mm. In the Lima river basin, the average annual precipitation is higher and varies between 1300 and 4200 mm, being over 3000 mm upstream and under 1900 mm downstream. The regime of the annual precipitation in this region is divided into three periods: a wet period which generally occurs between October and March (autumn and winter), a transition period throughout the spring months (April to June), and finally a drier period from July to September (summer months) [57,63,69].

2.2. Model Set-Up

In this modelling research, a single domain spherical coordinated grid was used that was previously created by [70], representing simultaneously the Minho and Lima estuaries (black rectangle shown in Figure 1), and the adjacent coastal region.

The bathymetry of the coastal region was constructed with General Bathymetric Chart of the Oceans (GEBCO) data [71]. The Minho estuary bathymetry was interpolated from Light Detection And Ranging (LIDAR) data with a spatial resolution of 2 m, collected in 2011 by the Portuguese Administration of Territory [72]. The Lima estuary numerical bathymetry was obtained from the interpolation of data provided by the Hidrodata company, with 5 m of spatial resolution. This data was measured in three different locations: the main port area (Hydrographic Survey 0255-ARN-1002), the port entrance (Hydrographic Survey 0286-ARN-0503A), and the area between the Eiffel bridge and A28 highway (Hydrographic Survey 0330-ARN-1203A), between October 2002 and December 2003 [73].

At the ocean open boundary, the astronomic forcing was applied with a spatial resolution of 0.25°, represented by the major 13 tide constituents obtained from the TOPEX/Poseidon solution [74]. The transport conditions were calculated from the reanalysis product IBI_MULTIYEAR_PHY_005_002 from the Copernicus Marine Environment Monitoring Service (CMEMS), using daily salinity and water temperature data from the surface down to 5500 m in depth with 9 km of horizontal resolution [75].

From the ERA-Interim atmospheric reanalysis model, with 3 h of temporal resolution, the data for the components of the heat flux model (relative humidity (%), air temperature (°C) and solar radiation (Jm⁻²s⁻¹)) and wind process (wind velocity (ms⁻¹) and direction (°)) was provided by the European Centre for Medium-Range Weather Forecasts (ECMWF) [76].

Two discharge points representing the fluvial open boundary were defined as outflows of the Minho and Lima rivers. For the Minho river, daily series of river discharge data at Frieira station was obtained from the Confederación Hidrográfica del Miño-Sil in Spain [77], while for the Lima river it was obtained the discharge data from the Portuguese National Information System of Water Resources at Ponte da Barca station [78].

The simulation period elapsed through fourteen months (between the years 2012 and 2013) with a time step of 30 s and 10 sigma layers (privileging the surface dynamics), and a bottom roughness with manning coefficients of 0.016 and 0.022 in the coastal area and inside the estuaries, respectively. The background horizontal eddy viscosity and diffusivity were set to 10 m²s⁻¹, and the background vertical eddy viscosity and diffusivity to 0 m²s⁻¹. The k-Epsilon was used as a turbulence model.

This model’s configuration was previously implemented and calibrated by [70]. The calibration of the hydrodynamic model was performed comparing sea surface elevation determined by the model with time-series data in six locations (five in the Minho river and one in the Lima river). The author found RMSE values between 8 and 10 cm, with relative errors of 4 and 6.5%, respectively. Lower RMSE values were found at stations closer to the estuary mouth, and higher values at stations further upstream. However, in all stations, a skill value of 0.99 was obtained. Regarding the transport module calibration, it was achieved comparing in situ water temperature and salinity data and model outputs in a total of 14 sampling sites. The analyses made by the author concluded that the model reproduces satisfactorily the patterns and seasonal evolution of water temperature and salinity in both estuaries. In the Minho estuary, the values of RMSE and skill are more homogeneous between the mouth and the head of the estuary, whereas in the Lima estuary the model results improve from downstream to upstream.

2.3. Climate Change Scenarios

Based on the IPCC synthesis report of 2014 [17], the future scenarios defined for this study were designated as RCP 4.5 and RCP 8.5, representing the respective Representative Concentration Pathways (RCPs) mentioned in the report. The RCPs describe four different 21st-century trajectories of land use, greenhouse gases and air pollutant emissions, and atmospheric concentrations.

To represent the actual situation, the calibration period (years 2012–2013) was defined as the Present scenario. The selection of this period was sustained on the higher confidence in the model results for these years as they were validated in the model’s calibration, and on the results of a statistical analysis of the main atmospheric, oceanic and river variables used as model inputs. For this analysis, a data collection was performed of air humidity and temperature, solar radiation, wind velocity, ocean salinity and temperature, and Minho and Lima outflows for the period 2010–2020 from the same sources stated in Section 2.2. These data were statistically analyzed and calculated the daily maximum, minimum and average values for each of these variables over 365 days, and the results were then compared with the calibration period (1 March 2012–1 March 2013).

For the future scenarios, changes were applied regarding the Present scenario in the model’s ocean, atmospheric, and river boundary conditions, namely, the mean sea level, the water temperature, the ocean salinity, the river discharge, the heat flux, and the wind process. In the model’s heat flux, climatic changes were made in the air temperature, air humidity, and solar radiation conditions, while in the model’s wind process, the wind velocity was the only adjustment.

By the end of the century, the IPCC report [2] exhibits a global sea-level rise between 0.32 and 0.63 m for scenario RCP 4.5, and between 0.45 and 0.82 m for scenario RCP 8.5, relative to the 1986–2005 period. In this work, a rise in sea level of 0.32 m was considered for RCP 4.5 and 0.82 m for RCP 8.5, representing the most optimistic and pessimistic cases, respectively. Regarding future oceanic water temperature and salinity, the implemented values were obtained from Phase 5 of the Coupled Model Intercomparison Project (CMIP) established by the World Climate Research Programme (WCRP) [79] averaging the predictions of four global models (CNRM-CM5, HadGEM2-ES, IPSL-CM5A-MR, and MPI-ESM-LR).

In riverine discharges, the adjustments for the future scenarios were acquired from the Swedish Meteorological and Hydrological Institute (SMHI) [80]. For the RCP 4.5 scenario, results from five global and five regional models with one hydrological model (E-HYPE2.5 model) were used to calculate the percentage difference to the reference period 1970–1999. For the RCP 8.5 scenario were considered the results from four global and regional models, and the same hydrological model for the same reference period.

Table 1. Changes at the ocean and river boundary inputs applied in the simulations of future scenarios.

	RCP 4.5	RCP 8.5
Sea Level (m)	+0.32	+0.82
Water Temperature (°C)	+1.55	+2.20
Ocean Salinity	−0.65	−1.55
River Discharge (%)	−5	−25

presents all the changes applied at the ocean and river model’s boundary.

Lastly, the changes applied in the model’s heat flux and wind system were based on Portal do Clima website data [81]. The variables withdrawn from the website for the model’s heat flux were air temperature, relative humidity, and solar radiation. For the wind system, only the wind speed was used. The extracted data consisted of monthly average anomalies regarding the reference period 1971–2000 and was calculated from five regional models (CLMcom-CCLM4-8-17, DHI-HIRHAM5, IPSL-INERIS-WRF331F, KNMI-RACMO22E, and SMHIRCA4) and four global models (CNRM-CERFACS-CNRM-CM5, ICHEC-EC-EARTH, IPSL-IPSL-CM5A-MR, and MPI-M-MPI-ESM-LR).

Table 2. Changes at the atmospheric boundary inputs applied in the simulations of future scenarios.

	Summer		Winter
	RCP 4.5	RCP 8.5	RCP 4.5	RCP 8.5
Air Temperature (°C)	+2.06	+4.73	+1.81	+3.75
Air Humidity (%)	−1.63	−5.93	+2.03	−0.42
Solar Radiation (Jm−2s−1)	+7.74	+13.42	−6.73	−3.80
Wind Velocity (ms−1)	−0.01	−0.02	+0.45	+0.11

presents the changes at the model’s atmospheric boundary.

All scenarios were divided into two periods (summer and winter), resulting in a total of six simulations (two for Present, two for RCP 4.5, and two for RCP 8.5). The future simulations for the summer season ran from 1 May to 1 September 2100, and for the winter season between 1 November 2100 and 1 March 2101. The first month of the simulations was used to spin up the model. Two additional scenarios derived from the RCP 8.5 scenario were simulated to assess the relative importance of changes in oceanic and river boundary conditions in defining future estuarine conditions. The RCP 8.5 SL scenario was simulated considering only the sea level rise, and the RCP 8.5 RD scenario only considering the changes in the river’s discharge.

3. Results

3.1. Present Scenario Analysis

As previously mentioned, the period from 1 March 2012–1 March 2013 was chosen to represent the Present scenario considering the higher confidence in the model results for this period, as well as the results of a statistical analysis performed on the main input variables of the model. The results of this analysis are exposed in Figure 2, showing that the climatological patterns for the decade 2010–2020 are very well represented by the calibration period. In fact, at no time were the maximum and minimum values of each variable exceeded during the calibration period. Except for rivers discharge, all the atmospheric and oceanic variables for the years 2012–2013 (represented by the grey dashed trend line) are very close to the climatological averages (represented by the green dashed trend line). Essentially, only the river flow values for 2012–2013 are below the climatological values, as the freshwater inflow into each estuary is dependent on human control carried out in the dams existing in these rivers.

3.2. Salt Intrusion

The saltwater movement upward a coastal system is designated as saline or salt intrusion. The higher salt content present in ocean water makes it denser than the brackish water that is found inside the estuary or lagoon. Consequently, when these two water bodies meet, the seawater sinks, leading to a salt intrusion near the bottom and in a wedge shape [22].

The Venice system is commonly used to classify a coastal system through its salinity distribution [82]. This system divides an estuary into six sections based on the longitudinal distribution of salinity. Near the estuary mouth exists the euhaline zone (with salinities over 30), followed by the polyhaline zone with a lower section (salinities from 25 to 30) and a middle section (salinities from 18 to 25). Further upstream stays the mesohaline zone with salinities ranging from 5 to 18, and the oligohaline zone with salinities between 0.5 and 5. Finally, there is the river section with salinities below 0.5 denominated as the limnetic zone.

In this study, the investigation of the saline maximum intrusion was performed averaging all the model results for the deepest layer’s under high tide conditions, for the summer (June, July, and August) and winter (December, January, and February) months, and considering that the saline intrusion extends until the 0.5 units of salinity. Thus, a map of the average salt intrusion at high tide during winter conditions is represented in Figure 3 for the five scenarios mentioned in Section 2.3.

From Figure 3 it was observed that both estuaries present similar longitudinal salinity distributions for all the scenarios, but with higher values in the Lima estuary. The mouths of both estuaries are characterized as euhaline, while the central zone of the Minho estuary is mesohaline and polyhaline and that of the Lima estuary is dominated by higher salinities (euhaline and polyhaline sections).

The upper estuaries present major variations in future projections. In the Minho case, there is a clear increase in the length of the oligohaline section from the Present scenario to the RCP 4.5 scenario, and from this last to the RCP 8.5. Consequently, the limnetic area with lower saline water may decrease for future scenarios. Considering Figure 3d,e, the oligohaline zone is visually longer in the first case, showing the major impact of the sea level rise regarding the river discharge. In the Lima estuary, compared to the Present scenario, is also found a propagation of the saltier sections upstream for RCP 4.5 and RCP 8.5 scenarios, however, the impact of the sea level rise in this estuary is not so significant as for the Minho estuary (Figure 3i,j).

Figure 4 represents the longitudinal distribution of salinity along the vertical sections represented in Figure 1, for the Present, the RCP 4.5, and the RCP 8.5 scenarios during the wet season. Both estuaries present vertical stratification in the first 5 km of the longitudinal section, although there is more evidence for the Lima estuary. In Minho estuary, the length of the salt intrusion is around 11 km from the estuary’s mouth for the Present scenario, about 12 km for the RCP 4.5 scenario, and near 13 km for the RCP 8.5 scenario. For the Lima estuary the length of the salinity intrusions is slightly higher (approximately 13 km for the Present and RCP 4.5 scenarios, and 14 km for the RCP 8.5).

From the analysis of subfigures of Figure 3 and Figure 4 corresponding to RCP 8.5 SL and RCP 8.5 RD scenarios, it is found that sea-level rise has a major impact on the higher saline intrusion observed in climate change scenarios.

During summer conditions, Figure 5 and Figure 6, both estuaries present typical estuarine longitudinal distributions, but with higher salinity values than during the wet season. The central region of the Minho estuary is dominated by euhaline and polyhaline zones, while the Lima estuary presents euhaline and lower polyhaline sections. These results from a reduction of the fluvial discharges during the summer months. In the mouth region, salinities varying between 30 and 36 are observed in both systems, while in the heads of the estuaries prevail salinities between 0.5 and 5, and between 18 and 25, in the Minho and Lima estuaries, respectively. During this dry season, the two systems are vertically homogeneous, and the salinity intrusion increases progressively within the different scenarios. In the Minho estuary, a salinity intrusion of nearly 14 km is observed in the Present scenario, 15 km in the RCP 4.5 scenario, and 17 km in the most critical scenario (RCP 8.5). As for the winter situation, the range of saltwater intrusion is larger for the Lima estuary, reaching 17 km, 18 km, and 20 km in the Preset, RCP 4.5, and RCP 8.5 scenarios, respectively.

The comparative analysis between the scenarios designed to assess the sea-level rise and reduced river flow impacts show results similar to that found during the winter months, evidencing the higher effect of the first.

3.3. Estuarine Plumes

Coastal low-density water masses carrying dissolved and suspended river materials with different properties from the surrounding ocean waters are designated as estuarine or riverine plumes [32]. Being that buoyancy is an inherent characteristic of these plumes, an average of all the top layer’s model results was performed during the ebb time for the summer (June, July, and August) and winter (December, January, and February) months, identifying the plume in the superficial layer.

In Figure 7 is presented the horizontal averaged salinity distribution after low tide during the wet season for the Present scenario, as well as the respective differences with the future scenarios. For Minho and Lima estuaries is observed the generation of a confined bulge in front of their mouths, however, in Minho this presents a larger dimension. This different pattern is essentially induced by the higher freshwater input compared to the Lima River, considering the wind conditions used in the simulations.

From the analysis of differences between the RCP 4.5/RCP 8.5 scenarios and the Present scenario, positive values are observed, resulting from the lower input of freshwater discharge in the future and therefore from higher salinities found within both estuaries. This pattern is similar for both climatic scenarios, however, it is more pronounced for the RCP 8.5 scenario.

A connection near the coast between Minho and Lima plumes is not observed in average maps, nevertheless, under high river flow and appropriate wind direction and velocity, the Lima plume can attach with the Minho plume, integrating the Western Iberian Buoyant Plume [70].

The analysis of the vertical salinity sections along the plumes/estuaries axis (Figure 8) shows that during winter months the Minho estuary has a slight stratification at its mouth, whereas in the Lima estuary the stratification extends along the first 6 km. The plume’s generation in front of the Minho estuary is more pronounced compared with the Lima case and is observed in the first 3 or 4 m of the water column. An increase in salinity at the mouth of the estuary is also observed for the future scenarios but is more relevant to the RCP 8.5 scenario. The vertical variations are negligible for the Lima estuary.

The surface salinity distribution under summer conditions is shown in Figure 9. Here the plume’s generation is less evident compared with the wet season. This difference is justified by the decrease in the average river discharge, and it is also corroborated by the higher differences in the central regions of the estuaries. Contrary to the winter situation, during the summer the plume size remains unchanged for all simulated scenarios.

The generation of the plume is not observed along the water column over the dry season (Figure 10). In the Minho estuary, the salinity vertical homogeneity is observed at the estuary mouth, with a slight increase between the future scenarios and the Present scenario. Neither vertical differences nor higher salinities for the future are noticeable for the Lima estuary.

4. Discussion

Channel topography, river outflow, and tidal characteristics are mainly responsible for the mixing processes in estuaries affecting directly the saline distribution [83]. Sometimes the influence of rainfall and evaporation is also important.

In this frame, human activities may strongly impact salt distribution in estuaries, especially considering that it is being pointed out as one of the main causes of climate change, particularly of global warming and consequent sea-level rise. Currently, the climate change impacts on coastal areas have been increasingly investigated due to the high socio-economic importance of these regions. Minho and Lima’s estuaries are coastal systems with different characteristics (such as size, shape and depth) that have been studied individually in past investigations. This study intends to apply a model configuration comprising both estuaries, capable of evaluating the estuarine hydrological response to climate change.

To achieve this, changes in the model’s oceanic, atmospheric and river boundaries were applied to create five different scenarios (Present, RCP 4.5, RCP 8.5, RCP 8.5 SL and RCP 8.5 RD). For each scenario, the estuarine response was assessed under summer and winter conditions, where the differences between both seasons were defined through different river discharge, air and water temperature, and wind velocity. Model results were analyzed considering longitudinal and vertical saline distribution, with special emphasis on the generation of the salt wedge and estuarine plumes.

Minho and Lima’s estuaries present a typically estuarine longitudinal salinity structure, with higher salinities near the mouth and lower saline water upstream. The vertical structure, during summer and winter conditions, shows a stratified Lima estuary and a homogeneous Minho estuary.

The salt intrusion was evaluated through the Venice system. During summer and winter seasons both estuaries, in central regions, are dominated by polyhaline and euhaline zones, however, the Lima estuary presents more saline waters. The upper areas have a different classification between summer and winter periods. Due to the higher river discharge during winter months, the upstream zones are identified as limnetic in Minho estuary and oligohaline in Lima estuary. In summer these estuaries become saltier, and therefore their classification changes: oligohaline and middle polyhaline in Minho and Lima, respectively. The highest salinity values found in the Lima estuary are caused by the fluvial effect. The Lima river has an annual average discharge about four times lower than the Minho river [46,67], turning the Lima estuary saltier. The salt wedge increment found in the RCP 8.5 scenario is a natural response of the system to the future sea-level rise and river flow reduction (Figure 3, Figure 4, Figure 5 and Figure 6). The rise in sea level pushes a higher volume of saltwater into the estuary prolonging the salt wedge, and these results are concordant with other Portuguese estuarine modelling studies [84,85,86], and other global locations like the United States and China [87,88,89].

Finally, it was analyzed the generation and propagation of the estuarine plumes. Model forecasts, during the dry season, shows small plumes (Lima’s plume being smaller than Minho’s plume) in front of the estuaries entrance. During this period, low-intensity winds (≈4 ms⁻¹) blew predominantly from the East, confining the plume near to the estuaries mouth. This result concurs with [46], which during identical wind conditions (East winds with small velocity 3 ms⁻¹) revealed a Minho plume with a similar shape. During winter conditions with the intensification of the rivers discharge, for the Present scenario, model predictions showed the generation of a robust plume in the Minho estuary mouth, and a smaller plume in Lima. The plume expands in both estuaries for the RCP 4.5 and RCP 8.5 scenarios, with a slight amplification in the last case. These plume patterns are induced by the quantity of freshwater runoff and wind conditions, being consistent with similar studies in other coastal systems [90,91,92].

5. Conclusions

The salt intrusion and estuarine plumes generation in Minho and Lima estuaries, and their modifications induced by climate change, have been evaluated using the Delft3D model, considering different scenarios simulated using a unique modelling application.

Due to climate change, future sea-level rise and river discharge reduction are the principal factors inducing the hydrological modifications in both estuaries, but the first was identified as the main feature prompting the predicted changes. For the future, it is expected the upstream prolongation of the salt wedge in both estuaries, especially under the most pessimistic scenario RCP 8.5, with a larger elongation in the Lima estuary. The estuarine plumes generated at the estuaries entrance are larger and more pronounced during the winter months due to the larger volume of freshwater discharged by the rivers. Vertically, the Minho estuary is classified as homogeneous in the summer and winter months, and the Lima estuary as homogeneous in the dry season and partially stratified in the wet season. For the future, a slight decrease in the size of the plumes is expected for both estuaries in the winter case, particularly for the RCP 8.5 scenario, while for the summer no changes in the size of the plumes are expected.

The present model configuration identified the major factors inducing the hydrological changes in both estuaries in contrasting scenarios, thus, it can be used as a solid tool to help in management policy decisions, especially related to both ecosystems health protection, and mitigation and climate resilience plans, crucial to coastal areas, and in line with the conclusions from the last IPCC report recommendations.

The use of numerical models in transitional ecosystems is therefore of the utmost importance, especially in changing scenarios resulting from increasing climatic changes. Their application to enhancement adaptive strategies to minimize the losses of the ecosystem services have been recognized as decisive by the environmental policy and the scientific community. In fact, recent modelling studies have highlighted the serious susceptibility of the marine coastal environment in South Europe to climatic modifications [93,94,95], being consistent with the results herein presented, reinforcing the importance of this study.