Effects of sea level rise on salinity and tidal flooding patterns in the Guadiana Estuary

Global mean sea level is rising at a rate of approximately 3.2 millimeters per year (Church et al. 2013). According to the most recent report from the Intergovernmental Panel on Climate Change (IPCC), the rate of global mean sea level rise has been increasing over the last two centuries and continues to accelerate. Sea level rise will affect coastal areas, which should be a concern considering 10% of the world's population live within 10 m elevation of the present sea level (Carrasco et al. 2016). Direct consequences of sea level rise on coastal areas include an increase in flooding area, an increase in erosion, and an increase in salinity and changes in ecosystems (Nicholls et al. 2011). Of relevance to the present study is the impact of sea level rise on estuaries, where rivers intersect sea and freshwater mixes with saltwater.

Salt intrusion is a direct consequence of sea level rise in estuaries (McLean et al. 2001). Estuarine circulation is mainly driven by freshwater flow, tides, and density differences (Garel et al. 2009). A study by Chua & Xu (2014) found a stronger longitudinal salinity gradient in estuaries due to sea level rise, which in turn drives a stronger gravitational circulation. The increase in salinity will cause the water to become denser, and thus increase the stratification of the water column. Changes in estuarine stratification and circulation will further cause oxygen depletion (Hong & Shen 2012). These effects are detrimental to ecosystem services and marine habitat as estuaries hold highly diverse and complex ecosystems (Sampath et al. 2015). The objective of this study is to evaluate the evolution of hydrodynamics and changes in salinity distribution in response to various sea level rise scenarios in a major estuary in the Iberian Peninsula, the Guadiana Estuary.

The Guadiana Estuary is formed at the interface of the Guadiana River and the Gulf of Cadiz. The head of the Guadiana River begins in Spain and extends 810 km south toward the Gulf of Cadiz (Delgado et al. 2012). From its mouth in front of Vila Real de Santo António, the Guadiana Estuary extends 80 km north to its tidal limit near Mértola, Portugal. Because of its narrow width and moderate depth, the Guadiana Estuary is considered a rock-bound estuary where the volume of water entering the estuary during the flood tide is larger than the freshwater discharge (Garel & D'Alimonte 2017). Tidal and riverine processes are the dominant forces in the estuary and the action of the waves is considered negligible (Garel et al. 2009). The Guadiana Estuary is characterized by a semidiurnal mesotidal regime (Sampath & Boski 2016) with an average neap tidal range of 1.28 m and an average spring tidal range of 2.56 m (Garel et al. 2009). The western margin of the estuary is characterized by a salt marsh sheltered by a littoral sand spit that drains into the estuarine channel (Garel et al. 2009). The eastern margin is composed of barrier islands and sand spits separated by extensive salt marshes that drain to the Gulf of Cadiz through a tidal inlet (Garel et al. 2009). The flow rate from the Guadiana River varies drastically from less than 10 to 4,660 m³/s. The construction of over 100 dams for water storage and irrigation since the 1950s has strongly reduced the freshwater flow rate (Garel et al. 2009). Of importance is the Alqueva Dam located 60 km upstream of the head of the estuary (Garel & D'Alimonte 2017). This large reservoir was completed in 2002 and since then the freshwater flow into the estuary has been reduced from a yearly average of 143 to 16 m³/s (Garel et al. 2009). The reduction in flow rate also has an impact on residence time, the time it takes a particle to exit the system (Oliveira et al. 2006). The residence time of salinity for a river discharge of 8 m³/s ranges between 6 and 60 days, but when the river flow is high, the residence time ranges from a few hours to 9 days (Oliveira et al. 2006). Thus, the discharge of the Guadiana River has a large effect on the horizontal distribution of salinity in the estuary. The Guadiana Estuary is well-mixed when the freshwater discharge from the Guadiana River is low and is partially stratified for higher discharges. When there is a higher tidal amplitude and lower discharge from the Guadiana River, tidal processes control the water circulation of the estuary and the estuary becomes well-mixed (Garel & D'Alimonte 2017). The estuary is weakly stratified when it is neap tide and the river discharge is low (Garel & D'Alimonte 2017). In the latter case, density-driven processes control the mixing of the estuary. The water column is stratified only under extreme conditions (Basos 2013).

Hydrodynamics and water properties of the Guadiana Estuary were simulated with the Modelo Hidrodinâmico (MOHID, hereafter) water modeling system in response to different scenarios of sea level rise. MOHID is programmed in ANSI FORTRAN 95 in order to produce object-oriented models integrating hydrodynamic processes for different marine systems (MARETEC 2017). The present study used a two-dimensional model in MOHID, as the Guadiana Estuary is classified as a well-mixed estuary (Garel & D'Alimonte 2017), meaning that the water column is vertically homogeneous. The Guadiana Estuary can become stratified, but only under freshwater flows of 1,000 m³/s (Fortunato et al. 2002). The present study ran simulations for freshwater flows of 10 and 100 m³/s as these are typical conditions of the system after the closure of the Alqueva Dam. This further justifies the use of a two-dimensional model, allowing for long (1–2 month) simulations, which would have had prohibitive computational costs with a three-dimensional model. MOHID solves Navier–Stokes equations with the hydrostatic approximation, using the finite volume method over a generic geometric grid. In this work, hydrodynamic and Eulerian transport modules were used. The finite volume method allows for the transport equations to be applied to the entire cell volume at specific points in discrete time (Neves et al. 2000). Of relevance are the results of a study by Mills et al. (2019) as the present study builds upon the methodologies used to simulate hydrodynamics and changes in salinity distribution with respect to sea level rise. These authors used MOHID in two dimensions (2D, hereafter) to simulate the hydrodynamics of the Guadiana Estuary considering different sea level rise scenarios up to the year 2100 along with varying freshwater flow rates of the Guadiana River (Mills et al. 2019). Results of the model demonstrated an overall increase in salinity in the estuary as well as flooding area around the estuary with respect to sea level rise. Sea level rise led to a reduction in water velocity in the main channel, most likely due to an increase in depth and water volume (Mills et al. 2019). The model considered only the present bathymetry as well as a tidal amplitude equivalent to the average tidal amplitudes. The horizontal distribution of salinity in estuaries is dependent on several factors. One major variable is the balance between spring and neap tidal cycles and freshwater flow (Vargas et al. 2017). It is, thus, the aim of this work to expand upon the methodologies of Mills et al. (2019) to determine how the horizontal distribution of salinity varies between spring and neap tide as well as across a bathymetry that will allow flooding of the surrounding marshes.

The present study implements the same general setup of the model by Mills et al. (2019) and uses the same Cartesian computational grid of 1,400 × 350 cells with a resolution of 30 m. This computational mesh was chosen as it provides the most appropriate spatial resolution without incurring excessive calculation time. The calculation time for each simulation is high, especially when the river discharge is low. Two months of simulation time is required due to the high residence time when the freshwater flow rate is low (Oliveira et al. 2006). The hydrodynamic MOHID model of the Guadiana Estuary has been previously validated, calibrated, and used in several studies (Lopes et al. 2003; Morais et al. 2012). The present model consists of two separate bathymetries: (1) a bathymetry in which coastal management strategies are implemented to keep the coastline as it is (Figure 1) and (2) a bathymetry that allows for geomorphological changes caused by sea level rise and thus allows flooding around the estuary (Figure 2). The future state of the bathymetry of the Guadiana Estuary could fall anywhere between these two extremes, and thus, the present paper aims to analyze the differences between an unchanged coastline and a coastline that allows flooding. The first bathymetry is the current bathymetry of the Guadiana Estuary and was computed by triangular interpolation of measured bathymetric data on the Cartesian grid of 1,400 × 350 cells (Mills et al. 2019). The second bathymetry was computed by Sampath et al. (2011) in a behavior-oriented model. This bathymetry allows the surrounding saltmarshes and low-lying areas of the estuary to be flooded. It was solved based on the rate of sea level change, accommodation space for the deposition of sediments, and vertical accretion of sediments dependent on inundation (Sampath et al. 2011).

Various sea level rise cases were chosen based on future sea level rise projections from the 5th assessment report by the IPCC. The 5th assessment report includes sea level rise forecasts up to the year 2100 for different Representative Concentration Pathways (RCP), more commonly known as greenhouse gas emissions (Church et al. 2013). The study by Mills et al. (2019) considered values forecasted by the IPCC for a case of low greenhouse gas emissions (RCP 4.5) and a case of very high greenhouse gas emissions (RCP 8.5). Values were chosen based on the average predictions for the years 2045–2065, 2081–2100, and 2100 for both RCP values. The values imposed in the model include a sea level rise of 0.24 m for 2040, 0.48 m for 2070, and 0.79 m for 2100. Two freshwater flow scenarios were chosen based on typical annual conditions of the freshwater flow rate from the Guadiana River and include: (1) a low discharge of 10 m³/s representing the summer months when there is little rain and (2) a high discharge of 100 m³/s for the remaining months when there is more rainfall. These two values for river discharge are appropriate as the freshwater flow has been strongly regulated since the closure of the Alqueva Dam, yielding a reduction in the seasonal variability of the river discharge (Quesada et al. 2019). Two tidal amplitudes were imposed in the model to examine how the dynamics of the estuary vary at spring tide and neap tide. Since this study did not aim at representing a specific period, a simple M2 tide was imposed. The amplitudes representing spring and neap tides are 1.28 and 0.64 m, respectively.

Areas of inundation were computed for the various scenarios of sea level rise and highest freshwater discharge scenario (100 m³/s) over the bathymetry allowing flooding. The methodologies follow those of Mills et al. (2019) who computed flooding area as a function of the number of hours of land submersion during one tidal cycle, but only for a freshwater discharge of 500 m³/s. The present study includes an analysis of flooding area for a river discharge of 100 m³/s as well as each spring and neap tide scenario. Ten classes of flooding were computed as a function of submersion time over 1 tidal cycle (12.5 h). For example, class 1 corresponds to a submersion time of 0–1.2 h, class 2 corresponds to a submersion time of 1.2–2.4 h, and so forth. Histograms and flood distribution maps were computed based on the percentage of time a land cell was covered by water in the grid during one tidal cycle.

Time series graphs of water velocity, water elevation, and salinity were produced along with several locations of the estuary. The locations and nomenclature used for the time series analysis are shown in Figure 3.

Time series graphs were produced at each location shown in Figure 3, allowing for an assessment of the evolution of salinity over two tidal cycles (approximately 24.48 h) for each sea level rise scenario. The average change in salinity every 30 years is summarized in Tables 1 and 2. The simulations for salinity distribution assume a salinity value of 0 for freshwater and a value of 36 for seawater.

The following section examines the horizontal distribution of salinity along with water velocity direction at a time instant 1 h before high tide. As can be seen in the time series results, changes in salinity throughout the scenarios of sea level rise vary for the present bathymetry, whereas results from the alternate bathymetry reveal a correlation between sea level rise and salinity. Thus, all horizontal distribution maps of salinity are shown for the present bathymetry, whereas only the present year compared with 2100 are shown for the bathymetry allowing flooding.

See Figures 4 and 5.

See Figures 6 and 7.

See Figures 8 and 9.

Salinity classes were computed based on the amount of area covered by specific salinity values. A salinity class of 1 corresponds to salinity values between 0 and 5, a salinity class of 2 corresponds to salinity values between 5 and 10, and so forth. Salinity classes were computed for all tidal scenarios and freshwater discharge scenarios to examine the differences between the present and 2100. The results are shown in Figures 10 and 11.

See Figures 12 and 13.

See Figures 14–17.

The results obtained from the MOHID model have shown the dynamics of the Guadiana Estuary to be complex, especially with respect to the tides. For the bathymetry allowing flooding, salinity increase is correlated to sea level rise for all cases of freshwater flow and for each tidal scenario. On the other hand, the present bathymetry that does not allow the coastline to be changed results in a complex dynamic where salinity values vary at different locations with respect to sea level rise.

The horizontal distribution of salinity intrusion depends on the freshwater flow and tidal scenario. When the freshwater flow from the Guadiana River is high the marshes east and west of the Guadiana Estuary increase in salinity. Areas further upstream of the estuary portray an increase in salinity when the river discharge is low. These results are consistent with those from Mills et al. (2019) who also found larger changes in water velocity and salinity in the marshes near Esteiro da Leziria (381_147) on the western margin of the estuary for higher freshwater discharges compared with changes in water velocity and salinity in the main channel. Likewise, results from Mills et al. (2019) found larger changes in water velocity and salinity upstream of the main channel for low freshwater discharges. The bathymetry allowing flooding demonstrates increases in salinity much further into the marshes, with some areas reaching salinity values close to 36 by 2100. Water velocity magnitude and direction also justify several of the salinity results. Especially for the present bathymetry, decreases in salinity coincide with decreases in water velocity. This result is attributed to the deepening of the channel from the increase in water volume, which has been known to occur in other shallow ebb-dominated estuaries (Friedrichs et al. 1990). Salinity increases with respect to mean sea level rise. Although the results from the present bathymetry may not portray an increase in salinity at the points chosen for the time series, salinity does increase on the Spanish side of the estuary as shown in Figures 4 and 8. The bathymetry allowing flooding reveals a stronger correlation between sea level rise and salinity at each of the time series locations. Results obtained from this model are consistent with results from other numerical models which have found increases in estuarine salinity with respect to sea level rise (Hong & Shen 2012; Chua & Xu 2014; Vargas et al. 2017). The results of this model indicate that bathymetry, freshwater flow, and spring-neap tide variability impact the horizontal distribution of salinity intrusion caused by sea level rise.

Results from the salinity classes indicate that as sea level rises, a greater area of the Guadiana Estuary is covered by higher salinity values.

In terms of land inundation, all varying hydrodynamic factors result in an increase in inundation due to mean sea level rise. The area of inundation at neap tide is much smaller compared with the area of inundation at spring tide, which reflects how much sea level and tide contribute to the amount of land submerged. Furthermore, larger amounts of land are covered by higher flooding classes with respect to sea level rise.

The two-dimensional MOHID water model produced maps demonstrating the hydrodynamics and salinity distribution to gain further insight into the future state of the Guadiana Estuary as it responds to climate change. The results of the model demonstrate that the estuary responds differently to sea level rise based on freshwater flow, bathymetry, and tidal amplitude. All results portray an increase in salinity in response to sea level rise. The spatial distribution of salinity is dependent on the bathymetry as well as freshwater flow and tidal amplitude. When the freshwater flow is low in the spring and summer months, areas located upstream of the estuary increase in salinity. When the freshwater flow is high due to rainfall in the winter, the marshes west and east of the Guadiana Estuary in Portugal and Spain increase in salinity. A limitation of this work is the use of a two-dimensional model instead of a three-dimensional model. This was justified because the estuary is well-mixed, but there may be some variability in the water column due to the complex dynamics of the system, especially at neap tide. Additionally, a simple M2 tide was used, excluding the effects of overtides identified by Garel & Cai (2018) and Quesada et al. (2019). Future studies should use a three-dimensional model with real tidal signals to allow for a more complete evaluation of the Guadiana Estuary and how it responds to climate change.

All relevant data are included in the paper or its Supplementary Information.