Impact of Relative Sea-Level Rise on Low-Lying Coastal Areas of Catalonia, NW Mediterranean, Spain

Low-lying coastal areas are high-risk areas regarding sea-level rise (SLR) due to the combination of low elevation, subsidence and currently low sediment supply, on the one hand, and their high natural and socio-economic values, on the other. Because of this, there is a growing need to assess the integrated impact of SLR taking into account the interactions between physical and ecological factors to make more informed decisions for the design of adaptation strategies. This study presents a methodology for improved SLR-induced flood-damage assessments in natural areas which is composed of (i) a pseudo-dynamic method coupling equilibrium-based coastal response and classical bathtub modeling approaches and (ii) a simple method to account for habitat conversion following inundation. This methodology is applied to low-lying coastal areas of Catalonia under different SLR scenarios, which show very different sensitivities depending on the configuration of the water-land border, topography, geomorphology, and degree of human impact on the floodplain. In terms of potential impact, natural habitats will be the most benefited from the likely conversion at the expense of agriculture lands. This can help coastal managers to adopt adaptation strategies where considering the inherent capacity of some coastal landscapes to adapt opens up new alternatives.


Introduction
Sea-level rise (SLR) will significantly alter coastal landscapes through inundation, erosion and salt-water intrusion of low-lying areas worldwide. Considering that 10% of the world's population inhabits areas less than 10 m above sea level [1], the occupation of which has led to the widespread conversion of natural areas into economically productive regions [2], the most dramatic and immediate effects of SLR will be the inundation of coastal lowland areas [3]. In these areas, such as deltas, where accelerated rates of SLR are exacerbated by natural subsidence due to sediment compaction, inundation is likely by the end of the century [4]. Therefore, SLR coupled with subsidence rates (called relative sea-level rise (RSLR)) will increase the vulnerability of coastal communities and economic sectors to flooding in the near future, causing both environmental and socioeconomic changes [5,6].
Within this context, the Mediterranean coast is especially vulnerable to the impact of RSLR due to the high concentration of sensitive low-lying areas, anthropogenic pressures, and natural hazards [7][8][9]. Moreover, impacts related to climate and environmental changes will be more severe than the expected global mean, with temperatures already reaching +1.5 • C above pre-industrial times [10]. As a consequence, there exist a large number of risk assessments to RSLR for low-lying environments in this region, such as for the Po delta and other Italian plains [11][12][13], the French Mediterranean coastline [14], and the Egyptian coast [15] among others.  Table 1). (The geographic coordinate system is ETRS89/UTM zone 31N).
Low-lying areas along the coast have high environmental value due to habitats that provide high ecosystem service values (e.g., [17]) including a very large percentage of environmentally protected areas [42]; 47% of all the designated wetlands in Catalonia are located in the three selected study areas of the Gulf of Roses, the Llobregat Delta, and the Ebro Delta (hereinafter GR, LD, and ED, respectively) ( Figure 1) [43]. An overview of the main geomorphologic and socioeconomic characteristics of these areas is given in Table 1. These areas include salt, brackish, and freshwater marshes; coastal lagoons; beaches; and sandy dune habitat, which are protected under different legislative provisions including RAMSAR sites, the Natura 2000 Network, and Natural Parks. The ED is one of the largest deltas in the Mediterranean, of which approximately 7800 ha are protected as a Natural Park. In GR, the Aiguamolls de l'Empordà National Park encompasses 4730 ha in two discontinuous regions separated by the Empuriabrava marina, where human activities are strictly regulated. In contrast, Baix Llobregat Agricultural Park is dominated by peri-urban agriculture [44], Figure 1. The Catalan coast, coastal comarcas and wetlands (shadowed purple areas). The insets are the low-lying study areas and their digital elevation models (see Table 1). (The geographic coordinate system is ETRS89/UTM zone 31N).
Low-lying areas along the coast have high environmental value due to habitats that provide high ecosystem service values (e.g., [17]) including a very large percentage of environmentally protected areas [42]; 47% of all the designated wetlands in Catalonia are located in the three selected study areas of the Gulf of Roses, the Llobregat Delta, and the Ebro Delta (hereinafter GR, LD, and ED, respectively) ( Figure 1) [43]. An overview of the main geomorphologic and socioeconomic characteristics of these areas is given in Table 1. These areas include salt, brackish, and freshwater marshes; coastal lagoons; beaches; and sandy dune habitat, which are protected under different legislative provisions including RAMSAR sites, the Natura 2000 Network, and Natural Parks. The ED is one of the largest deltas in the Mediterranean, of which approximately 7800 ha are protected as a Natural Park. In GR, the Aiguamolls de l'Empordà National Park encompasses 4730 ha in two discontinuous regions separated by the Empuriabrava marina, where human activities are strictly regulated. In contrast, Baix Llobregat Agricultural Park is dominated by peri-urban agriculture [44], containing 2900 ha of Water 2020, 12, 3252 4 of 28 fruit and vegetable crops with the objective of promoting the integration of agricultural activities with the natural environment.
In the ED, approximately 210 km 2 of the coastal plain are devoted to rice production, generating approximately 98% and 13% of the total yield in Catalonia and Spain, respectively [45]. The land is cultivated under continuous paddy inundation, which requires a constant water supply from the Ebro River distributed via an extensive network of irrigation channels. In the LD area, the main crops are vegetables and fruits, most of which are consumed in the metropolitan area of Barcelona. Agricultural production in the GR is mainly cereals, fruit trees, vines, and olives. Table 1. Main characteristics of the study sites (see Figure 1).

Gulf of Roses (GR)
Llobregat Delta (LD) The main factor controlling the inundation of these low-lying areas is their topography (Table 1). Thus, ED is most vulnerable to SLR as approximately 53% of its land area is less than 0.5 m above the mean sea level (MSL) [18]. Anthropogenic influences also modulate the natural inundation of these areas. For example, the dense network of irrigation and drainage channels that crisscross the ED plain acts to extend the area of inundation across the deltaic plain [26].
The configuration of the coastline is also important since, given an absence of flood protection infrastructure, natural landforms serve as a natural barrier to inundation. The LD and GR have a coastal fringe formed by sandy beaches, with coastal dunes of moderate height [46] and some inlets and creeks. These areas are fronted by active shorelines able to respond to SLR and maintain some level of protection assuming they have enough accommodation space. On the other hand, the ED has an active sandy outer coastline with small areas of dunes and a passive muddy coastline along two semi-enclosed bays.

Low-Lying Areas
The topography of the study sites was characterized using a digital elevation model (DEM) with a grid resolution of 2 × 2 m obtained from Light Detection and Ranging (LiDAR) data from the Cartographic and Geologic Institute of Catalonia [47]. Land use and habitats were characterized using two databases: (i) a land-cover map of Catalonia developed by the Ecological and Forestry Applications Research Centre (4th version) and (ii) the habitat distribution maps produced by the Department of Environment of the Catalan Government (2nd version). The former is a high-resolution Water 2020, 12, 3252 5 of 28 thematic map obtained by photo-interpretation analysis with a scale of 1:2500 and a pixel resolution of 0.25 m [48]. The latter details habitats in general including those of interest at a European Union (EU) level compiled from aerial orthophotos (1:5000) between 2008 and 2012 based on the interpretation and adaptation of the EU CORINE classification [49,50]. The official Catalonia Wetland Inventory [51] was also used.

Sea-Level Rise
Sea-level rise projections were based on the Intergovernmental Panel on Climate Change (IPCC) AR5 RCP4.5 and RCP8.5 scenarios (50% probability level) [52], with used values representative of regional estimates as showed in [53]. In addition to this, we have also included a high-impact (H+) scenario to take into account the uncertainties associated with polar ice-sheet processes, which is given by the projection of sea level at 95% probability of the RCP8.5 steric component [54]. The inclusion of this scenario has been done from a risk management perspective to include very adverse conditions with a potentially high impact despite their low probability [55]. Relative to 2010, these three scenarios yielded 2100 SLR values of 0.49 m, 0.70 m, and 1.70 m, respectively ( Figure 2).
Water 2020, 12, x FOR PEER REVIEW 5 of 28 resolution thematic map obtained by photo-interpretation analysis with a scale of 1:2500 and a pixel resolution of 0.25 m [48]. The latter details habitats in general including those of interest at a European Union (EU) level compiled from aerial orthophotos (1:5000) between 2008 and 2012 based on the interpretation and adaptation of the EU CORINE classification [49,50]. The official Catalonia Wetland Inventory [51] was also used.

Sea-Level Rise
Sea-level rise projections were based on the Intergovernmental Panel on Climate Change (IPCC) AR5 RCP4.5 and RCP8.5 scenarios (50% probability level) [52], with used values representative of regional estimates as showed in [53]. In addition to this, we have also included a high-impact (H+) scenario to take into account the uncertainties associated with polar ice-sheet processes, which is given by the projection of sea level at 95% probability of the RCP8.5 steric component [54]. The inclusion of this scenario has been done from a risk management perspective to include very adverse conditions with a potentially high impact despite their low probability [55].  The three study sites are among the areas susceptible to subsidence along the Catalan coast [56]. The ED is most affected by subsidence, with reported rates ranging from 1.75 mm/yr [57] to 3 mm/yr [58,59]. The subsidence in ED under current conditions has been analyzed by [60] using DInSAr (Differential Interferometry Synthetic Aperture Radar). They found varying rates across the deltaic plain averaging about 3 mm/yr along the deltaic front, which we have selected as a representative rate. Subsidence rates for the LD and GR are poorly studied but lower than the ED, with reported average values of 1.25 mm/yr [61] and 0.8 mm/yr [62], respectively.

Methodology
Our methodology consisted of two steps: (i) the delineation of inundation-prone areas due to RSLR and (ii) the assessment of the impact of inundation on the affected habitats. In both cases, we adopt a pseudo-dynamic methodology that, in the first instance, accounts for the capacity of active shorelines to respond dynamically to RSLR and, in the second instance, considered the capacity for habitat conversion.

Delineation of Inundation-Prone Areas
The extent of the area susceptible to inundation under the SLR scenarios was delineated using a pseudo-dynamic bathtub approach. For this, areas below a given water level (i.e., the RSLR of interest) and hydraulically connected to the sea were assumed to be potentially inundated [22,23]. Hydraulic connection was defined using the "eight-side rule", where the diagonal and cardinal neighbors are used to determine if a model cell is flooded and to remove isolated low-lying inland areas [22]. However, this approach is valid only when the coastal zone is passive, that is, it does not The three study sites are among the areas susceptible to subsidence along the Catalan coast [56]. The ED is most affected by subsidence, with reported rates ranging from 1.75 mm/yr [57] to 3 mm/yr [58,59]. The subsidence in ED under current conditions has been analyzed by [60] using DInSAr (Differential Interferometry Synthetic Aperture Radar). They found varying rates across the deltaic plain averaging about 3 mm/yr along the deltaic front, which we have selected as a representative rate. Subsidence rates for the LD and GR are poorly studied but lower than the ED, with reported average values of 1.25 mm/yr [61] and 0.8 mm/yr [62], respectively.

Methodology
Our methodology consisted of two steps: (i) the delineation of inundation-prone areas due to RSLR and (ii) the assessment of the impact of inundation on the affected habitats. In both cases, we adopt a pseudo-dynamic methodology that, in the first instance, accounts for the capacity of active shorelines to respond dynamically to RSLR and, in the second instance, considered the capacity for habitat conversion.

Delineation of Inundation-Prone Areas
The extent of the area susceptible to inundation under the SLR scenarios was delineated using a pseudo-dynamic bathtub approach. For this, areas below a given water level (i.e., the RSLR of interest) and hydraulically connected to the sea were assumed to be potentially inundated [22,23]. Hydraulic connection was defined using the "eight-side rule", where the diagonal and cardinal neighbors are used to determine if a model cell is flooded and to remove isolated low-lying inland areas [22]. However, this approach is valid only when the coastal zone is passive, that is, it does not actively respond to RSLR. This includes areas of resistant geology and inner coastlines that do not experience direct wave action. In contrast, exposed low-lying areas fringed by sandy coastlines will react dynamically to changing sea-level conditions (e.g., [3]). Here, we assumed an equilibrium-based response of sandy coastlines to RSLR, as depicted by the Bruun model [63], i.e., an upward and landward translation of the active profile in-pace with rising sea level and maintaining the shape of the equilibrium profile. The validity (and uncertainty) associated to the use of Bruun or, in fact, of any model predicting SLR induced morphological changes has been covered in the literature by different authors [64,65], and it is still an open question [66,67]. In the absence of a universally accepted model, here we have used the Bruun model to assess the magnitude of the RSLR-induced shoreline retreat in the study sites. Thus, we followed the approach of Jiménez et al. [38] when analyzing SLR-induced erosion in Catalonia, whereby the predicted landward response is modulated or prevented by the existence or lack of accommodation space in the hinterland, respectively. Under this assumption, a sandy beach-dune system protecting a low-lying area under current conditions can migrate landwards, while maintaining its relative elevation (and thus protective function) under RSLR provided accommodation space (and sand) is available ( Figure 3). Following [38], to reduce the uncertainty in the selection of the closure depth (see [68]), the Bruun rule was applied using the characteristic values of the inner shelf slope for each study site, calculated from the shoreline to 10 m of depth and thus extending deeper than the medium-term closure of depth along the Catalan coast, which is about 7 m [69] ( Table 2). This permits to account for the expected increase of the limit of the active profile under increasing time scales as those corresponding to RSLR [70].
Water 2020, 12, x FOR PEER REVIEW 6 of 28 actively respond to RSLR. This includes areas of resistant geology and inner coastlines that do not experience direct wave action. In contrast, exposed low-lying areas fringed by sandy coastlines will react dynamically to changing sea-level conditions (e.g., [3]). Here, we assumed an equilibrium-based response of sandy coastlines to RSLR, as depicted by the Bruun model [63], i.e., an upward and landward translation of the active profile in-pace with rising sea level and maintaining the shape of the equilibrium profile. The validity (and uncertainty) associated to the use of Bruun or, in fact, of any model predicting SLR induced morphological changes has been covered in the literature by different authors [64,65], and it is still an open question [66,67]. In the absence of a universally accepted model, here we have used the Bruun model to assess the magnitude of the RSLR-induced shoreline retreat in the study sites. Thus, we followed the approach of Jiménez et al. [38] when analyzing SLR-induced erosion in Catalonia, whereby the predicted landward response is modulated or prevented by the existence or lack of accommodation space in the hinterland, respectively. Under this assumption, a sandy beach-dune system protecting a low-lying area under current conditions can migrate landwards, while maintaining its relative elevation (and thus protective function) under RSLR provided accommodation space (and sand) is available ( Figure 3). Following [38], to reduce the uncertainty in the selection of the closure depth (see [68]), the Bruun rule was applied using the characteristic values of the inner shelf slope for each study site, calculated from the shoreline to 10 m of depth and thus extending deeper than the medium-term closure of depth along the Catalan coast, which is about 7 m [69] ( Table 2). This permits to account for the expected increase of the limit of the active profile under increasing time scales as those corresponding to RSLR [70]. The landward extension of the active fringe adapting to RSLR was given as the reach of overwash-induced transport, which depends on the wave and water-level climates, and permits the beach to rebuild landwards while the shoreline erodes. This dynamic adaptation approach has previously been used to simulate the long-term (rollover) behavior of the Trabucador barrier system [71] and to explain the survival of the Tortosa barrier despite experiencing extreme erosion rates (approximately 40 m/year) in the ED [72]. The role of overwash transport in transferring material towards the hinterland is modulated by coastal morphology, particularly dune/beach elevation and beach slope (e.g., [73]).
To estimate the critical width, the extension of overwash deposits was identified and measured using aerial photos to obtain an averaged representative value for each study site (Table 2). Thus, in The landward extension of the active fringe adapting to RSLR was given as the reach of overwash-induced transport, which depends on the wave and water-level climates, and permits the beach to rebuild landwards while the shoreline erodes. This dynamic adaptation approach has previously been used to simulate the long-term (rollover) behavior of the Trabucador barrier system [71] and to explain the survival of the Tortosa barrier despite experiencing extreme erosion rates (approximately 40 m/year) in the ED [72]. The role of overwash transport in transferring material towards the hinterland is modulated by coastal morphology, particularly dune/beach elevation and beach slope (e.g., [73]).
To estimate the critical width, the extension of overwash deposits was identified and measured using aerial photos to obtain an averaged representative value for each study site (Table 2). Thus, in each case, beach profiles were assumed to adjust to RSLR as predicted by the Bruun rule provided the hinterland is wider than the projected shoreline erosion plus the critical width. If the existing accommodation space is bounded by an inner fixed boundary, this was assumed to be the innermost limit for beach migration. In the latter case, the elevation of the physical boundary controls inundation following the bathtub approach. Table 2. Representative values to calculate coastal fringe re-adaptation to RSLR following Bruun rule for each study site.

Potential Flood Damage
Once the areas susceptible to inundation under a given RSLR scenario were delineated, the resulting damage was estimated. One common approach is to determine the area of existing land use and habitats and assume they will be lost if inundated (e.g., [26]), which we refer to as the "total damage approach". However, when an inundated area is composed of natural habitats, some habitats have the capacity to change/adapt to new sea-level conditions, which we call the "conversion approach". The conversion approach is similar to that used in coastal landscape models where a given habitat affected by inundation has an assumed capacity to change/adapt depending on its properties and the new conditions to which it is exposed [27,28]. The evolution of habitats within an inundated area is evaluated by applying a series of conversion rules linking their current spatial distribution and their expected capacity to respond to each RSLR scenario.
First, coastal habitats were re-classified into generalized types based on their distinctive responses to SLR according to geomorphology, ecology, and level of development. These categories and general summaries of the land-cover types and main habitats found in each study site are provided in Table A1 (in Appendix A). The differences in habitats within a given category correspond to site specificities that depend on local characteristics such as geomorphology, freshwater availability, and land use.
The first variable affecting a habitat's susceptibility to inundation is vertical distribution and distance with respect to mean sea level. Here, we assumed that specific habitats thrive within a certain elevation range according to their salinity tolerance and resilience to flooding and physical (wave) disturbance. To this end, habitat and topographic data were jointly analyzed to evaluate their vertical distribution, with resulting distribution ranges varying among the study sites ( Figure 4), due to their local characteristics, particularly topography and water flux across the coastal plain (see also [74]). Considering the vertical distributions of existing habitats ( Figure 4 and Table A1), a series of basic transition rules were applied (Table 3). This followed a simplified version of previously adopted dynamic transitions including the SLAMM [30][31][32][33] and is similar to other simplified approaches [28].  An important aspect to consider in habitat response is the effect of hydraulic connectivity between the affected area and the sea. This was represented by delineating areas of total and partial inundation, which are inextricably linked to vegetation distribution [75]. To this end, we included tidal influences using the following three vertical levels: (i) the mean sea level (MSL), given by the target RSLR; (ii) the high-sea level (HSL); and (iii) the low-sea level (LSL). The HSL and LSL bound the intertidal zone around the future MSL, which corresponds to an average astronomical tidal amplitude of 0.20 m. Thus, when analyzing the inundation extent of each RSLR scenario, the pseudodynamic bathtub approach was applied as delimited by the MSL. Two additional inundation simulations were then performed, the first to derive the area below the HSL and the second to derive the area below the LSL. The area below the LSL was taken as the area of permanent inundation, and the intertidal zone between the HSL and LSL was taken as the area of temporary inundation. Importantly, to delineate the area below the HSL, hydraulic connectivity criteria were applied. This meant that areas within this elevation range but which are isolated by the presence of a physical barrier/infrastructure were not inundated. Habitat conversion under a given RSLR scenario was controlled by its location with respect to the MSL ( Table 3). The transition rules did not consider the time required for habitats to accommodate new conditions, but as a general rule, the longer they have to reach new conditions, the greater the capacity to change [76]. This time factor is further discussed in Section 5.2.  An important aspect to consider in habitat response is the effect of hydraulic connectivity between the affected area and the sea. This was represented by delineating areas of total and partial inundation, which are inextricably linked to vegetation distribution [75]. To this end, we included tidal influences using the following three vertical levels: (i) the mean sea level (MSL), given by the target RSLR; (ii) the high-sea level (HSL); and (iii) the low-sea level (LSL). The HSL and LSL bound the intertidal zone around the future MSL, which corresponds to an average astronomical tidal amplitude of 0.20 m. Thus, when analyzing the inundation extent of each RSLR scenario, the pseudo-dynamic bathtub approach was applied as delimited by the MSL. Two additional inundation simulations were then performed, the first to derive the area below the HSL and the second to derive the area below the LSL. The area below the LSL was taken as the area of permanent inundation, and the intertidal zone between the HSL and LSL was taken as the area of temporary inundation. Importantly, to delineate the area below the HSL, hydraulic connectivity criteria were applied. This meant that areas within this elevation range but which are isolated by the presence of a physical barrier/infrastructure were not inundated. Habitat conversion under a given RSLR scenario was controlled by its location with respect to the MSL ( Table 3). The transition rules did not consider the time required for habitats to accommodate new conditions, but as a general rule, the longer they have to reach new conditions, the greater the capacity to change [76]. This time factor is further discussed in Section 5.2. The conversion rules shown in Table 3 were applied assuming two boundary conditions: (i) no modification to existing structures in the floodplain (e.g., dikes, levees, channels, etc.) and (ii) negligible sediment input to the floodplain. The first condition could be changed to simulate the effect of floodgates or barriers as adaptation measures to RSLR. However, as our primary objective was to assess impact, we chose to evaluate the effect of inundation without anthropogenic influence other than the existing infrastructure. The second condition limits the inherent capacity of wetlands to adapt in response to SLR [75,76]. However, this assumption reflects current conditions in the region and especially in the ED (e.g., [58,77]). Thus, instead of considering vertical accretion rates to control habitat shifts [32,78], the applied rules were based on relative elevation, i.e., sea-level criteria and current vertical distribution of each habitat type.
Finally, the effect of saltwater intrusion on cropland areas was also considered, as this is a major threat to agricultural lands in the region [79,80]. Under this scenario, impacts are expected to differ between the study areas given that their main crop types are not the same (see [81] for the ED, [82] for the GR, and [83] for the LD), meaning variable salinity thresholds with respect to impacts on yield (e.g., [84,85]). Furthermore, local geomorphology dictates the suitability of an area for agriculture and crop type ( Figure 4). For example, given the low or very low tolerance of existing crops to saltwater intrusion, croplands are sparse in low-lying areas close to the water level in GR and LD (<1% of the total surface). Therefore, a minimum elevation threshold of 0.5 m with respect to the simulated future water levels was applied for the assumed maintenance of croplands; agricultural land below this level under a given SLR scenario was assumed to be abandoned due to salt intolerance and the likely rapid colonization by other salt-tolerant vegetation (see [86]). Such conversion of abandoned agricultural land to halophytic vegetation communities (or barren land) follows ecosystems connectivity, helping to preserve linkage between them [87,88]. The application of this agricultural transition rule was specific to GR and LD, as agriculture in the ED is devoted to rice production in paddies at all elevations along the delta plain (Figure 4), which are maintained by active freshwater inundation, which productivity also depends on elevation (e.g., [80]).

Flood Analysis and Potentially Inundated Areas
For comparison, Figure 5 shows the ED inundation extent for the entire range of RSLR scenarios based on both the classical (i.e., passive inundation) and pseudo-dynamic bathtub approaches. By 2050, including natural coastline adaptation to SLR decreased the inundated area by approximately 34.7% under RCP4.5 and RCP8.5, and by 21.1% under the H+ scenario. By 2100, the protection provided by the morphodynamic coastal response decreases, with the inundated area being just 11.3%, 5.1%, and 2% lower than that calculated by using the passive bathtub approach for RCP4.5, RCP8.5, and H+ scenarios, respectively. Therefore, the inclusion of the capacity for active sandy shorelines to respond dynamically to SLR yielded a smaller inundation extent relative to a passive inundation approach, and this effect was most pronounced under the low-medium RSLR scenarios. However, as RSLR increases, the natural protective effect decreased, reflecting site-specific modulation via anthropogenic modification of the deltaic plain ( Figure 6). Under lower RSLR conditions, most of the coast is protected by beaches, with inundated areas occurring (i) along the inner bays where there is no beach protection; (ii) in the inner part of the spits; and (iii) in areas with very flat beaches (<RSLR) (indicated by the light-blue areas in Figure 6). However, if a beach can maintain its relative elevation with respect to MSL-as the pseudo-dynamic approach assumed-the hinterland should be protected. This effect can only be detectable provided there are very low sandy stretches (i.e., below the target level); if the entire coast is higher than the simulated RSLR, the computed inundation extent of both methods will be the same.
The capacity for such a dynamic protective response is lost in some locations due to the lack of accommodation space, which is constrained by the existence of infrastructure in the hinterland. In the ED, this was mainly identified in the northern region, where some narrow beaches are backed by levees causing coastal squeeze (see the red square in Figure 6c), with the elevation of this infrastructure (with no adaptation capability to SLR) being the main control on inundation. Accordingly, approximately 1435 ha was flooded assuming 40 cm of SLR primarily as a result of the beach breaching in the absence of active adjustment. The distribution of floodwater across the deltaic plain is another element controlling the magnitude of inundation. This will largely be controlled by the topography of an area, but in the case of the ED, the existing network of channels crisscrossing the area extend the inundation from adjoining areas (Figure 6d). Water 2020, 12, x FOR PEER REVIEW 11 of 28 Light-blue shading indicates additional areas of inundation not accounting for dynamic beach response (i.e., computed using the static bathtub approach but not the pseudodynamic method). Active beach adaptation to RSLR is prevented when accommodation space is insufficient (e.g., the red square in (c) where narrow beaches are backed by infrastructure). Due to the network of channels crisscrossing the plain, some areas otherwise protected by a beach can be inundated by water from adjoining areas (red arrows in (d)). (The geographic coordinate system is ETRS89/UTM zone 31N).
The computed inundation-prone areas for the three study sites under the considered scenarios are shown in Table 4 and Figure 7. The results reflect the influence of geomorphology and relief, with the ED being the most vulnerable site, and GR and LD only significantly affected under high RSLR conditions due to relatively high beach profiles protecting the hinterland and their higher elevation. This higher topography also implies that differences in the inundation extents determined using the static and modified bathtub approaches will be relatively low in the GR and LD because, as previously noted, this is only detectable when beach/dune heights are lower than the RSLR scenario. (c) 40 cm; and (d) 50 cm. Light-blue shading indicates additional areas of inundation not accounting for dynamic beach response (i.e., computed using the static bathtub approach but not the pseudo-dynamic method). Active beach adaptation to RSLR is prevented when accommodation space is insufficient (e.g., the red square in (c) where narrow beaches are backed by infrastructure). Due to the network of channels crisscrossing the plain, some areas otherwise protected by a beach can be inundated by water from adjoining areas (red arrows in (d)). (The geographic coordinate system is ETRS89/UTM zone 31N).
The computed inundation-prone areas for the three study sites under the considered scenarios are shown in Table 4 and Figure 7. The results reflect the influence of geomorphology and relief, with the ED being the most vulnerable site, and GR and LD only significantly affected under high RSLR conditions due to relatively high beach profiles protecting the hinterland and their higher elevation. This higher topography also implies that differences in the inundation extents determined using the static and modified bathtub approaches will be relatively low in the GR and LD because, as previously noted, this is only detectable when beach/dune heights are lower than the RSLR scenario. Water 2020, 12, x FOR PEER REVIEW 12 of 28   The difference in vulnerability is reflected in the percentage of the affected surface. Thus, by 2050, the inundated surface in the GR and LD represents less than 1% and 2% of the total area under both RCP and H+ scenarios, respectively. In comparison, the proportion of the ED inundated under these scenarios increase by 9% and 27%, respectively (Table 4). Except for the ED, by 2100, the increase in the magnitude of RSLR under RCP4.5 and RCP8.5 did not result in a significant increase in the inundated area. Specifically, the long-term projection of the inundation extent in the GR and LD increased by less than 7% under RCP8.5. In comparison, in the ED, the affected area covered approximately 50% of the deltaic plain under RCP8.5. As expected, the vulnerability of these systems to the H+ scenario significantly increased, with approximately 35% of the GR and LD, and 80% of the ED, being susceptible to inundation. Figure 8 shows the relative percentage of habitats (by area) affected by inundation at the three study sites under the different RSLR scenarios. Although the evolution of affected surfaces with RSLR is similar, significant differences exist in the extent of inundation and the shape of the curve for each habitat type, largely reflecting their vertical distributions (Figure 4). The affected surface area for each habitat by 2050 and 2100 under the tested scenarios are given in Tables S1-S3 for GR, LD, and ED, respectively. The difference in vulnerability is reflected in the percentage of the affected surface. Thus, by 2050, the inundated surface in the GR and LD represents less than 1% and 2% of the total area under both RCP and H+ scenarios, respectively. In comparison, the proportion of the ED inundated under these scenarios increase by 9% and 27%, respectively (Table 4). Except for the ED, by 2100, the increase in the magnitude of RSLR under RCP4.5 and RCP8.5 did not result in a significant increase in the inundated area. Specifically, the long-term projection of the inundation extent in the GR and LD increased by less than 7% under RCP8.5. In comparison, in the ED, the affected area covered approximately 50% of the deltaic plain under RCP8.5. As expected, the vulnerability of these systems to the H+ scenario significantly increased, with approximately 35% of the GR and LD, and 80% of the ED, being susceptible to inundation. Figure 8 shows the relative percentage of habitats (by area) affected by inundation at the three study sites under the different RSLR scenarios. Although the evolution of affected surfaces with RSLR is similar, significant differences exist in the extent of inundation and the shape of the curve for each habitat type, largely reflecting their vertical distributions (Figure 4). The affected surface area for each habitat by 2050 and 2100 under the tested scenarios are given in Tables S1-S3 for GR, LD, and ED, respectively. The most affected habitats are those found at the lowest elevations, namely coastal lagoons and wetlands. Although this was common between the three sites, the greatest susceptibility was simulated in the ED, where more than 90% of the wetland area is inundated assuming a RSLR of 70 cm (broadly equating to 2100 under the RCP8.5 scenario). Under the same conditions, 38% and 16% of the existing wetlands in the GR and LD are inundated, respectively. This general trend also applies The most affected habitats are those found at the lowest elevations, namely coastal lagoons and wetlands. Although this was common between the three sites, the greatest susceptibility was simulated in the ED, where more than 90% of the wetland area is inundated assuming a RSLR of 70 cm (broadly equating to 2100 under the RCP8.5 scenario). Under the same conditions, 38% and 16% of the existing wetlands in the GR and LD are inundated, respectively. This general trend also applies to lagoons, although differences between the GR and LD are higher (Figure 8). This difference is even more pronounced for agricultural lands, where a 70 cm RSLR resulted in the submergence of approximately 50% in the ED compared to <1% in the GR and LD. Under the higher RSLR scenarios, agricultural lands in the GR and LD were increasingly affected but to a much lesser extent than in the ED.

Flood Damage Analysis
According to our methodology, habitats occupying inundated areas are not necessarily lost as they have some capacity to adapt to the new conditions. Tables S1-S3 show the evolution of habitat surfaces at each site under the different scenarios by 2050 and 2100. Natural habitats, which are located at the lowest elevations and are most directly affected by inundation, are also more likely to evolve or migrate in response to changing sea levels, and they are coastal lagoons, coastal vegetation, and wetlands. Figure 9 shows the projected habitat evolution in the three study sites by 2100 under the RCP8.5 scenario.
Water 2020, 12, x FOR PEER REVIEW 14 of 28 to lagoons, although differences between the GR and LD are higher (Figure 8). This difference is even more pronounced for agricultural lands, where a 70 cm RSLR resulted in the submergence of approximately 50% in the ED compared to <1% in the GR and LD. Under the higher RSLR scenarios, agricultural lands in the GR and LD were increasingly affected but to a much lesser extent than in the ED. According to our methodology, habitats occupying inundated areas are not necessarily lost as they have some capacity to adapt to the new conditions. Tables S1-S3 show the evolution of habitat surfaces at each site under the different scenarios by 2050 and 2100. Natural habitats, which are located at the lowest elevations and are most directly affected by inundation, are also more likely to evolve or migrate in response to changing sea levels, and they are coastal lagoons, coastal vegetation, and wetlands. Figure 9 shows the projected habitat evolution in the three study sites by 2100 under the RCP8.5 scenario.  With the adopted rules, the largest habitat variation was estimated for coastal lagoons whose surfaces increase with time and with the magnitude of SLR, and they occupy almost all of the projected inundated surface. By far, the largest increase occurs in the ED, where a significant portion of the deltaic plain was inundated, generating a large, shallow waterbody partially protected by the sandy coastal fringe, with different eco-geomorphological characteristics depending on the site. Along the inner northern and southern bays, saltwater open lagoons occur due to their passive shorelines. On the other hand, along the seaward coast, brackish-saltwater leaky lagoons occur, partially protected by a narrow sandy barrier ( Figure 9C). In the LD, the new lagoon surfaces would have a similar typology to the existing ones, mostly chocked lagoons elongated perpendicular to the coast ( Figure 9B). This is due to the relatively high coastal profile that results in inundation progressing through existing channels that connect lagoons to the open sea. Lastly, the simulated future lagoon surface in the GR has a combination of elongated chocked brackish lagoons at the northern end of the bay and restricted/leaky lagoons in the central part ( Figure 9A).
For wetlands, the largest conversion was also projected in the ED. At relatively low RSLR values (RCP8.5 by 2050), an increase of approximately 25% was predicted with losses due to full submergence (i.e., conversion to lagoons) compensated by the conversion of rice fields and areas of coastal vegetation in the intertidal zone to transitional wetlands. Importantly, approximately half of this projected area would not be captured without considering the dynamic adaptation of habitats (i.e., the "total damage approach" in Table S3). Under higher RSLR values (RCP8.5 by 2100), the wetland surface decreases to approximately 79% of the original values. In the LD and GR, the simulated variations are different from the ED as topography significantly differs ( Figure 9). Thus, at low RSLR values (RCP8.5 by 2050), change in the wetland surface is negligible. With a higher rate of change (RCP8.5 by 2100), the wetland surface area in the GR and LD decreased to 71% and 90%, respectively (Tables S1 and S2).
Coastal vegetation shows a similar variation to wetlands, where a new fringe of halophyte vegetation developed close to the sea level affected by high tides (the HSL-MSL range). The largest variation was found in the ED, where a significant increase of up to 240% was simulated relative to the current state under low SLR scenarios (RCP4.5 and RCP8.5 by 2050). However, under the higher SLR scenario (e.g., RCP8.5 by 2100), although the area occupied by coastal vegetation is larger than under current conditions, it progressively decreased due to the effects of existing infrastructure. This change in areas of coastal vegetation was also observed in the GR and LD but was modulated according to their topographic and spatial characteristics. These factors determine the extension of halophytic vegetation in comparison with losses via direct inundation and landward beach migration, whereby current coastal vegetation could be buried by overwash deposits (Tables S1 and S2). The simulations also indicate that the composition of new areas of wetland and coastal vegetation would differ from the current communities via large-scale conversion to transitional wetlands, and halophytic communities would occur at the expense of freshwater vegetation and agricultural land.
Beach and dune habitats were simulated to decrease at all sites, with spatial patterns driven by the effect of existing infrastructure that limits landward migration. Thus, considering the expected change by 2050 under the RCP8.5 scenario, beaches decrease to approximately 96%, 89%, and 86% of their current area for GR, LD, and ED, respectively. As expected, the loss of beaches significantly increases beyond 2050 due to SLR acceleration and the lack of accommodation space for natural adaptation (Tables S1-S3).
Finally, agricultural lands were simulated as being the most negatively affected habitat, with estimated losses larger than those resulting from direct inundation given no capacity to adapt and because the fringes closest to the projected water levels are developed preventing their future use for agriculture. Thus, croplands currently occupying the future intertidal zone under a given RSLR scenario are likely to be replaced by halophytic vegetation and transitional wetlands. Moreover, in the GR and LD, agricultural land located at elevations below +0.5 m with respect to future MSL is too saline to support the current cropping systems. For example, under a relatively low RSLR scenario (2050 under RCP8.5), losses of cropland due to habitat conversion are 1.5%, 2.6%, and 5.1% above those from direct inundation in the GR, LD, and ED, respectively. Under the same scenario, by 2100 the relative increases in cropland losses are 8.9%, 26.1%, and 3.3% for GR, LD, and ED, respectively. This reflects the potential abandonment of the orchards and crops on the higher land in the GR and LD due to the saltwater intrusion.

Methodological Aspects
We developed a pseudo-dynamic method to account for active shoreline responses to RSLR in predictions of the inundation of low-lying areas of sandy coasts. This approach is based on the coupling of an equilibrium-based coastal response and classical bathtub modeling approaches. Although simple, it permits the movement beyond a passive inundation model, which ignores the dynamic response of coastal environments to RSLR (e.g., [21]).
The dynamic adaptation of sandy shorelines to RSLR has been simulated by using the Bruun model, which was implemented considering the availability of accommodation space in the hinterland. When no obstacles exist, beaches adjust by migrating upward and landward so that their relative elevation is maintained and, consequently, afford the same level of protection against inundation. However, where a physical barrier prevents this landward migration, beaches will progressively erode and lose their dynamic protective capacity, and inundation will ultimately be controlled by the height of the barrier. Thus, a lack of accommodation space depends not only on the existence of a barrier but also on the rate of SLR, which determines the velocity of beach landward migration to reach the barrier. Our approach can be adapted for other models of coastal response to SLR (e.g., [68,89,90]), with the protection afforded by a beach ultimately controlled by the dynamic response to future conditions. It has to be noted that we have not included shoreline changes other than RSLR-induced, and this should be equivalent to "isolate" the RSLR component in the long-term behavior of these areas. However, other factors such as river sediment supplies and longshore and cross-shore sediment transport patterns would also contribute to their long-term evolution (e.g., [59,71]). Integrating all these components into a single long-term morphodynamic model is an issue that is far from being solved at present (e.g., [66]). Importantly, our approach is closer to reality for assessing the inundation extent in sandy coasts, being relevant for low-lying beaches lower than the target RSLR. Thus, when we compared the computed inundated extents for the three study areas with those obtained by using the conventional bathtub approach, significant differences were only found for the ED (a <35% in the inundated area by 2050 under the RCP8.5 scenario). Although this trend should be maintained for higher RSLR provided accommodation space is available, both methods tended to produce similar results, with a reduction in the inundation extent of just about 5% by 2100. However, this is largely attributed to the hydraulic connectivity of the ED plain, where seawater entering through passive areas is distributed across the plain via the channel network (see also [26]).
We also implemented a simple method to account for the possibility of habitat conversion following inundation. Thus, SLR-induced damage projects are modulated and, as a consequence, may affect decision-making and the development of adaptation strategies. In most projections, an assumption is made that all inundated areas of land are lost. While this can be acceptable for highly modified environments, this may lead to unrealistic damage estimations for those natural systems able to adapt to future conditions.
There are some advanced tools to simulate habitat conversion, such as the SLAMM, which was specifically used to characterize wetland response under future sea-level scenarios [31][32][33][34]. This model allows for the dynamic assessment of SLR, transitioning away from simple and static inundation models, although concerns regarding application in different geomorphological, ecological, and hydrological settings have been raised [91,92]. Nevertheless, the use of complex landscape models provides detailed information on habitat responses despite the constraints of time-consuming computational effort and model calibration [35]. Here, we followed a much simpler approach based on the same underlying concept-that habitats have a capacity to adapt to new conditions-formulated as a series of rules based on ecological succession, tolerance to inundation and/or salinity, vertical distribution, and location with respect to the MSL. This approach considers the intertidal zone associated with the future MSL as a buffer area of potential high environmental value. When located within the flood plain outside the sandy fringe, areas of typical vegetation will naturally evolve, with habitat shifts being controlled by hydraulic connectivity to the sea and salinity tolerance (e.g., [93]). This is similar to the assumptions of some existing models for wetlands that predict vegetation distribution based on dominant hydrodynamic conditions, such as water depth and duration of inundation [94,95].
Although the capacity for wetland adaptation is highly linked to sediment availability and vertical accretion [27,96,97], this was not considered here as sediment input to the study areas is practically zero. For example, in the case of the ED, significant retention of riverine sediments behind dams in the upper Ebro River catchment (up to 99%) dramatically reduces the capacity for the vertical accretion of the delta plain in response to SLR [58,98]. Indeed, the suspended sediment load is currently <0.01 g/L [99], meaning the delta is unable to grow and suffers from intense coastal reshaping by wave action [59,100]. This anthropogenic impact on hydrology and sediment budgets can disrupt the eco-geomorphological feedbacks needed to adjust to SLR and hinder the natural capacity of systems to survive [101,102]. Therefore, the conversion rules we applied in our simulations are appropriate for the target study sites but should be adapted for application elsewhere if conditions permit vertical accretion. Overall, although a broadly simple approach, our methodology improves the impact assessment of the RSLR-induced inundation in low-lying natural areas both with respect to inundation extent and potential habitat loss. Importantly, the approach can be easily implemented in a GIS environment and, thus, used as an additional element in the decision-making process to design adaptation strategies in this type of environment.

Inundation-Driven Impacts on Study Sites
While the study sites are considered the most sensitive areas of the Catalan coast to SLR-induced permanent inundation, they show very different sensitivities depending on the configuration of the sea-land border, topography, geomorphology, and degree of human impact on the floodplain. The ED was found to be most vulnerable to SLR inundation due to (i) a very long passive coastline unable to protect the floodplain behind along the semi-enclosed bays, where floodwaters affect a large area of the hinterland; (ii) a very low relief, which exposes a large proportion of the floodplain to inundation under relatively low RSLR scenarios (approximately 73% of the plain is <1 m above MSL); and (iii) a dense network of channels crisscrossing the plain that, although may act as barriers, if not managed appropriately, can facilitate the transfer of floodwaters across the plain. In contrast, the LD and GR show more resilient configurations due to (i) coastlines formed by relatively wide sandy beaches with dunes able to adapt to changing sea levels while maintaining a good level of protection; (ii) limited exposure of the lowest part of the plains due to their topography (approximately 15% of their plains are <1 m above MSL); and (iii) the main areas of water entrance onto the plain are existing outlets connecting lagoons and rivers with the open sea, which restricts inundation in these areas.
In this context, it is important to highlight the need for reliable estimates of local subsidence rates, as these will determine the local acceleration of SLR and, consequently, modulate possible differences in RSLR under the same climatic scenario, as those reported by [103,104] for the Mediterranean coast. In the analyzed areas, used subsidence rates are based on existing local estimates, with GR and LD being the sites with a lower coverage and, consequently, with a larger associated uncertainty. Differences in vulnerability among study sites are also reflected in the expected impacts of SLR-inundation on representative habitats. For croplands, the ED is considered the only site to be significantly impacted by inundation based on both the relative and absolute magnitudes of the affected area. This is due to the extension and covered range of elevations of the deltaic plain surface, which indicated a direct impact (by inundation) on cultivated land of between 10% and 52% by 2050 and 2100, respectively, under RCP8.5. Existing studies on the impact of SLR on rice production due to salt intrusion have found a similar pattern [81] due to the relationship between soil salinity and surface elevation. Considering these results, it is expected that the productivity of agriculture land above the simulated inundation levels will also decrease due to progressive soil salinization.
The simulated impact of inundation on agriculture is much smaller in the LD and GR compared to the ED, both in absolute and relative terms, with approximately 26% and 10% of agricultural land affected, respectively, by 2100 under RCP8.5. It is important to note that these projections only account for direct inundation and do not reflect areas affected by saltwater intrusion, which would lead to higher losses in the GR and LD. For the LD and GR, as most of the agriculture lands affected by SLR would become disconnected from water, cultivation will likely be abandoned. This implies a net loss of ecosystem services provided by agricultural lands (e.g., [105]), including cultural services [82,83]. In the ED, the inundation of croplands maintained hydraulic connectivity in our simulations, driving the development of a shallow lagoon along the entire coast ( Figure 9). Thus, the combination of a wide and gently sloping topography fringed with coastal vegetation, lagoons, and wetlands at the seaward boundaries of rice fields and protection from direct wave action afforded by an active sandy coastline provides ideal conditions for marsh development and migration. Indeed, in addition to vertical accretion, transgression into adjacent uplands is a primary mechanism for marsh survival, where gently sloping uplands favor marsh migration (e.g., [106]).
The current configuration of the ED plain appears to be the most important factor determining differences in the degree of inundation as well as the induced damage. On the one hand, the inundation extent is strongly controlled by the dense channel network crisscrossing the delta plain [25]. On the other hand, this network moderates habitat migration and the formation of new natural areas. To enhance this landscape shift from agriculture to natural habitats formed by lagoons and marshes, existing barriers in the plain (e.g., channels and other minor infrastructure) should be removed, as they function as the main obstacle for marsh expansion by upland migration (for example, [107,108]). This boundary effect has also been considered by [36] who analyzed the effects of SLR on the pristine and anthropic configurations of the ED.
The largest expansion of natural habitats was predicted for lagoons. In the case of the LD and GR, lagoons increased in area but maintained their typology by expanding from their current configuration over the surrounding wetlands according to local topography. This topographic control means that the LD is expected to be the site experiencing the lowest degree of lagoon expansion-approximately 184% by 2100 under RCP8.5 compared to 600% for the GR. In absolute terms, the simulated lagoon surface in the ED under this scenario will be approximately 0.16 times the current one. In contrast, once the critical rate of SLR at which wetlands drown is exceeded, the extensive inundation of the ED plain would result in a significant increase in lagoon surface area (approximately 565% by 2100 under RCP8.5) and significantly altered typology. A large proportion of this new surface would be occupied by saltwater open shallow lagoons along the bayside shorelines and brackish-saltwater leaky lagoons along the seaward coast. Under current conditions, freshwater inputs to lagoons occur via irrigation channels, a practice that can be maintained under future conditions. Importantly, the ecological characteristics of converted natural areas may differ from the existing ones. Projected wetland areas would be dominated by saltwater ecosystems, whereas under current conditions, they are a mix of saltwater, brackish, and freshwater environments. In this context, salinity is generally one of the most important factors determining coastal wetland habitat type (e.g., [109]), and for the three study areas, SLR-related inundation will drive a shift in plant communities to saltand flood-tolerant wetlands [110,111], having implications for wetland management [109].
The time available to adapt becomes crucial when considering the ability of habitats to naturally respond to changing sea levels. We assessed the impact of inundation as a function of projected water levels but did not consider the time taken to reach these levels. Indeed, the rate of rising will strongly affect the capacity of these systems to respond. Thus, the longer the time taken to reach new conditions, the higher the likelihood that a system can respond via ecological succession. This temporal effect was examined for coastal marshes by [27], demonstrating threshold SLR rates that lead to marsh submergence and irreversible conversion to unvegetated areas. Such thresholds are site-specific and depend on factors such as sediment availability and tidal range (e.g., [106,112]). This potential effect is illustrated for the ED in Figure 10, which shows the inundation extent for three elevation ranges The inclusion of time as a variable in assessments of potential inundation damage is important because it has direct consequences in terms of (1) flood management, by constraining the design, development, and implementation of adaptation strategies, and (2) habitat response, by increasing the probability of adaptation and ecological succession.

Implications for Designing Adaptation Strategies
The approach adopted to assess the impact of RSLR in low-lying areas provides an alternative vision to considering the entire inundated surface as a total loss of value. The classical total damage approach results in management strategies based primarily on flood control. This is the case of the ED, where rice producers have developed an adaptation plan that includes building a dike along the bayside shorelines and at the back of beaches to prevent inundation of the deltaic plain and, thus, the inundation of rice fields [113]. However, even where such a barrier is advantageous, rice fields within the potential inundation area are typically below future projected water levels and, consequently, will be severely affected by saltwater intrusion (e.g., [81]). Thus, to maintain their productivity, additional investment to adapt the existing hydrological infrastructure is also required, such as pumping systems to remove water from low-lying fields. In the case of the LD and GR, this strategy is not likely plausible, given that the extent of habitats of economic interest is very small, unless the H+ scenario is considered. However, even in this extreme case, the damage is expected to be relatively low.
The adoption of an adaptation strategy based on protecting the anthropogenic system would artificially bind existing natural areas, such as lagoons and wetlands, disconnecting them from the sea, and further limiting accommodation space along the coastline. This implies a strategy based on poldering, although alternative nature-based strategies are increasingly being considered in areas that traditionally use this approach, including removing existing barriers to enhancing natural values (e.g., [114,115]). In fact, given there is no "silver bullet" solution to adaptation, a transformative

Implications for Designing Adaptation Strategies
The approach adopted to assess the impact of RSLR in low-lying areas provides an alternative vision to considering the entire inundated surface as a total loss of value. The classical total damage approach results in management strategies based primarily on flood control. This is the case of the ED, where rice producers have developed an adaptation plan that includes building a dike along the bayside shorelines and at the back of beaches to prevent inundation of the deltaic plain and, thus, the inundation of rice fields [113]. However, even where such a barrier is advantageous, rice fields within the potential inundation area are typically below future projected water levels and, consequently, will be severely affected by saltwater intrusion (e.g., [81]). Thus, to maintain their productivity, additional investment to adapt the existing hydrological infrastructure is also required, such as pumping systems to remove water from low-lying fields. In the case of the LD and GR, this strategy is not likely plausible, given that the extent of habitats of economic interest is very small, unless the H+ scenario is considered. However, even in this extreme case, the damage is expected to be relatively low.
The adoption of an adaptation strategy based on protecting the anthropogenic system would artificially bind existing natural areas, such as lagoons and wetlands, disconnecting them from the sea, and further limiting accommodation space along the coastline. This implies a strategy based on poldering, although alternative nature-based strategies are increasingly being considered in areas that traditionally use this approach, including removing existing barriers to enhancing natural values (e.g., [114,115]). In fact, given there is no "silver bullet" solution to adaptation, a transformative change in current policies is needed that goes beyond the existing portfolio of accepted approaches [116].
Addressing maladaptation in coastal areas may have devastating effects and, sometimes, may cause more damage than a "do nothing" approach in which nature is allowed to take its course [117]. In areas of high environmental value, ensuring accommodation space for natural habitats to adapt may be a feasible strategy to cope with the effects of SLR [118]. However, urban sprawl and the existence of flood-prevention structures can act as barriers that inhibit such dynamic ecosystem responses [119]. Promoting natural protection by creating natural buffer areas must be accompanied by the progressive removal of anthropogenic infrastructure that prevents migration, thus favoring ecosystem connectivity. This may be feasible in the ED because of its low elevation and flat topography, since landward migration of wetlands is largely controlled by the topographic slope of the adjacent upland [106]. This could be achieved by reclaiming stretches of agricultural land located at the lowest elevations, which are the first to be affected by inundation and saline intrusion, to be transformed into wetlands and, in turn, restrict urban infrastructure and rice fields to higher elevations on the plain. This approach to managing inundations will result in a dynamic coastal landscape where habitats would shift according to future conditions and where the existing channel network could be used to control salinity. This strategy is not usually well received by local stakeholders set to be most economically affected [120][121][122] in general, and in the ED in particular [123,124]. In the GR, potential social conflicts exist between different stakeholders in the context of environmental protection, largely as a result of perceived threats to agriculture [125]. The existence of multiple stakeholders with different perspectives is always a challenge in adaptation planning, with "winners" favoring adaptation and "losers" opposing it [126]. Nevertheless, the creation of new wetland areas may support alternative economic activities focused on environmental conservation, which could be more sustainable given the high environmental value of the area (e.g., birdwatching, cycling, eco-tourism, etc.) [43,127]. Further understanding of the perceptions of all stakeholders is required during the design of sustainable adaptation actions in response to climate change [128]. When addressing SLR-induced flooding in highly vulnerable areas, such as lowlands, quantifying affected areas as direct losses is justified assuming no change in the status quo; however, considering the inherent capacity of some coastal landscapes to adapt, SLR may present new opportunities. In Catalonia, wherever possible, promoting habitat creation may help increase coastal resilience in the face of SLR and provide additional ecosystem services in environmental hotspots, even where these are progressively degraded by human actions.

Conclusions
In this work, we have presented a methodology for improved SLR-induced inundation-damage assessments in natural coastal areas. To improve the delineation of inundation-prone areas, the classical bathtub approach was integrated with an equilibrium-based coastal response to RSLR to account for the dynamic adaptation capacity of sandy shorelines. To improve the assessment of induced damage in inundated areas, the likely ability of coastal habitats to adapt to changing conditions was simulated by linking their spatial distribution with their expected ability to respond to changing sea-level conditions. This methodology was applied to assess the impact of different scenarios of RSLR in the most vulnerable low-lying areas of Catalonia, which also have the highest natural values along the coast.
Obtained results showed a very different susceptibility of the study sites to inundation, being the ED the most affected one under all considered scenarios. In terms of associated impacts, the most affected habitat was croplands, which is unable to adapt to new conditions, and it will convert to natural habitats (wetlands/lagoons) depending on local hydraulic connectivity. In this sense, existing man-made infrastructures in the ED plain play a double role. On the one hand, the existing network of channels controls the inundation extent on the plain, and on the other hand, they act as physical barriers limiting the horizontal habitat migration and the formation of new natural areas. The adopted approach permits to move beyond traditional inundation-damage assessments where any affected area is considered as a total loss of value. In areas of high natural values, this may lead to consider SLR not only as a threat but also as an opportunity for a change in their management that allows a range of adaptation strategies different from classic protection measures.
Supplementary Materials: The following are available online at http://www.mdpi.com/2073-4441/12/11/3252/s1, Table S1: Projected surface (Ha) under the TDA and CA for different SLR scenarios at 2050 and 2100 in the GR. Table S2: Projected surface (Ha) under the TDA and CA for different SLR scenarios at 2050 and 2100 in the LD. Table S3: Projected surface (Ha) under the TDA and CA for different SLR scenarios at 2050 and 2100 in the ED.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.