Journal of Computational Science

Understanding the impact of sea level rise on coastal areas is crucial as a large percentage of the population live on the coast. This study uses computational tools to examine how two major consequences of sea level rise: salt intrusion and an increase in water volume affect the hydrodynamics and ﬂooding areas of a major estuary in the Iberian Peninsula. A 2D numerical model created with the software MOHID was used to simulate the Guadiana Estuary in different scenarios of sea level rise combined with different freshwater ﬂow rates considering varying tidal amplitudes. An increase in salinity was found in response to an increase in mean sea level in both high and low freshwater ﬂow rates at all areas around the estuary. An increase in ﬂooding areas around the estuary was also positively correlated with an increase in mean sea level.© 2020 The Authors. Published by Elsevier B.V. This is an open access article under the CC BY license


Introduction
A rise in mean sea level is a global concern as 10 % of the world's population live within 10 m elevation of the current sea level [3]. Based on projections by the Intergovernmental Panel on Climate Change (IPCC), the rate of sea level rise is accelerating. An acceleration in global sea level rise has been reported throughout the 20th century based on satellite altimeter data [6]. Global mean sea level rise could be as high as 1 m by the year 2100 if greenhouse gas emissions continue to be very high [5]. As the number of people accommodating coastal areas increases, an acceleration in the rate of sea level rise will severely impact society and the economy. Impacts of sea level rise on coastal areas include submergence of land, increased flooding, increased erosion, changes in ecosystems and an increase in salinity [16]. Not yet fully explored in the literature is how the dynamics of salinity in estuaries will evolve due to sea level rise. It is thus the aim of this work to quantify a relationship between sea level rise and salinity content in estuaries, as an increase in estuarine salinity has the potential to damage estuarine environments.
From a physical standpoint, an estuary is defined as a semienclosed body of water where river meets sea, extending into the river as far as the upper limit of the tides. Estuaries are productive environments, hosting highly diverse and complex ecosystems [20]. Fisheries, aquaculture, ecotourism, and port facilities rely on estuaries to contribute revenues to national, regional, and local economies.
These transitional areas where fresh water derived from land mixes with saltwater from the sea are sensitive to changes in climate such as sea level rise. Estuaries adapt to sea level rise by deepening their channels. This increases the accommodation space of sediments, further enhancing ebb and flood asymmetry [20].
A water column is considered stratified when less dense water lies over denser water. Thus, the stratification of an estuary is determined by the balance between buoyancy and mixing of denser saltwater from the sea and less dense fresh water from land. An increase in sea level allows more saltwater to enter the estuary [19]. Thus, a major consequence of sea level rise in estuaries is the intrusion of salt from sea water into fresh ground water, which can further impact stratification [15].
The intrusion of salt will have a profound effect on the physical properties of estuaries as the dynamic mixing of salt water and fresh water is the driving factor regulating stratification [19]. Changes in oceanic salinity combined with varying flow rates of freshwater discharge into the estuary will alter the entire estuarine circulation [19]. Hong and Shen [11] found an increase in stratification in the Chesapeake Bay as a result of an increase in salinity based on the results of a three-dimensional numerical model. This change in horizontal salinity gradients will further alter estuarine circulation and cause oxygen depletion [11]. Alterations in estuarine circulation will be detrimental to marine ecosystems that cannot tolerate high salinity content [4]. Furthermore, sea level rise increases salinity upstream and impacts tidal currents [11]. A study on the impact of an increase in salinity on Louisiana estuaries by Wiseman et al. [22] found both positive and negative trends attributed to the variability of freshwater flow rates. These authors found that an increase in salinity results in the death or decline in the productivity of marshes, which eventually leads to land loss [22]. Further impacts of salt intrusion include contamination of water supplies for both consumption and for the industry [10].
The aim of this article is to present the results of a numerical model simulating salt intrusion in the Guadiana Estuary, a major estuary on the Iberian Peninsula. The model simulates different sea level rise scenarios as projected by the IPCC along with different river discharge flow rates and various tidal amplitudes. As changes in sea level inevitably alter the topography of estuaries [17], the simulations of the present study were run over two possible future bathymetries of the Guadiana Estuary to determine if action must be taken to maintain the coastline or to allow flooding of surrounding areas of the estuary.

Study area
The Guadiana Estuary connects the Guadiana River to the Gulf of Cadiz. The Guadiana River is the fourth longest river on the Iberian Peninsula with a length of 810 km, 200 of which form a border between Portugal and Spain [7]. Beginning at its mouth in front of Vila Real de Santo António in the Algarve region of Portugal and Ayamonte, Spain, the estuary extends about 80 km north to its tidal limit at Mértola. The estuary drains a total area of 66,960 km 2 [9] in front of these highly populated regions. The Guadiana Estuary is at its widest in front of Vila Real de Santo António with a width of 800 m and is most narrow near Mértola with a width of 70 m. The average depth of the estuary is 5 m, but in some cases, it can reach 10 m.
Downstream, the estuary is in the form of a submerged delta where there is moderate wave energy [9]. The estuary is classified as a rock-bound estuary due to its narrow and relatively deep channel, located on a passive margin where the volume of water entering the estuary during the flood tide is larger than the freshwater discharge from the Guadiana River [9]. The Guadiana Estuary is characterized by a semi-diurnal meso-tidal regime with tidal influence 80 km upstream from where it meets the Gulf of Cadiz [7]. The average neap tidal range is 1.28 m and the average spring tidal range is 2.56 m with a maximum spring tidal range up to 3.44 m [9].
Tides and freshwater input control the vertical mixing and stratification of the estuary, which ranges from very well-mixed to partially stratified [9]. The mixing of the estuary varies by season, with a well-mixed water column when there is low discharge from the Guadiana River due to a lack of rainfall and a partially stratified water column when there is more discharge from the river. It is only under extreme conditions that the water column can become stratified [9]. The ebb currents of the estuary are generally faster than flood currents, even in low river conditions, an effect most likely due to the large hydraulic depth of the estuarine channel. The ebb-dominance becomes more pronounced as the river discharge increases [9].
The freshwater discharge of the Guadiana River can range between 10 m 3 /s and 4660 m 3 /s. The construction of over 100 dams since the 1950s has controlled 70 % of the drainage basin by the end of 2000 [9]. Of relevance is the Alqueva Dam, located 60 km from the head of the estuary. It is the largest reservoir in southern and western Europe storing 4150 hm 3 of water [9]. Its closure in 2002 has significantly reduced the mean river flow from an average of 143 m 3 /s over the last 26 years to 16 m 3 /s after 2002 [9]. Downstream the dam, the river has ephemeral flows and is dry for 40 % of the year [8]. The Guadiana Estuary currently faces a reduction in river flow due to scarce rainfall and an increase in freshwater storage. A map of the study area is shown in Fig. 1.

Hydrodynamic model setup and simulated scenarios
The MOHID water modelling system was used to simulate the hydrodynamics of the Guadiana Estuary and the effects of different sea level rise scenarios on salinity distribution and transport. MOHID is a three-dimensional hydrodynamical modelling system that integrates diverse numerical models in an objectoriented programming approach [13]. MOHID simulates physical and biogeochemical processes in the water column and sediments [13]. MOHID has been used to simulate the hydrodynamics of many coastal and estuarine systems, mainly the major estuaries of Portugal, but also in the Netherlands, France and Ireland [18]. The model has been previously validated and calibrated in the Guadiana Estuary by comparing in-situ measurements of salinity with model results [12]. A recent study by Quesada et al. [18] used MOHID in 2D barotropic mode to simulate the effects of freshwater discharge and tidal forcing on tidal propagation in the Guadiana Estuary. This model accurately produced phase values of tidal elevations and currents as confirmed by a harmonic analysis [18]. The model was able to demonstrate how freshwater discharge has a stronger effect on tidal wave amplitude upstream the estuary, whereas tidal forcing controls tidal wave amplitude downstream the estuary [18].
The model solves the Navier-Stokes equations considering the hydrostatic and Boussinesq approximations. u i is the velocity component in each direction; f is the Coriolis parameter, = 0 +' is the density decomposed in its reference value plus anomaly, and A i is the turbulence viscosity in each direction. The equations are discretized using the finite volume method. The solution scheme is semi-implicit of the Alternating Direction Implicit (ADI) type. The S21 scheme of Abbott et al. [1] for shallow water equations was used and adapted to 3D according to Martins et al. [14]. The equations are as follows: The same discretization and computational solution were used to transport salinity, which is assumed to be a passive tracer.
The calculation time to stabilize the salinity field is very high when the river flow is low due to the high residence time of the estuary in low river discharge conditions. For example, with a river flow of 10 m 3 /s it takes two months of simulated time for the salinity field to stabilize. The computational mesh was chosen to provide a spatial resolution appropriate for the study without incurring excessive calculation time. Thus, each simulation used a Cartesian computational mesh of 1400 × 350 cells with a constant space step of 30 m to evaluate salinity distribution and transport as well as coastal flooding at the mouth of the Guadiana Estuary. The simulations were run in an Intel Xeon Gold 6138 server and the code was parallelized using OpenMP and distributed over 8 cores. With this setup, two months of simulated time takes more than fifteen days. The model incorporates various sea level rise scenarios combined with different river discharge scenarios typical of the Guadiana River over two separate bathymetries. The first bathymetry was computed by triangular interpolation of measured bathymetric data on the mesh of 1400 × 350 cells with a space step of 30 m. This bathymetry is indicative of the present bathymetry of the estuary and corresponds to a case in which maximum human intervention is taken to keep the coastline unchanged. The second bathymetry represents a situation in which there is no human intervention, thus allowing maximum flooding of the saltmarshes and low-lying areas. Process-based models such as MOHID cannot evaluate long-term geomorphological dynamics so the results from a behavior-oriented model by Sampath et al. [21] were used for this specific bathymetry. The bathymetry was computed based on the rate of sea level rise, accommodation space for sediment deposition and vertical accretion of sediments due to inundation [21]. The two different bathymetries are shown in Fig. 2.
Two different freshwater discharges were imposed in the model representing: (1) a low freshwater discharge of 10 m 3 /s when there is little rain and (2) a high freshwater discharge of 100 m 3 /s when there is more rain. The model further explores three varying tidal amplitudes of type M2 representing: (1) typical system conditions equivalent to the average tide conditions, resulting in an amplitude of 1 m; (2) spring tide equivalent to an amplitude of 1.28 m and (3) neap tide equivalent to an amplitude of 0.64 m. To simulate the salinity transport in the Guadiana system, the salinity of water from the Guadiana River was considered null and the salinity of the coastal region adjacent to the estuary was set to a constant value of 36. These values are typical for this region. There are no units for salinity, in agreement with recommendations by UNESCO.
The two different river discharge scenarios mentioned above were combined with different sea level rise scenarios derived from sea level rise projections provided by the 5th IPCC report for Representative Concentration Pathway (RCP) 4.5 and RCP 8.5. Representative Concentration Pathway refers to greenhouse gas concentration [5]. The 5th report of the IPCC projected several possible values in global mean sea level rise for these two RCP values over the years 2046-2065, 2081-2100 and 2100. Due to the high computational cost and the fact that sea level rise scenarios would be combined with two different river discharge scenarios as well as three different tidal amplitudes, it was not possible to simulate all sea level rise forecasts from the IPCC. It was decided to consider three cases of sea level rise that covered all possibilities within all scenarios from the IPCC. The sea level rise projections correspond to a sea level rise of 0.24 m, 0.48 m, and 0.79 m corresponding to the years of 2040, 2070 and 2100, respectively. The values were imposed in the model by changing the reference level of the tidal components.
Three different sea level rise scenarios along with the current sea level were combined with the two river discharge scenarios along two possible bathymetries of the estuary. Additionally, in order to assess the future state of the Guadiana Estuary with respect to specific tidal events, spring and neap tide scenarios were also simulated. These boundary conditions on the model produce an evaluation of the previously mentioned scenarios resulting in projections, not forecasts of the future state of the Guadiana Estuary. While it is expected for the tides and river fluxes to change due to sea level rise, the aim of this study is to analyze the trends produced by sea level rise on spring and neap tides. Therefore, the results provide an evaluation of the sensibility of the Guadiana Estuary to changes in tidal amplitude combined with sea level rise.
Each sea level rise, river discharge and tidal scenario were imposed in the model to produce several time series and horizontal distribution maps. The locations for each time series can be seen in Fig. 1, corresponding to the intersection of the Guadiana Estuary and Esteiro da Leziria (381 147), in the middle of Esteiro da Leziria (391 96) and at the westernmost region of Esteiro Da Leziria (392 23).
The time series were constructed over a period of 2 tidal cycles (24.83 h) and represent the temporal evolution of average values of salinity, current velocity, and water elevation along with their  respective minimum/maximum values observed within the simulated scenarios. Average, maximum, and minimum values were obtained by grouping all three tidal amplitudes combined with both freshwater flow rates to examine the differences between the present and 2100. Scenarios based on the present bathymetry were analyzed separately from those simulated on the flood-permitting bathymetry.
Flooding areas were computed for ten different classes of submersion time during one tidal cycle over the bathymetry allowing maximum flooding. Table 1 outlines the different flooding classes.
For each class in Table 1 the area of flooding was computed, and a histogram was constructed for the various scenarios. Flood distribution maps were created to show the amount of land submerged during a tidal cycle according to flooding class.

Results
The results of the model are displayed as follows: section 3.1 includes the results of three time series locations on the western margin of the estuary in the saltmarshes of Castro Marim, section 3.2 includes the horizontal distribution of salinity as well as changes in salinity intrusion length in response to sea level rise, section 3.3 explores salinity distribution at spring and neap tide, and section 3.4 examines flooding areas due to sea level rise.

Horizontal distribution of salinity transport
The salinity distribution maps for the present bathymetry in high water conditions for the various scenarios of sea level rise can be seen in Fig. 6 and for low water condition in Fig. 7. Only salinity maps for a river discharge of 10 m 3 /s and 100 m 3 /s along with average M2 tidal amplitude of 1 m are shown as they are most representative of typical system conditions.
For scenarios with a low freshwater discharge (10 m 3 /s) and at high tide, the salinity progresses upstream in response to sea level rise in the main channel for the average M2 tide. At the inter-  section of the main channel and Carresqueira and Lezíria where the salinity has a value of 20, salinity progresses about 3100 m upstream from the present situation to the year 2040. The salinity intrusion after the year 2040 continues to progress upstream, but at a much slower rate. The evolution of salinity at the mouth of the estuary is less dramatic as mentioned in the time series analysis, but salinity intrusion can evolve by tens or hundreds of meters Fig. 8.
The distance of salinity intrusion at spring and neap tide for the two freshwater flow rates was also computed for both high water and low water conditions. Table 2 shows the increase in salt intrusion length from now to 2100 for both types of bathymetry.
All cases of tide, freshwater flow, and bathymetry result in an increase in salt intrusion length due to sea level rise. The increase is much more significant for the bathymetry allowing flooding compared to the increase in salinity for the present bathymetry. The increase in salinity intrusion length is also higher for the low freshwater discharge. The furthest extension of salinity intrusion is at spring tide with a freshwater flow of 10 m 3 /s where the salinity intrusion length increases from 26.31 km to 28.14 km upstream the main channel from now to 2100.

Salinity distribution at spring and neap tide
Salinity classes for spring and neap tide scenarios were computed based on the amount of area containing salinity values within the intervals as shown in Table 3.
The results demonstrate an overall increase in area covered by higher salinity values with respect to sea level rise for all cases of freshwater flow and tide. More explicitly, areas move up in salin-     ity class between the present and 2100. The differences are more pronounced at high tide and for lower freshwater flows.

Areas of flooding
Flooding area was computed as a function of the number of hours of submersion over one tidal cycle in the various scenarios of sea level rise as shown in Table 4. The analysis was carried out for the scenario of the highest freshwater discharge (100 m 3 /s).
The total intertidal submersion area increases from a value of 10 km 2 for the present year to a value of 14 km 2 in the "2100 scenario. The relationship between mean sea level and flooding area is nearly linear. Furthermore, an increase in mean sea level corresponds to a higher submersion time. Inundated areas also move  up in flooding class, meaning a higher amount of land will be submerged for a longer period with respect to sea level rise. The flood distribution maps (  Fig. 9) demonstrate that areas located further downstream the estuary are not as affected by mean sea level rise compared to areas located upstream. The evolution of the inundated area is more pronounced in the north margin of the inlet not confined by Vila Real de Santo António. The marshes near the terminal section of Ribeira do Beliche also suffer high inundation because of sea level rise. The present study does not have bathymetry data for locations further upstream this point, so no conclusions can be drawn for the flooding areas upstream. Since the river has a narrower upstream profile the areas of flooding are expected to be lower.

Discussion
The model simulations demonstrated that the mixing and transport of salinity are the result of a complex dynamic in which the lower portion of the Guadiana Estuary is dependent on both freshwater flow from the Guadiana River and the mean sea level of the adjacent coastal waters. The Guadiana Estuary alternates between purely fresh water when there is a high discharge from the Guadiana River and water of varying salinities for low flows from the river. In the latter situation, the tide dominates the salinity with a semi-diurnal pattern.
For lower freshwater discharge rates, salinity increased in response to an increase in mean sea level with a more pronounced change in salinity in the main channel. For higher discharges, the increase in water volume due to sea level rise induced a higher penetration of fresh water from the outlet and therefore resulted in lower salinity values at the mouth of the estuary. The velocity time series showed a reduction in advective transport at the mouth of the estuary due to mean sea level rise. Salt is transported mainly by diffusion, in the innermost regions where the signal of the tide is not detected. This result is consistent with other values of velocity obtained for this region. Time series results for locations further upstream the estuary revealed an increase in salinity of about 4 in response to sea level rise. The increase in salinity in response to an increase in mean sea level is more pronounced closer to the coastline than it is upstream the main channel. The increase is also much more pronounced for the bathymetry allowing flooding, with a nearly linear relationship between sea level rise and salinity content in the marsh area of Castro Marim. This is due to the extension of salinity further upstream the estuary when the simulations were run over the present bathymetry. Changes in salinity for the flood-permitting bathymetry show the impact the topography has on the hydrodynamics and salt transport of the estuary. The results from the specific locations used in the time series analysis can be used by coastal managers to determine if the increase in salinity at Castro Marim will impact the ecosystem, society, and the economy.
The results of the flood distribution areas for all types of tide demonstrated an increase in submerged land area due to mean sea level rise with an estimated increase of 4 km 2 for the entire intertidal area of the present situation. There is most significant flooding at Esteiro da Carrasqueira and the marshland surrounding Ribeira de Beliche, as shown inFig. 9. Results of the areas of flooding can be utilized by coastal managers and landowners to determine what action should be taken to manage the inevitable flooding due to sea level rise.

Conclusion
This study used a numerical model to analyze several impacts of sea level rise in the Guadiana Estuary. Sea level rise projections obtained from the 5th report of the IPCC for RCP 4.5 and RCP 8.5 scenarios for the years 2040, 2070 and 2100 were imposed in the model to simulate three different scenarios of sea level rise. The parameters of the simulations allowed an analysis of the velocity modulus, salinity distribution, and flooding areas of the Guadiana Estuary. Several variables such as tide, freshwater flow, mean sea level and changes in topography of the estuary were all considered for future outcomes of the estuary due to sea level rise. The high number of simulations allowed for an analysis of several possibilities of the future state of the Guadiana Estuary.
The model demonstrated that an increase in mean sea level generally results in a decrease in the velocity modulus of the main channel of the estuary. This is most likely due to the increase in depth of the main channel as well as an increase in water volume associated with the rise in mean sea level. Increases in salinity are most significant for cases of a low freshwater flow and for a more likely bathymetry that allows flooding of the surrounding areas.
The results portray an overall increase in both salinity and flooding area in the Guadiana Estuary with respect to an increase in mean sea level. As these two consequences are detrimental to estuarine ecosystems, society and the economy, further research must be done to further analyze the impacts of sea level rise in the Guadiana Estuary. One limitation of this study is that a two-dimensional model was used as opposed to a three-dimensional model. It is representative of most cases of the Guadiana Estuary because it is generally well-mixed, but as previously mentioned the estuary can become partially stratified in high freshwater flow conditions. A three-dimensional model would provide a more complete analysis of the hydrodynamics of the estuary. Additional numerical models should be simulated to portray a more quantifiable relationship between salinity and sea level rise and to further examine how the salinity changes with respect to other parameters, such as overtides.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.