Characterization of hypersaline zones in salt marshes

Salt pans are commonly found in coastal marshes and play a vital role in the marsh plant zonation. However, the correlation between these hypersaline zones and the marsh hydrological conditions have barely been characterized. This study numerically investigates the effects of evaporation rate, tidal amplitude, and marsh platform slope on salt pan formation, and found that salt pans can hardly grow in the intertidal zone due to regular tidal flushing, while tend to form in the lower supratidal zone, where evaporation is sustained. The accumulated salts create an upward salinity gradient that trigger downward unstable flow. The decreases of potential evaporation rate, tidal amplitude and/or marsh platform slope strengthen the hydraulic connection between the marsh surface and the underlying watertable, the key to sustaining evaporation, and therefore result in thickener and wider salt pans. These findings offer a deeper insight into the marsh eco-hydrology and guidance for their degradation prevention.


Introduction
Salt marshes at the ocean-land interface are one of the most productive ecosystems worldwide, since they provide numerous eco-functions, e.g. maintaining coastal biodiversity (Adam 1990, Gopi et al 2019 and protecting inland area from the threat of sea level rise (Kirwan andMegonigal 2013, Fagherazzi et al 2020). These eco-functions are closely linked to marsh plants, which have been found to exhibit striking zonation patterns across elevational gradients in many salt marshes (Fariña et al 2018). Various studies were carried out to understand the underlying mechanisms and major drivers of plant zonation (e.g. Silvestri et al 2005, Moffett et al 2010a, Wilson et al 2015, Carol et al 2019, Xie et al 2020. Generally, it is believed that physical stress, rather than biotic interactions, greatly affects the distribution pattern of marsh plants (Bertness and Pennings 2000, Emery et al 2001, Feng et al 2018. For example, He et al (2012) found that, due to the harsh physical conditions of high estuarine marshes in the Yellow River Delta, the invasion of Spartina is limited to low estuarine marshes only.
Salinity is an important physical stressor controlling plant zonation, because different species vary greatly in their salt tolerance levels. For low-latitude salt marshes, salinity may determine the distributional limit of Spartina (Pennings et al 2005, Qi et al 2017. While salinity changes were rated by Moffett et al (2010a) as the most effective metrics for distinguishing species habitats and halophyte zones in the marsh site they investigated. The hypothesized link between salinity variations and zonation has driven various studies in the forms of field measurements, laboratory experiments and numerical modeling. Laboratory experiments under controllable conditions have demonstrated that excess salinity is stressful to angiosperms (Pennings and Moore 2001, Noe 2002, Veldhuis et al 2019. Field and numerical studies showed that soil surface salinity in salt marshes increases with marsh surface elevation toward the landward and reaches the maxima just above the mean high sea level (MHSL), beyond which the soil surface salinity starts decreasing (Wang et al 2007, Carol et al 2019. Such a consistent cross-shore salt distribution pattern is due to that longer evaporation periods at higher elevations lead to salts rather concentrated in surface soil, while the barely flooded areas above the MHSL receive much less salt load. These results can partially explain the zonation by linking topographic elevation and the distribution of halophytes, of which the responses to salinity are quite species-dependent. In some salt marshes, however, instead of decreasing, the surface soil salinity beyond MHSL may continue rising to a high level and result in salt precipitation on the marsh surface. Consequently, a hypersaline and non-vegetative zone called salt pans (alternatively known as salt barren or salt flat, figures 1(a) and (b)) will form (Wilson et al 2009, Goudie 2013. The non-vegetative feature makes salt pans rather discernible from satellite or aerial images. Due to the inhibitory effect on plant growth, salt pans have received widely attention from many studies, with an attempt to understand their formation mechanisms. Based on field data of salt pans in a Juncus marsh of northwestern Florida, Hsieh (2004) developed a theory of salt pan formation, where salt pans tend to be found in tide-dominated salt marshes and the mean high water level dictates the position of a salt pan. Hsieh (2004) called for mechanismbased modeling studies to test the theory and unravel important insights into the formation of salt pans. Wang et al (2007) modified a salt and water balance model to simulate the effects of soil, tides, topography and climate on the salt dynamics of the Atlantic and Gulf coastal marshes. Their simulation results indicate that the mean higher high water level determines the position of the salinity plateau, while the tidal irregularity affects the width of the salinity plateau. More importantly, Wang et al (2007) found that salt pans may form once the maximum salinity reaches a threshold under the influences of hydraulic conductivity, temperature, evaporation and tidal salinity. Nonetheless, the underlying mechanism of salt pan formation remains unclear, since their model only simulates salinity variations without considering salt accumulation and precipitation in marsh soils and on marsh surfaces. More recently, Shen et al (2018) numerically investigated the salt dynamics in coastal marshes with a focus on salt pans. Their simulations highlighted the importance of the hydraulic connection between the watertable and marsh surface, which sustains evaporation and supports salt pan formation in the supratidal zone. In situations of high potential evaporation rates or low marsh soil permeability, the hydraulic connection cannot be maintained and evaporation will be disrupted, thereby impeding the formation of salt pans. These findings provide an important insight into the mechanism of salt pan development.
While previous studies have attempted to reveal the formation mechanism of salt pans, the characteristics of these hypersaline zones under different conditions remain poorly understood, and such characterization is important. Since salt marshes with well-developed salt pans are widely distributed over the world, their vegetation distributions patterns can be highly different. Therefore, characterizing salt pans would have significant implications for better understanding plant zonation. To fill the research gap, this study numerically characterized salt pans in tidal marshes by simulating and analyzing salinity distribution patterns in a 2D creek-normal section representative of natural salt marshes. Three key controlling variables were examined, including slope of marsh platform, potential evaporation rate and tidal amplitude.

Methods
This study conducted the numerical simulations using SUTRASET (Zhang et al 2014, Shen et al 2018. SUTRASET was modified based on the variablesaturation, variable-density groundwater model SUTRA (Voss and Provost 2008), by incorporating an evaporation module to simulate the salt accumulation/precipitation process in porous medium. More details of the SUTRASET model are provided by Zhang et al (2014) and Shen et al (2018).

Model setup
This study considered a 2D soil section perpendicular to a creek in a tide-dominated salt marsh (figure 1(c)). By neglecting inland freshwater input, the inland boundary AB was set to be no-flow. Such a boundary condition has been used in several numerical studies on the subsurface hydrodynamics in salt marshes (e.g. Wilson and Gardner 2006, Xin et al 2010, and can be valid when assuming dry seasons with little inland freshwater input. The marsh bottom AF and the hydraulic divide EF were also set as no-flow boundaries. Three types of boundary conditions were applied to the marsh platform BC and tidal creek CDE, according to the tides and local porewater saturation conditions: (1) boundary nodes underneath the tides were assigned a hydrostatic pressure calculated based on the tidal level; (2) boundary nodes above the tides were treated as a seepage node with an atmospheric pressure if they were saturated in the previous time step (Wilson and Gardner 2006); and (3) evaporation (e.g. a specified-flux boundary condition) was applied to the exposed unsaturated section (e.g. boundary nodes above the seepage face). Evaporation of seawater is neglected in the model as it affects little the seawater salinity. Readers may refer to Shen et al (2018) for more information of the boundary conditions. being the tidal period and was set to 12 h. The initial hydraulic heads in the entire model domain were set as hydrostatic based on the MSL. Also, the initial porewater salinity in the whole domain was set to that of seawater salinity, namely 35 ppt (parts per thousand), a value commonly used in numerical studies on subsurface flow in salt marshes (Xin et al 2009(Xin et al , 2010. Hsieh (2004) identified that the physical factors involving surface topography, tidal regime and local climate are critical to the formation of salt pans in salt marshes. Therefore, this study examined three major controlling variables: tidal amplitude (A), potential evaporation rate (EVP), and slope of marsh platform (S p ). Tidal amplitudes were varied between 0.7 and 1.6 m with an interval of 0.3 m, and these values were within the range adopted by other relevant numerical studies (e.g. Wilson and Morris 2012, Wilson et al 2015, Xin et al 2017. Following Zhang et al (2014), to cover different weather conditions, the potential evaporation rates were set to 2 mm d −1 , 4 mm d −1 , 6 mm d −1 and 8 mm d −1 . Three marsh platform slopes were chosen (0.01, 0.02 and 0.03) to represent mildly, moderately and highly sloping marsh surfaces, respectively (Chapman 1974). All these values led to a total of 48 simulation cases, allowing us to characterize evaporation-induced salt dynamics in salt marshes under a wide range of scenarios (e.g. different inundation periods, different potential evaporation rates). The marsh top soil was assumed to be homogeneous and sandy loam, a commonly found soil type in coastal marshes due to tides that brought in varying ranges of sediments and vegetation that accumulate them around the root zones (Chapman 1974). With stronger water retaining capacity than sands which are critical for evapoconcentration, and higher permeability than clays which are important to replenish evaporated porewater, sandy loam has a unique characteristic that promotes the formation of salt pans. While early studies revealed significant impact of hydraulic conductivity on evaporation and thus the possible formation of salt pans (Wang et al 2007, Xin et al 2017, this study only considered sandy loam without varying the hydraulic conductivity of the marsh soil. This is because the current study aimed to characterize well-developed salt pans rather than identify the critical condition of salt pan formation.

Simulation cases and parameter ranges
The same mesh discretization scheme with 36 491 nodes and 36 000 elements and a 60 s time step size were used in all the simulations. To avoid numerical oscillations, mesh discretization must satisfy the stability criterion of grid Péclet number P e ≈ ∆L/α L ⩽ 4, with ∆L being the transport distance [L] between two sides of an element along subsurface flow direction and α L being the longitudinal dispersivity [L] (Hughes and Sanford 2004). The maximum P e under such a discretization scheme was 1.0, thereby meeting the stability criterion. Tests on the dependence of results on mesh resolution were conducted to ensure the numerical convergence of results presented in this paper. We compared the salinity distributions (in the marsh soils and on the marsh surface), hydraulic heads and boundary fluxes for four mesh resolutions: 18 361 nodes, 28 431 nodes, 36 491 nodes and 48 521 nodes. Results under the two highest resolutions were almost identical, so the resolution of 36 491 nodes was used.
The model was allowed to simulate salt dynamics driven by tide and evaporation over 5 years. Tidedriven groundwater circulations transfer low-salinity pore water from the intertidal zone to the supratidal zone, where evaporation may sustain. Readers may refer to Shen et al (2018) for details of the groundwater circulation patterns between the intertidal zone and the supratidal zone. As evaporation takes away freshwater only, salt would accumulate. When the pore water salinity exceeds the solubility, salt precipitation in the supratidal zone would occur, forming salt pans. In this study, the simulations by the elapsed time of 5 years did not reach a quasi-steady state (i.e. the total amount of salt in the marsh soil does not change over time) because evaporationinduced salt accumulation would continue, but our tests showed rather slow variations of surface salinity and thickness of precipitated salt if the simulation kept running. Since this study primarily focused on the evaporation-induced salt accumulation on the marsh surface, the results of the elapsed time were used for analysis. Note that the simulation period of 5 years is exclusive to the model setup (e.g. domain dimensions, evaporation rate, and tidal amplitude) in this study, and it may vary to reach the state of slow changes of surface salinity and salt precipitation thickness, if the model setup is highly different. Figure 2 shows the time-dependent variations of evaporation rate, saturation, salt concentration and thickness on the soil surface, and salinity distributions in the marsh soil. The liquid water saturation does not change over time, indicating that the flow condition reaches a steady condition in less than a month. On the marsh platform starting from the creek bank, the surface saturation increases with elevation and reaches a maximum value at the middle of the intertidal zone. Such rising saturation pattern was formed by porewater circulation, where saltwater infiltrates the unsaturated marsh soil from the marsh platform, and discharges from the creek bank. The porewater circulation weakens as the distance to the creek increases, resulting in higher saturation. Further landward to the saturation maxima, the surface saturation declines with elevation due to less inundation period and hence limited recharge. Within the intertidal zone, salinity remains near 35 ppt and salt precipitation is absent, as the tides can flush the salts accumulated by intermittent evaporation. Due to abundant moisture and no salt stress, the evaporation rate in the intertidal zone is maintained close to the potential evaporation rate. Despite never being inundated, the supratidal zone between 60 m < x < 130 m maintains a surface saturation above 0.5, as the local water table is rather close to the surface. The high surface moisture availability leads to persistent evaporation, which takes water away with salt left behind, resulting in the formation of a salt pan on the surface. In turn, salt accumulation reduces the evaporation rate evidently in the supratidal zone. Such an inhibitory effect of salt accumulation on evaporation rate has been previously reported by America et al (2020). At the very early stage (e.g. 0.07 year, figure 2(a)), salinity near the marsh surface gradually increases but does not reach a high level. As evaporation continues, salts underneath the marsh surface become more concentrated (e.g. 1.03 years, figure 2(b)), creating significant density gradients vertically. Consequently, the salt plume sinks through the unsaturated zone (above the watertable represented by the white dashed line), penetrating and forming salt fingering below the watertable. These salt fingers merge as they move toward the marsh bottom, finally forming a sand-clock shaped plume that spreads laterally (figures 2(c)-(f)). Further inland to the salt pan (where x < 60 m), surface saturation reduced to nearly zero, due to the disruption of the hydraulic connection between the surface and water table disrupts. As a consequence, evaporation and excess salts are absent in that region.

Results and discussions
Salt fingering in coastal marshes have been reported by several studies. In their investigations on the drivers of groundwater flow at a marsh-back barrier island transect, Ledous et al (2013) have pointed out the potential occurrence of salt fingering induced by high salinity, stating that 'high temperatures can cause a high amount of evaporation across the marsh and salinities in saltpans can reach 60 ppt during parts of the year, potentially inducing salt fingering' . Shen et al (2015) numerically discovered that, when the salinity of surface water is higher than that of subsurface water, in certain situations, unstable fingering would occur in salt marshes, greatly modifying the porewater regime, particularly for the marsh interior near the inland boundary. In addition, the study of Xin et al (2017) regarding the impact of evaporation, tidal fluctuations and rainfall on marsh soil conditions also reported the formation of salt fingering as salts become very concentrated on the marsh surface. However, both Shen et al (2015) and Xin et al (2017) did not consider the process of evaporation-induced salt precipitation on the marsh surface. This study demonstrates how salt accumulation may result in salt fingering under the effect of evaporation, thereby advocating the statement of Ledous et al (2013).
Figures 3(a1)-(d3) compare the distributions of marsh surface salinity estimated by different cases. Clearly, surface salinity variations in all these cases exhibit a similar spatial pattern: a slight increase with the surface elevation in the intertidal zone, and a sharp increase beyond the MHSL (vertical gray dashed line) to reach the plateau before falling back to 35 ppt in the supratidal zone. Generally, with the tidal amplitude and marsh platform slope unchanged, a higher potential evaporation rate leads to higher surface salinity (i.e. figure 3(a1)), since more freshwater is evaporated from the salty porewater with salt left behind. Meanwhile, the rise of overall surface salinity is more prominent when the potential evaporation rate increases from a relatively small (e.g. 2 mm d −1 ) to a moderate value (e.g. 4 mm d −1 ). Therefore, for salt marshes located in temperate area, the phenomenon of plant zonation may be more prominent, due to the greater possibility of salt pan formation. Moreover, under the same tidal amplitude and marsh platform slope, the salinity plateau is narrower when the potential evaporation rate is higher, with the seaward end remaining at the MHSL while the landward end shifting further seaward. This is because, as the topographic elevation increases landward ( figure 1(c)), the distance between marsh surface and watertable generally increases accordingly. In this regard, at higher elevations, the hydraulic connection between the evaporating surface and underlying watertable, critical for sustaining evaporation, as pointed out by Shen et al (2018), is disrupted more easily by a high potential evaporation rate. Correspondingly, the narrower salt pans (when potential evaporation rate is relatively high) would exert less influences on marsh plant distributions, in comparison to wider salt pans (when potential evaporation is relatively low).
The comparison in figures 3(a1)-(d3) further demonstrates that, the salinity plateau becomes narrower when the marsh surface is steeper (e.g. figures 3(a1)-(a3)). Such a trend is also linked to the above-mentioned hydraulic connection, which dictates the location of the right end of the salinity plateau. Therefore, the steeper the marsh surface is, the closer the location of such a critical distance is to the boundary between the supratidal zone and the intertidal zone. Another notable trend from figure 3 is that, under the same marsh platform slope and potential evaporation rate, the increase of tidal amplitude narrows the salinity plateau. This trend can be manifested in the comparison between figures 3(a1), (b1), (c1) and (d1). When tidal amplitude increases, the boundary between supratidal and intertidal zones shifts landward, but the tidal signal attenuates more quickly as it transports landward. As a result, the overall distance between marsh surface and watertable in the supratidal zone becomes greater, resulting in weaker hydraulic connections that maintains evaporation and subsequently narrower salt pans. Moreover, when the tidal forcing is stronger, the location of the salinity plateau shifts landward (e.g. figure 3(a1) vs figure 3(b1)). According to Hsieh (2004), in event of sea level rise, the upper boundary of salt pan tends to invade into the old high marsh while the corresponding lower boundary may invade into the old salt pan, which are consistent with our results, as the effect of increasing tidal amplitude (with MSL unchanged) on modifying the width and location of salinity plateau is equivalent to that of sea level rise (with tidal amplitude unchanged).
Apart from surface salinity distributions, the present study further compared the maximum surface salinity in different cases. Figures 3(e1)-(e3) show that, for given potential evaporation rate and marsh surface slope, the highest salinity decreases monotonously with the increase of tidal amplitude. Moreover, the maximum salinity increases by a relatively great extent as the potential evaporation rate rises from 2 mm d −1 to 4 mm d −1 . However, the increase of maximum salinity is much less sensitive when the potential evaporation rate further increases (e.g. figure 3(e1)). Among all the cases, the maximum salinity is the greatest when tidal amplitude is 0.7 m, potential evaporation rate is 8 mm d −1 , marsh platform slope is 0.01 ( figure 3(e1)). This situation possibly fosters an optimal hydraulic connection between the marsh surface and the watertable for sustaining evaporation.
Based on the exponential ratio between solid and solute salts (used in the SUTRASET model) and simulated surface salinities, this study calculated and plotted the distribution of precipitated salt thickness along the marsh surface (figures 4(a1)-(d3)). The section where the solid salt thickness is nonzero represents the salt pan. Apparently, the landward and seaward boundaries of salt pans formed in all the cases are within the corresponding salinity plateaus shown in figure 3. This is consistent with the findings of Hsieh (2004) that salt pans are confined to the salinity plateau. It is clear from figures 4(a1)-(d3) that, in all the simulated cases, the thickness of precipitated salt in the supratidal zone increases first and then decreases toward the inland boundary. Another common feature is that, as the potential evaporation rate increases (tidal amplitude and slope of marsh platform remain the same), the salt precipitation on the marsh surface will be thicker (e.g. figure 4(a1)). The comparison among figures 4(a1)-(d3) indicates that, when a salt pan is formed in the supratidal zone of a salt marsh, it tends to be thicker when the potential evaporation rate is high (e.g. EVP = 8 mm d −1 ), the tidal amplitude is small (e.g. A = 0.7 m), and the marsh surface is mildly sloping (e.g. S p = 0.01).
Using the results of precipitated salt thickness, we calculated the salt pan widths in different cases, based on the thickness criterion of 1 × 10-4 m. Although the use of another criterion may somewhat change the calculated salt pan width, the trend would keep the same. It is evident from figures 4(e1)-(e3) that, when tidal amplitude increases while the other two controlling variables stay the same, the variation of salt pan width almost exhibits a linear decreasing trend. More importantly, in accordance with figure 3, despite the larger surface salinity, the salt pan is narrower when the potential evaporation rate is higher (e.g. figure 4(e1)). Also, clearly, the width of salt pan also greatly depends on the slope of marsh platform. When the marsh surface becomes steeper while the other factors (e.g. potential evaporation rate, tidal amplitude) remain unchanged, the salt pan will be much narrower (e.g. figure 4(e1) vs figure 4(e2)). From figures 4(e1)-(e3), it can be inferred that large  areas of salt pans are more likely to be found in microtidal marshes that are mildly sloping and locating at temperate climate zones.

Concluding remarks
This study investigated the characteristics of salt pans in coastal marshes based on a 2D creek-perpendicular cross section. Major findings from this study are: (1) Under the continuous effect of evaporation, salinity near the marsh surface of the supratidal zone keeps rising, eventually leading to salt precipitation and the subsequent formation of a salt pan. Meanwhile, the accumulated salt in the marsh sediments create an upward salinity gradient that triggers unstable flow.
(2) The overall salt mass and the width of salinity plateau tend to decrease with the decline of potential evaporation rate, and the increase of either tidal amplitude or marsh platform slope.
(3) Consistent with the conclusions of previous studies, the locations of salt pans are within the salinity plateaus. Overall, the thickness of precipitated salt on the marsh surface will decrease when the potential evaporation rate declines, the tidal amplitude increases, or the marsh surface becomes steeper. Moreover, the salt pan is wider when the potential evaporation rate is lower, the tidal amplitude is smaller, and the marsh surface is mildly sloping.
Despite the above findings, this study used a simple semi-diurnal solar tidal signal, whereas the multi-constitute tides (e.g. spring-neap tides) may complicate the salt dynamics in coastal marshes. Moreover, our study considered a 2D normal-creek section and neglected the 3D topography, which has been found to significantly affect the porewater flow and salinity distribution patterns (Moffett et al 2010b, Xin et al 2011. In addition, the present study set a constant potential evaporation rate in each case, while in reality it may vary remarkably even during a day. These factors should be examined in future research to better characterize salt pans in coastal marshes. Despite these limitations, this study still provides a deep insight into the features of hypersaline zones in salt marshes, thereby promoting a better understanding of the eco-functions of these wetlands. A systematic and quantitative examination of the salt pan characteristics in coastal marshes would not only advance our understanding of the marsh ecosystems but also provide guidance for developing effective measures to prevent them from degradation, causing unexpected release of blue carbon.

Data availability statement
The data that support the findings of this study are available upon reasonable request from the authors.