Nitrate fate in coastal unconfined aquifers influenced by preferential flows

This study examined the influence of preferential flow on pore water flows and marine nitrogen transport reaction in variable saturation and variable density coastal aquifers. The 2-D unconfined aquifer model established was based on the software COMSOL by coupling the dynamic and chemical processes together. The results showed that preferential flow affects groundwater flow and salinity distribution, leading to a more complicated mixing process. The preferential flow resulted in an increase in mixing zone area and the upper saline plume area of 10.33 and 2.62 m², respectively, a decrease in saltwater wedge area of 7.22 m², and an increase in nitrate (NO₃^-) removal efficiency from 7.9% to 8.97%. The NO₃^- removal efficiency increases progressively with the depth (h) and quantity (n) of preferential flows; however, it decreases after a certain quantity. Further quantitative analysis revealed an increase in the intensity of nitrification and dissolved oxygen inflow flux with preferential flow depth and quantity increase. This phenomenon usually occurs on coasts where biological caves are abundant. The results also offer significant implications for designing engineering measures to mitigate saltwater intrusion and are significant to prevent groundwater quality deterioration in coastal zones.

1 Introduction

Groundwater is an important freshwater resource for industrial and agricultural purposes in coastal areas (Lu and Werner, 2013; Zhang et al., 2019). Nitrate pollution has deteriorated the quality of groundwater around the world due to an increase in human population and urban development (Anwar et al., 2014; Sun et al., 2021; Gao et al., 2024). Saltwater intrusion (SWI), as a global problem, seriously threatens freshwater resources and coastal productivity (Chang et al., 2019, 2020; Zheng et al., 2020, 2021, 2022; Chang et al., 2023, 2024). The NO₃^- concentration in drinking water derived from groundwater exceeds the standard permissible limit of 50 mg L^-1 set by the World Health Organization (WHO) in most parts of the world (Kringel et al., 2016; Radfard et al., 2018). In China, approximately 34.1% of groundwater resources are contaminated with NO₃^- (Zhang et al., 1996; Lu et al., 2019).

Groundwater flow and solute transport in aquifers under the tidal influence are very complex (Figure 1) (Robinson et al., 2007a; Fang et al., 2022a, b, 2023). Periodic tidal fluctuations can lead to the formation of two coexisting saltwater plumes in aquifers: (i) a saltwater wedge (SW) due to density-driven saltwater recirculation and (ii) an upper saline plume (USP) due to tide-driven saltwater recirculation (Robinson et al., 2007b; Kuan et al., 2019). Land-derived groundwater is generally transported between the SW and USP and discharged into the ocean near the low tide mark, which is an area where saltwater and freshwater mix between two different saltwater plumes. Several surveys have also reported a large number of macro-pores (such as crab burrows and invertebrate nests) in aquifers (Xin et al., 2009; Guimond et al., 2020; Xu et al., 2021; Pan et al., 2022; Xiao et al., 2022) and found that macro-pores are mainly concentrated in the upper supratidal zone and intertidal zone (Fanjul et al., 2008). The presence of crab holes and macro-pores increases the overall surface infiltration rate to a range of 0.1 to 1 m d^-1, which is one to two orders of magnitude higher than the matrix hydraulic conductivity (Hughes, 1998). They are well recognized to act as preferential flow paths and enhance the infiltration rate of surface water (Xiao et al., 2019), impact the distribution of salinity (Williams et al., 2014; Edith et al., 2015; Enrique et al., 2018), and change solute transport by increasing the connectivity of otherwise impermeable muddy soils (Xiao et al., 2019; Guimond et al., 2020), thus affecting biogeochemical processes in salt marshes (Xu et al., 2021; Pan et al., 2022; Xiao et al., 2022). Xin et al. (2009) studied a 3-D model simulated in the marsh with a two-layer soil configuration at the Chongming Dongtan wetland. Their results suggested that preferential flow can act as drains for the surrounding soil during the falling tide, which increases water exchange between marsh soils and tidal creeks. Xiao et al. (2019) constructed a model, and the results of their analysis suggested that preferential flow can promote soil permeability and facilitate solute transport in salt marshes. Additionally, preferential flow increases the connectivity between the soil surface and groundwater (Xu et al., 2021). It also intensifies the heterogeneity of the aquifer, leading to greater hydraulic conductivity (Kreyns et al., 2020). Although Gao et al. (2023) demonstrated that large pores increase USP and reduce SW, the mechanism is unclear. All combined, these studies confirmed that preferential flow or tidal-induced complex pore water flow and preferential flow have significant influences on solute transport. However, the effects of preferential flow distributed in the tidal zone on groundwater flow and solute transport in coastal unconfined aquifers under tidal conditions are not clear.

Figure 1

Figure 1 (A) Conceptual diagram of the seawater level, groundwater flow, and salinity distribution in a near-shore aquifer under tide action. Still water level (SWL), high tide level (HTL), low tide level (LTL), the black area is the preferential flow, the blue arrows indicate the direction of preferential flow. The upper saline plume (USP) formed by tide-driven recirculation (TDR), the saltwater wedge formed by density-driven recirculation (DDR), and the saltwater wedge toe (TOE) are shown. The freshwater–saltwater mixing zone associated with the saline plums are also depicted. (B) Numerical model domain, parameters, and boundary condition, the dashed box in (B) illustrates the area depicted in (A).

Saltwater transports a large quantity of dissolved oxygen (DO), NO₃^-, ammonium ion (NH₄⁺), and dissolved organic carbon (DOC) into the aquifer. Several chemical reactions such as DOC degradation, aerobic respiration, nitrification (NH₄⁺ transformation to NO₃^-), and denitrification (NO₃^- transformation to N₂) can take place (Meile et al., 2009; Sun et al., 2019). These processes are affected by tidal force, which can change the recycling rate and chemical composition of saltwater (Rocha, 2013; Zheng et al., 2020; Gao et al., 2022). Shuai et al. (2017) established a 2-D estuarine subsurface flow model based on the numerical model, assuming that DO, NO₃^-, NH₄⁺, and DOC primarily originate from rivers. Nitrification occurs in the shallow layers, while denitrification occurs in the deeper anaerobic layers. These reactions, which are both controlled by redox conditions and DOC concentrations, play an important role in determining NO₃^- in the aquifer. Preferential flow can affect nitrification and denitrification reactions by directly or indirectly changing soil properties (Pan et al., 2022). These reactions are also influenced by the redox conditions, aeration, and particle composition of the solution (Xin et al., 2009). Xiao et al. (2019) conducted field observations and flow modeling to assess how crab burrows drive carbon exchange in an intertidal marsh in South Carolina. They found that the concentrations of dissolved inorganic and organic carbon in crab burrow pore water differ from that in the surrounding soil matrix, and the gas-phase concentrations of CO₂ in crab burrows were approximately six times greater than that in ambient air. These results suggested that crab burrows increase the reaction of carbon. Heiss (2020) evaluated the effect of deep whale burial on NO₃^- removal efficiencies in aquifers and found that organic carbon sources affected NO₃^- removal since DOC not only provided the reactant for denitrification but also provided an anaerobic environment (DOC oxidation). Cheng et al. (2020) investigated the effect of preferential flows on nitrification and denitrification in aquifers, and the results showed that preferential flows promoted sediment nitrification, and the degree of promotion was significant with the increase of density. Li et al. (2021) conducted laboratory experiments and found that crab burrows significantly affected the concentrations of NH₄⁺ and NO₃^- and significantly promoted the nitrogen exchange flux at the sediment–water interface. It can be seen that the preferential flow not only effectively changes the pore water flow and solute transport rate but also increases the oxygen content of the aquifer and the contact area between the sediment and O₂, thus interfering with the reaction of NO₃^-.

So far, the effects of preferential flow distributed in the tidal zone under tidal conditions on pore water flow and nitrogen–sea source transport reaction in coastal unconfined aquifers remain unclear. In this study, numerical simulations were conducted to address the effects of preferential flow number and depth on the pore water flow, nitrification, and denitrification in sea sources under tidal action. The results might provide a theoretical reference for designing engineering measures to mitigate saltwater intrusion, groundwater quality management, and the ecological restoration of coastal zones.

2 Methods

The groundwater flow and solute transport reaction models are used to study the influence of preferential flow on salinity change, pore water flow, and nitrogen transport reaction in aquifers.

2.1 Groundwater flow and solute transport model

A 2-D coastal wetland model was developed to study the salinity change and solute transport of aquifers with variable saturation and density using the COMSOL finite element software. The simulation of flow with variable saturation and density was based on Richard’s equation (Equation 1). Solute transport was based on the advection–dispersion–reaction equation (Equation 2) in porous media (Shuai et al., 2017):

$\begin{array}{ll}\rho (S_{e}S+\frac{w_m}{\rho g})\frac{\partial p}{\partial t}+\nabla ·\rho \left(\frac{-k_s}{\mu }{k}_r(\nabla P+\rho g\nabla Z)\right)=Q_m& (1)\end{array}$

$\begin{array}{ll}\frac{\partial (\theta C_i)}{\partial t}+\nabla ·\rho (C_{i}u)=\nabla ·(D\nabla C_i)+R_i& (2)\end{array}$

where ρ is fluid density (kg m^-3), C_i is the concentration of species i (mM), P is the pressure (Pa), and Z is the elevation head (m). The value of k_s is the saturated hydraulic permeability (Equation 3), while k_r stands for the relative permeability (m s^-1) (Equation 4), u is the Darcy velocity (m s^-1), θ is the water content (-) (Equation 6), g is the acceleration of gravity (9.81m s^-2), Q_m in this equation represents a stress source term (kg m^-3 s^-1), w_m is the specific moisture capacity (m^-1) (Equation 5), D is the hydrodynamic dispersion coefficient (m² s^-1), R_i is the reaction rate for species i (mM s^-1), S_e is the effective saturation (-) (Equations 7, 8), and S is water storage (P_a^-1). The constitutive relationship between relative permeability and pore pressure obeys the van Genuchten model (van Genuchten, 1980):

$\begin{array}{ll}{k}_s=\frac{\mu K}{\rho g}& (3)\end{array}$

$\begin{array}{ll}{k}_r=S_e^{0.5}{(1-{[1-S_e^{1/m}]}^m)}^2& (4)\end{array}$

$\begin{array}{ll}{w}_m=\frac{\alpha m}{1-m}({\theta }_s-{\theta }_r)S_e^{1/m}{[1-S_e^{1/m}]}^m& (5)\end{array}$

$\begin{array}{ll}\theta ={\theta }_r+S_e({\theta }_s-{\theta }_r)& (6)\end{array}$

$\begin{array}{ll}{S}_e=\frac{1}{{(1+{\left|\alpha P\right|}^{n_i})}^m}& (7)\end{array}$

$\begin{array}{ll}m=1-\frac{1}{n_i}& (8)\end{array}$

where θ_s and θ_r is the saturated water content and relative water content (-), respectively. K is the hydraulic conductivity (m s^-1), where n_i, which is related to m in Equation 8, and a is the fitting parameter which describes the shape of both the moisture (m^-1) and relative permeability functions (-) obtained by van Genuchten (van Genuchten, 1980).

The model is a 2-D beach aquifer cross-section (Figure 1), where the left and bottom boundary (BAE) and the upper boundary (CD) are both set as no flow and zero solute flux boundaries (Xin et al., 2010; Shen et al., 2018). The right boundary (DE) is set as the inland boundary, while the ocean boundary (BC) is the permeable layer that is used to represent the interface between the beach surface and seawater or atmosphere and assigned a semi-pervious layer or seepage face boundary. The process is realized by the following formula (Equations 9–11):

$\begin{array}{ll}-n·\rho \frac{-k_s}{\mu }{k}_r(\nabla P+\rho g\nabla Z)=\rho R_b(H-H_b)& (9)\end{array}$

$\begin{array}{ll}{R}_b=\frac{k_s}{L}& (10)\end{array}$

$\begin{array}{ll}\begin{array}{l}C=C_{sea}n·u<0\\ -n·D\nabla C=0n·u\ge 0\end{array}& (11)\end{array}$

where R_b is the conductance term (s^-1), defined as the ratio of the saturated hydraulic conductivity (k_s) with a coupling length scale (L, m). It was set at a high value allowing water to readily move in and out of the interface. H_b is the external head representing sea level (m), H is the total head (m), and C_sea is seawater salinity (ppt).

In order to deal with the boundary conditions reasonably, two conditions have been considered. When the aquifer is saturated and the external head (H_b) of the sea level is higher than the elevation head, the head difference is determined by the difference between the overlying seawater and the beach interface. However, if the total head (H) is higher than the sea level (H_b), and the pressure head is equal to or higher than the atmospheric pressure. The boundary is a seepage surface. If the value of R_b is 0 when the aquifer is unsaturated, the boundary is a no-flow boundary. The salinity is determined by the groundwater velocity at the node. If the velocity is inward (to the aquifer), the salinity is estimated by seawater salinity (C_sea = 35 ppt). If the groundwater flows from an aquifer, the zero concentration gradient was specified at the nodes. Previous studies have demonstrated that this modeling approach is useful in describing the water level movement and residence times within the transition zone (Xin et al., 2010; Shen et al., 2019).

A time-varying head (H(t)) Equation 12 was implemented at selected cells in this zone by:

$\begin{array}{ll}H(t)=H_{msl}+A\sin (\omega t)& (12)\end{array}$

where H(t) is tidal head (m), determined by time; H_msl is the mean sea level (10 m); A is tidal amplitude (0.5 m); $\omega =2\pi /T$ is the angular frequency; and T is the tidal cycle (semi-diurnal tide, 12 h). These parameters were used in many studies (Anwar et al., 2014; Shuai et al., 2017).

The model domain represented a homogeneous and isotropic coastal aquifer with a thickness of 12 m and a sloping beach boundary (slope of 0.1; Figure 1), which was comparable to a typical sandy coastal aquifer system (Carsel and Parrish, 1988) with a hydraulic conductivity of 15 m d^-1, longitudinal dispersion coefficient (a_L) of 0.2 m, transverse dispersion coefficient (a_T) of 0.02 m (Anwar et al., 2014), and porosity of 0.45. Moreover, seawater salt concentration and density were set to 35 ppt (mass fraction) and 1,025 kg m^-3, respectively. The freshwater salt concentration of 0 ppt and a density of 1,000 kg m^-3 (see Table 1) were both adopted from Robinson et al. (2014) and Shen et al. (2019).

Table 1

Table 1 Model parameters, boundary condition, and kinetic parameter.

The preferential flow had a width of 0.16 m and was assumed to be a homogeneous medium with high permeability and porosity (the permeability was 100 times higher than that of the matrix and the porosity was set as 1). This technology enabled simulations of the preferential flow rapid responses to tidal-water-level fluctuations without affecting the simulation results (Xin et al., 2009). Similar techniques and parameters have been used in groundwater models to simulate pore water flow and solute transport induced by macro-porosity (Xin et al., 2009; Xiao et al., 2019). The preferential flow was set upstream of the intertidal zone, and this was adopted from other studies where it was reported that preferential flow is mainly concentrated in the upper supratidal zone and intertidal zone (Fanjul et al., 2008; Ying, 2021).

For all simulations, the mesh density increased in the preferential flow zones and at the beach. The values for the Courant number and numerical Péclet number did not exceed 1 and 4, respectively, which satisfied the stability criterion and avoided numerical oscillation. During the simulation, a constant temperature was maintained, resulting in constant dynamic viscosity.

The simulation first conducted a no-preferential flow model under the tidal action. Subsequently, to investigate the influence of depth and quantity of preferential flow, we set different depths (2.3, 3.3, and 4.3 m) and quantities (from 1 to 3) of preferential flow conditions under the tide action (five groups of experiments were conducted). The specific parameters are listed in Table 2.

Table 2

Table 2 Model parameters and TOE and average salinity.

2.2 Solute reaction model

To study the transformation of marine nutrients, a reaction network of four reactive species was applied to the reaction model. The model considered nitrification, denitrification, aerobic respiration, and DOC degradation.

The reaction model simulated the first-order reaction that required mixing for the inflow of DO, NO₃^-, NH₄⁺, and DOC. The study omitted the production of NH₄⁺ from DOC degradation and nitrification and only focused on the transformation of NO₃^- from marine-derived sources. The reactions and kinetic rate expressions are shown in Table 3. Kim et al. (2017) validated this reaction network by comparing the simulation results with field data from the Cape Shores beach adjacent to Delaware Bay.

Table 3

Table 3 Reaction and kinetic rate expressions.

Based on the studies of Bardini et al. (2012) and Spiteri et al. (2008), the parameters of the reactions are shown in Table 1. This study determined the influence of preferential flow under tide influence on the chemical transformations in a nearshore aquifer. The sensitivity analyses on the specific kinetic parameter values adopted were not performed, although these parameter values vary in different coastal settings. Similar to salt transport, the nutrient concentration was set along with the aquifer–ocean interface for the reactive solutes. The boundary conditions of the nutrients are presented in Figure 1, and the specific nutrient values are shown in Table 1. The TOE (the saltwater wedge toe) was used to evaluate the length of the saltwater wedge, and the NO₃^- removal efficiency (R_N) Equation 13 was quantified by:

$\begin{array}{ll}{R}_N=\frac{{\displaystyle {\int }_0^t{\displaystyle {\int }_{\text{Ω}}{R}_{{\text{NO}}_3^{-}}\theta d\Omega dt}}}{{\displaystyle {\int }_0^t{\displaystyle {\int }_l{c}_{0{\text{NO}}_3^{-}}{f}_{b}dl}}}& (13)\end{array}$

where Ω is the domain area, t is the tidal cycle time (12 h), f_b is the boundary flux, $R_{{\text{NO}}_3^{-}}$ is the rate of NO₃^- reaction, and l is the length of the boundary layer.

We performed four steps: (i) groundwater flow and salinity distributions for all cases at the static water level were run to steady state, (ii) the hydraulic gradient and salinity distributions obtained at the static water level were used as the initial conditions under tidal action, (iii) the salinity distributions obtained in step (ii) were used as the initial conditions for the chemical solute transport and were run to steady state; (iv) and the results of the previous step were run as the initial conditions for the reaction network until they stabilized. The simulations were run for 700 days to allow the reactive solutes to reach equilibrium.

3 Results

3.1 Model calibration

To verify the reliability of the model, we developed a numerical model consistent with the lab experiment of Xie et al. (2023), with an aquifer length of 4 m and width of 0.7 m and porosity of 0.45; the parameters a and n_i were set to 5.9 m^-1 and 2.68, respectively, and the longitudinal dispersion coefficient (a_L) was 0.005 m and transverse dispersion coefficient (a_T) was 0.0005 m [for the specific parameters, please refer to Xie et al. (2023) for details]. We consider three case conditions: non-preferential flow, preferential flow at x = 0.8 m, z = 0.3 m, and preferential flow at x = 1.8 m, z = 0.3 m, and the preferential flow scales and parameter settings are also consistent with those of Xie et al. (2023).

From the experimental results in Xie et al. (2023) and numerical simulation results (Figure 2), the simulated results matched the experimental in all cases well (Figure 2). We can find that tidal forcing formed a USP and SW in the intertidal zone where saltwater infiltrated near the high tide mark and was discharged near the low tide mark. A fresh groundwater discharge zone (FDZ) separated the circulation cell from the lower interface located at the base of the beach. As expected, under the impact of preferential flow, when the preferential flow location is at x = 0.8 m, it has some effect on the salinity distribution and reduces the saltwater wedge deed extension. When the preferential flow location was at x =1.8 m, it has little effect on the salinity distribution, which is consistent with the salinity distribution pattern of the no-preferential flow condition.

Figure 2

Figure 2 Experimental (A1), (B1), (C1) (Xie et al., 2023) and simulated (A2), (B2), (C2) results of non-preferential flow and preferential flow. The horizontal dotted lines is tidal range; in (A1), (B1), and (C1), the solid black lines indicate the numerical results for 50% salinity contour isohalines by Xie et al. (2023), and in (A2), (B2), and (C2), the solid black lines indicate 50% salinity contour isohalines based on COMSOL. PF is preferential flow.

The experimental and numerical results consistently demonstrated that the preferential flow could affect the salinity distribution in terrestrial environments, and the 2-D numerical model could be used to predict the water–salt exchange processes in aquifers.

3.2 Salinity and pore water flow

The salinity distribution under two different conditions is presented in Figure 3. Preferential flow modified the salinity distribution under tidal influence. The average salinity (average value of the whole aquifer) in the aquifer decreased from 7.21 to 7.03 ppt, while the TOE decreased from 24.2 to 23.09 m. Under the tidal condition (Figure 3), the area of the mixing zone (MZA) (the area between 10% and 90% salinity contour) was 53.87 m², the area of USP (USPA) was 10.37 m², and the area of SW (SWA) was 156.0 m². However, under the preferential flow, the mixing zone and USP were increased, and SW was decreased, i.e., the MZA is 64.20 m², the USPA is 12.98 m², and the SWA is 148.79 m². This indicated that the preferential flow enhanced freshwater–saltwater mixing and USP (Figure 3E); the preferential flow weakened the length of SW in the aquifer and reduced the salinity.

Figure 3

Figure 3 Simulation results: salinity distribution of tidal action in (A, B), the pink, black, and white solid lines show the 10%, 50%, and 90% salinity contour, respectively. Saltwater wedge (SW), upper salinity plume (USP), saltwater wedge toe (TOE), and fresh water discharge area (FDZ). TOE is the intersection of 50% salinity contour and x-axis. The horizontal black dotted line shows the average sea level and tidal range. In (C, D), the white solid line and data are particle path and travel time (unit: day). The particles start from the inland boundary and the beach, respectively. The gray vertical line indicates the preferential flow (h = 3.3 m). (E) Salinity difference graph (preferential flow - no preferential flow).

To determine the effect of preferential flow on groundwater flow, saltwater circulation, and saltwater–freshwater exchange under tidal conditions, the travel time for specified particles in the aquifer was calculated [readers may refer to Xin et al. (2010) for more details]. Three particles were released both near the inland boundary and the beach surface (refer to Figure 3, Table 4 for the specific release location).

Table 4

Table 4 Particle travel times and starting position.

Under tidal action (Figure 3), a USP developed in the intertidal zone, with particles released from inland along the moving path of freshwater bypassing the USP and SW discharge to the sea. The deeper particles took longer than those in the shallow, e.g., the one starting from the shallow (x = 75 m, z = 9 m) took 257.1 d, while that from the deep inland (x = 75 m, z = 3 m) took 267.8 d (Figure 3), increasing by 4.2%. The difference in travel time among the three particles at the seaside varied dramatically as the particles released to the deeper beach moved further landward and took longer paths, e.g., the particle starting from x = 13 m, z = 8.8 m took 176.7 d, while that from x = 3  m, z = 7.42 m took 1,026.9 d (Figure 3, Table 4).

Preferential flow leads to pore water flow and an increase in USP. Therefore, the particles released from the inland shallow layer experienced USP resistance across the preferential flow, with the travel time increased slightly, e.g., it took 257.1 d for particles released from x = 75 m and z = 9 m and 258.2 d for preferential flow (Figure 3). The travel time of particles released from the beach was consistent with the travel time under non-preferential flow, decreasing with an increase in USP, e.g., the particle released from x = 3 m and z = 7.42 m needed 1,026.9 d with the non-preferential flow but only 908.4 d in the presence of preferential flow.

To further discuss the influence of preferential flow on groundwater flow and circulation, Figure 4 presents the velocity distribution of rising tide, high tide, falling tide, and low tide without and with preferential flow. Only the velocity distribution in the area of preferential flow is shown for comparison. As shown in Figure 4, the velocity distribution differed with the seawater level. When the tide rose, salt water flowed into the aquifer through the intersection of the tidal level and the aquifer. When the tide fell, the tidal water retreated, and the pore water flowed to the ocean (Figure 4). Under the preferential flow condition at the rising tide, a relatively large pore water velocity appeared and rapidly flowed into the aquifer, thereby providing a preferential gateway for the pore water. When the flood completely flowed over the marsh platform (t = 3 h), the saltwater level reached the maximum (10.5 m), and a sharp rise in the pore water velocity near the preferential flow model was observed [Figure 4—(2B)]. The direction of the pore water flow under preferential flow was vertically downward, and the magnitude increased considerably. This indicated that the over-topping water rapidly enters preferential flow and diffuses into the surrounding marsh soil. Furthermore, preferential flow enhances the mixing and dilution of saltwater–freshwater, thus reducing the salinity of groundwater and increasing the quality of drinking water.

Figure 4

Figure 4 Simulation results: flow velocity at rising tide (t = 1 h), high tide (t = 3 h), falling tide (t = 6 h), and low tide (t = 9 h), (1A) - (4A), (1B) - (4B) under without and with preferential flow, respectively. (Note: for comparison, only the velocity distribution in the area where the preferential flow).

3.3 Distribution of nutrients

The distributions of DO, NO₃^-, NH₄⁺, and DOC under tidal [Figure 5—(1A), (2A), (3A), and (4A)] and preferential flow conditions [Figure 5—(1B), (2B), (3B), and (4B)] are presented in this section. There were differences in the distribution of nutrients. DO was mainly concentrated in the shallow layers of the aquifer and NO₃^- was in the shallow and medium layers of the aquifer (Note: we used shallow, middle, and deep zones to differentiate the characteristics of the distribution area of DO, NO₃^-, NH₄⁺, and DOC. In this study, these three intervals are a relative region varying with the preferential flow, so we are only distinguishing between the same model) [Figure 5—(1A), (1B), (2A), and (2B)], while NH₄⁺ [Figure 5—(3A) and (3B)] and DOC [Figure 5—(4A) and (4B)] were distributed in the entire aquifer with a relatively high concentration in the shallow layer and in the middle and shallow layers, respectively, e.g., the average concentration of DOC in the region z > 4 m is 0.058 mM larger than that at z< 4 m (0.025 mM), and the average concentration of NH₄⁺ in the region z > 4 m is 0.012 mM larger than that at z< 4 m (0.011 mM). The nutrient concentration increased along preferential flow in the reaction process [Figure 5—(1B), (2B), (3B), and (4B)].

Figure 5

Figure 5 Simulation results: nutrient concentration distribution without preferential flow action (A) and with preferential flow action (B) (1), (2), (3), and (4) is DO, NO₃^-, NH₄⁺, and DOC distribution after reaction, respectively. The black solid line is the isoline of denitrification, and the gray vertical line indicates the preferential flow; note that only the intertidal zone is shown.

3.4 Reaction rate distribution

The reaction rate distributions of nitrification [Figure 6—(1A) and (1B)], denitrification [Figure 6—(2A) and (2B)], respiration [Figure 6—(3A) and (3B)], and DOC degradation [Figure 6—(4A) and (4B)] are described in this section.

Figure 6

Figure 6 Simulation results: reaction distribution (mM s^-1). (1C)–(4C) is nitrification, denitrification, aerobic respiration, and DOC degradation reaction zone comparison under tidal and preferential flow. The gray vertical line indicates the preferential flow, PF is preferential flow.

Nitrification and respiration occurred primarily in the surface layers of the aquifer, while denitrification was dominant under the influence of anaerobic bacteria in the middle and deep layers where oxygen is depleted (Figure 6). DOC degradation occurred in the USP and SW of the whole aquifer under tidal action and preferential flow. However, the reaction rate in the middle and upper layers was relatively high [Figure 6—(4A) and (4B)], e.g., the average rate of DOC degradation in the region z > 4 m is 1.49 × 10^-5 mM d^-1 larger than that at z< 4 m (6.54 × 10^-6 mM d^-1). Compared to the non-preferential flow condition, preferential flow resulted in the extension of nitrification, respiration, and DOC degradation along with the preferential flow in the reaction process. However, denitrification consistently took place in the deeper USP.

We quantitatively analyzed the NO₃^- removal efficiency (R_N). Under non-preferential flow condition, the R_N is 7.9%. When we insert the preferential flow, the NO₃^- removal efficiency (R_N) was modified, e.g., under the preferential flow action, the R_N is 8.97%. These results suggest that preferential flow increases the NO₃^- removal efficiency, thus decreasing the salinity distribution.

3.5 Sensitivity analysis

3.5.1 Preferential flow depth

The mixing zone and USP increased vertically with preferential flow depth, the retreat of SW, and the decrease in salinity in the aquifer with an increase in preferential flow depth (Figure 7), e.g., when the preferential flow depth was 2.3 m, MZA was 59.04 m², USP was 12.50 m², SW was 152.50 m², salinity was 7.19 ppt, and TOE was 23.10 m, and when the preferential flow depth was increased to 4.3 m, MZA was 69.11 m², USP was 13.01 m², SW was 147.22 m², salinity was 7.00 ppt, and TOE was 22.90 m (Figure 7).

Figure 7

Figure 7 Simulation results: the mixing zone area (MZA) (A), the USP area (USPA) (B), the SW area (SWA) (C), the salinity (D), TOE (E), and the NO₃^- removal rate (R_N) (F) with a different preferential flow depth.

The travel time of particles released from deep inland increased with preferential flow depth. With the persistent vertical increase of USP, the travel time and path of the particle passing the USP became longer, e.g., the particle released from x = 75 m and z = 3 m took 273 d and 280.6 d with depth that was 2.3 and 4.3 m, respectively. The travel time of the particle released from the beach was consistent with the previous section, decreasing with the preferential flow depth, e.g., the particle released from x = 3 m and z = 7.42 m took 1,007.6 d with depth that was 2.3 m, while under the preferential flow with depth at 4.3 m, it only required 927.1 d. This might be due to the retreat of SW with the change in preferential flow depth (Figure 8, Table 4).

Figure 8

Figure 8 (A–C) Simulation results: rate of change in particle transport time. Positive values indicate an increase, while negative values indicate a decrease, and the colors and the data in parentheses are the coordinates of the particle release point.

The distribution of nutrients was consistent with the discussions in Sections 3.4 (Supplementary Figure S1). The significance of preferential flow depth on the distribution of nitrification, denitrification, respiration, and DOC degradation is shown in Supplementary Figure S2. All reactions were consistent with the discussions in Section 3.5. Nitrification and respiration consumption of deeper DO have their long residence time in the aquifer. Nitrification and respiration primarily occurred in the surface layers (Supplementary Figure S2). Denitrification was dominant in the middle and deep layers of the aquifer (Supplementary Figure S2), and DOC degradation occurred throughout the aquifer with a faster reaction rate in the middle and upper layers (Supplementary Figure S2). However, each reaction at the USP acted differently with preferential flow depth due to the rapid inflow of nutrients through the preferential flow, e.g., when the preferential flow depth was 2.3 m, the NO₃^- produced by nitrification was 3.184 g d^-1, and the NO₃^- removal efficiency was 8.72%. However, under a preferential flow depth of 4.3 m, the NO₃^- produced by nitrification was 3.298 g d^-1, and the NO₃^- removal efficiency was 9.04%.

3.5.2 Preferential flow quantity

The USP increased with preferential flow quantity, MZ increases and then decreases, SW decreases and then increases and decrease in salinity in the aquifer with an increase in preferential flow quantity (Figure 9), e.g., when the preferential flow was 1, MZA was 69.11 m², USP was 13.01 m², SW was 147.22 m², salinity was 7.00 ppt, and TOE was 22.9 m, and when the preferential flow quantity was increased to 3, MZA was 66.40 m², USP was 13.63 m², SW was 145.04 m², salinity was 6.94 ppt, and TOE was 22.50 m (Figure 9).

Figure 9

Figure 9 Simulation results: the mixing zone area (MZA) (A), the USP area (MZA) (B), the SW area (MZA) (C), the salinity (D), TOE (E), and the NO₃^- removal rate (F) with a different number (n) of preferential flow.

The travel time of particles released from deep inland increased with preferential flow quantity. With the persistent vertical increase of USP, the travel time and path of the particle passing the USP became longer, e.g., the particle released from x = 75 m and z = 3 m took 280.6 d and 281.5 d at n = 1 and n = 3, respectively. The travel time of the particle released from the beach was consistent with the previous section, decreasing with the preferential flow quantity, e.g., the particle released from x = 3 m and z = 7.42 m took 927.1 d with n that was 1, while under the preferential flow with n that was 3, it only required 872.2 d (Figure 10, Table 4).

Figure 10

Figure 10 (A, B) Simulation results: rate of change in particle transport time. Positive values indicate an increase, while negative values indicate a decrease. The colors and the data in parentheses are the coordinates of the particle release point.

The significance of preferential flow quantity on the distribution of nitrification, denitrification, respiration, and DOC degradation is shown in Supplementary Figure S3. Nitrification and respiration consumption of deeper DO have their long residence time in the aquifer. Nitrification and respiration primarily occurred in the surface layers (Supplementary Figures S3, S4), while denitrification was dominant in the middle and deep layers of the aquifer (Supplementary Figure S3), and DOC degradation occurred throughout the aquifer with a faster reaction rate in the middle and upper layers (Supplementary Figure S4). However, each reaction at the USP acted differently with preferential flow quantity due to the rapid inflow of nutrients through the preferential flow, e.g., when the preferential flow quantity was 1, the NO₃^- produced by nitrification was 3.298 g d^-1, and the NO₃^- removal efficiency was 9.03%; when the preferential flow quantity was 2, the NO₃^- produced by nitrification was 3.347 g d^-1, and the NO₃^- removal efficiency was 9.04%. However, under a preferential flow quantity that was 3, the NO₃^- produced by nitrification was 3.34 g d^-1, and the NO₃^- removal efficiency was 8.73%.

4 Discussion

4.1 Implication of preferential flow on pore water flow

A USP promotes the retreat of SW under tidal action, which renders the SW closer to the ocean and limits saltwater intrusion. Preferential flow enhanced (Figure 3) this effect, driving SW further seaward, and this impact gradually increases with the depth of preferential flows (Figures 7, 9). These results are consistent with those of Gao et al. (2023), but they consider the under spring-neap tide action. The high hydraulic conductivity of the preferential flow increased the salt transport, resulting in a larger USP than the non-preferential flow; this could control SWI to some extent. Several previous studies have been performed to investigate saltwater–freshwater mixing dynamics—for example, Xie et al. (2023) investigated the effect of fracture characteristics on salinity distribution and groundwater flow through experiments and numerical simulations. They reported that the vertical fractures had a limited impact on most SWI properties. This seems to be different from our results. This is mainly due to the effect of fracture on SWI properties depending on their relative position to the saltwater. In our study, we focused on the longitudinal rift located inside the USP and connected to it. In addition, they noted that fractures can increase the mixing zone, which is consistent with our results. Preferential flow acted as a drain for the shallow layer during the falling tide and for a recharge well during the rising tide. This was the main reason for the increase in USPs and the decrease in the salinity of the aquifer. Additionally, Xiao et al. (2019) found that the salinity in crab burrows in the intertidal zone at high tide was slightly higher than in the soil matrix, e.g., the salinity of crab burrows reached 34 PSU and larger than the nearby soil matrix (33.4 PSU). This suggested that the presence of crab burrows can greatly enhance salt transport in salt marshes.

4.2 Implication of preferential flow on nitrification and denitrification

The nutrients from the sea were initially input in the numerical model with the salinity simultaneously. Nitrification primarily occurred in the surface layer (Figure 6, Supplementary Figures S2, S4) since nitrification, respiration consumption of deeper DO, and denitrification were dominant in the middle and deep layers of the aquifer (Figure 6, Supplementary Figures S2, S4). Shuai et al. (2017) had the same report, and they assumed that DO, NO₃^-, NH₄⁺, and DOC primarily originate from rivers. They noted that nitrification occurs in the shallow layers, while denitrification occurs in the deeper anaerobic layers. In the presence of preferential flow, nitrification mainly occurred in the surface layer of the aquifer and extended along with the preferential flow (Figure 6, Supplementary Figures S2, S4). Nitrification and respiration consumption of deeper DO (Figure 5, Supplementary Figures S1, 3) and the preferential flow increase in the area of sediment air/water interface increased the exchange between surface water and groundwater and the concentration of DO in the aquifer. This resulted in the extension of nitrification and respiration along the preferential flow (Figure 6, Supplementary Figures S2, S4) and decreased the carbon stock of the sediment. Preferential flow increased the concentration of nutrients at the USP by dramatically enhancing nutrient transport and reducing SW as the USP pushed the SW seaward. However, total nitrification increased under preferential flow conditions and gradually increased with depth and amount of preferential flow; however, it decreases after a certain quantity.

In the presence of tidal action, the distribution of nutrients in the aquifer formed a USP and a SW. The travel and residence times were longer in the deeper layers of the aquifer for nutrients, nitrification, and respiration consumption of deeper DO (Figure 5, Supplementary Figures S1, S3), creating anaerobic conditions for denitrification. Subsequently, the NO₃^- produced was used up during denitrification, resulting in the reduction of NO₃^- in the deep layer (Figure 5—(2B), NO₃^- distribution). Denitrification occurred in the deeper USP as nitrification in the upper USP not only provided the reactant for denitrification but also offered an anaerobic environment. The NO₃^- removal efficiency increased compared to that in non-preferential flow, and this effect increases progressively with the depth and number of preferential flows; however, it decreases after a certain quantity. It can be observed that the existence of macro-porosity (such as crab burrows and invertebrate nests) facilitates the removal of NO₃^- from sea sources. Some previous studies demonstrated that denitrification also increased with the size of the mixing zone (Spiteri et al., 2008; Heiss et al., 2017; Gao et al., 2023). Heiss et al. (2017) reported that saltwater–freshwater mixing could promote denitrification in the intertidal zone.

In addition, Gao et al. (2023) investigated the effect of macropores on the salinity distribution and denitrification of the aquifer under spring-neap tide action. They pointed out that macropores can lead to NO₃^- removal efficiency increase. They investigated the reaction between DOC from the sea source and NO₃^- from the land source. In this study, we considered the relationship between DOC from the sea source and NO₃^-, which is the main reason for the difference, and it is also noteworthy that they also suggested that macropores can increase the mixing of saltwater and freshwater. This is consistent with our results. Thus, our study provides new insights into the relationship between macropores and salinity distribution and denitrification in the aquifers.

4.3 Knowledge gaps and research needs

This study elucidated the effect of preferential flow on salinity change and solute transformation under tidal conditions. It determined the importance of the depth and quantity of preferential flow, which guide the mitigation of saltwater intrusion and contaminant removal. However, the location and density (Fanjul et al., 2008; Ying, 2021) distribution of preferential flow should be further studied. The preferential flows are based on the actual and generalized, and the real aquifer preferential flow distribution pattern and specifications are extremely complex. Future research can focus on the real aquifer preferential flows on the saltwater intrusion characteristics as well as on the impact of solute transport transformation. Waves drive seawater recirculation (Anwar et al., 2014), and seasonal groundwater level drives groundwater fluctuations (Michael et al., 2005; Liu et al., 2016). This leads to the more complex nearshore mixing dynamics of freshwater and seawater, which, in turn, complicates NO₃^- transformation. In addition, other factors may also impact the biogeochemical processes in aquifers, including precipitation, evaporation, and long fluctuating tides (e.g., spring-neap tides), and sediment heterogeneity (Heiss and Hichael, 2014; Geng and Boufadel, 2017; Kreyns et al., 2020; Gao et al., 2023; Zheng et al., 2023), which need to be considered in the future work.

5 Conclusions

This study examined the combined effect of preferential flow and tide on pore water flow and marine nitrogen transport reaction in coastal aquifers. The following conclusions might be drawn:

(1) Under tidal action, preferential flow increases the hydraulic conductivity of the aquifer and accelerates the pore water flow and solute transport. The preferential flow results in the vertical increase of the USP and further retreat of SW. In the presence of preferential flow, the NO₃^- removal efficiency is increased.

(2) Nitrification mainly occurs in the surface layer of the aquifer, while denitrification dominates the middle and deep layers. The nitrification effect increases with the increase in the depth and quantity of preferential flow, and the NO₃^- removal efficiency increases progressively with the depth and number of preferential flows; however, it decreases after a certain quantity.

Data availability statement

The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding authors. The model calibration experimental data used to produce Figure 2, A1, B1 and C1 are compiled from Xie et al. 2023.

Author contributions

CG: Conceptualization, Investigation, Methodology, Writing – original draft, Writing – review & editing. JK: Funding acquisition, Resources, Supervision, Visualization, Writing – review & editing. JW: Writing – review & editing. WC: Supervision, Writing – review & editing.

Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This research was supported by the National Natural Science Foundation of China (51979095) and acknowledges the Postgraduate Research and Practice Innovation Program of Jiangsu Province (KYCX22_0658).

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmars.2024.1369869/full#supplementary-material

References

Anwar N., Robinson C., Barry D. (2014). Influence of tides and waves on the fate of nutrients in a nearshore aquifer: numerical simulations. Adv. Water Resour. 73, 203–213. doi: 10.1016/j.advwatres.2014.08.015

CrossRef Full Text | Google Scholar

Bardini L., Boano F., Cardenas M., Revelli R., Ridolfi L. (2012). Nutrient cycling in bed-form induced hyporheic zones. Geochim. Cosmochim. Acta 84, 47–61. doi: 10.1016/j.gca.2012.01.025

CrossRef Full Text | Google Scholar

Carsel R., Parrish R. (1988). Developing joint probability distributions of soil water retention characteristics. Water Resour. Res. 24, 755–769. doi: 10.1029/WR024i005p00755

CrossRef Full Text | Google Scholar

Chang Q., Zheng T., Chen Y., Zheng X., Walther M. (2020). Investigation of the elevation of saltwater wedge due to subsurface dams. Hydrol. Process 34, 4251–4261. doi: 10.1002/hyp.v34.2210.1002/hyp.13863

CrossRef Full Text | Google Scholar

Chang Q., Zheng T., Gao C., Zheng X., Lin Y., Song X., et al. (2024). Hydrodynamic behavior of freshwater-saltwater mixing zone in the context of subsurface physical barriers. J. Environ. Manage. 353, 120080. doi: 10.1016/j.jenvman.2024.120080

PubMed Abstract | CrossRef Full Text | Google Scholar

Chang Q., Zheng T., Zheng X., Gao C., Song X., Walther M. (2023). Repulsion driven by groundwater level difference around cutoff walls on seawater intrusion in unconfined aquifers. Sci. Total Environ. 874, 162535. doi: 10.1016/j.scitotenv.2023.162535

PubMed Abstract | CrossRef Full Text | Google Scholar

Chang Q., Zheng T., Zheng X., Zhang B., Sun Q., Walther M. (2019). Effect of subsurface dams on saltwater intrusion and fresh groundwater discharge. J. Hydrol. 576, 508–519. doi: 10.1016/j.jhydrol.2019.06.060

CrossRef Full Text | Google Scholar

Cheng H., Jiang Z., Ma X., Wang Y. (2020). Nitrogen dynamics in the mangrove sediments affected by crabs in the intertidal regions. Ecotoxicology 29, 669–675. doi: 10.1007/s10646-020-02212-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Edith A., Claudia S., Lange S., George T., Pierre L., Li J. (2015). Potential efficiency of riparian vegetated buffer strips in intercepting soluble compounds in the presence of subsurface preferential flows. PloS One 10, e0131840. doi: 10.1371/journal.pone.0131840

PubMed Abstract | CrossRef Full Text | Google Scholar

Enrique O., Muoz-Carpena R., Gao B., Fox G. (2018). Riparian vadose zone PF: review of concepts, limitations, and perspectives. Vadose Zone J. 17, 1–20. doi: 10.2136/vzj2018.02.003

CrossRef Full Text | Google Scholar

Fang Y., Qian J., Zheng T., Zheng X., Walther M. (2023). Submarine groundwater discharge in response to the construction of subsurface physical barriers in coastal aquifers. J. Hydrol. 617, 129010. doi: 10.1016/j.jhydrol.2022.129010

CrossRef Full Text | Google Scholar

Fang Y., Zheng T., Wang H., Guan R., Zheng X., Walther M. (2022b). Experimental and numerical evidence on the influence of tidal activity on the effectiveness of subsurface dams. J. Hydrol. 603, 127149. doi: 10.1016/j.jhydrol.2021.127149

CrossRef Full Text | Google Scholar

Fang Y., Zheng T., Wang H., Zheng X., Walther M. (2022a). Nitrate transport behavior behind subsurface dams under varyinghydrological conditions. Sci. Total Environ. 838, 155903. doi: 10.1016/j.scitotenv.2022.155903

PubMed Abstract | CrossRef Full Text | Google Scholar

Fanjul E., Grela M., Canepuccia A., Iribarne O. (2008). The southwest atlantic intertidal burrowing crab Neohelice granulata modifies nutrient loads of phreatic waters entering coastal area. Estuar. Coast. Shelf Sci. 79, 300–306. doi: 10.1016/j.ecss.2008.04.005

CrossRef Full Text | Google Scholar

Gao C., Kong J., Zhou L., Shen C., Wang J. (2023). Macropores and burial of dissolved organic matter affect nitrate removal in intertidal aquifers. J. Hydrol. 617, 129011. doi: 10.1016/j.jhydrol.2022.129011

CrossRef Full Text | Google Scholar

Gao S., Zheng T., Zhang B., Fang Y., Zheng X. (2024). Combined effects of aquifer heterogeneity and subsurface dam on nitrate contamination in coastal aquifers. J. Environ. Manage. 351, 119740. doi: 10.1016/j.jenvman.2023.119740

PubMed Abstract | CrossRef Full Text | Google Scholar

Gao S., Zheng T., Zheng X., Walther M. (2022). Influence of layered heterogeneity on nitrate enrichment induced by cut-off walls in coastal aquifers. J. Hydrol. 609, 127722. doi: 10.1016/j.jhydrol.2022.127722

CrossRef Full Text | Google Scholar

Geng X., Boufadel M. C. (2017). The influence of evaporation and rainfall on supratidal groundwater dynamics and salinity structure in a sandy beach. Water Resour. Res. 53, 6218–6238. doi: 10.1002/2016WR020344

CrossRef Full Text | Google Scholar

Guimond J., Seyfferth A., Moffett K., Michael H. (2020). A physical-biogeochemical mechanism for negative feedback between marsh crabs and carbon storage. Environ. Res. Lett. 15, 034024. doi: 10.1088/1748-9326/ab60e2

CrossRef Full Text | Google Scholar

Heiss J. (2020). Whale burial and organic matter impacts on biogeochemical cycling in beach aquifers and leachate fluxes to the nearshore zone - sciencedirect. J. Contam. Hydrol. 233, 103656.1–16. doi: 10.1016/j.jconhyd.2020.103656

CrossRef Full Text | Google Scholar

Heiss J., Hichael H. (2014). Saltwater-freshwater mixing dynamics in a sandy beach aquifer over tidal, spring-neap, and seasonal cycles. Water Resour. Res. 50, 6747–6766. doi: 10.1002/2014WR015574

CrossRef Full Text | Google Scholar

Heiss J., Post V., Laattoe T., Russoniello C., Michael H. (2017). Physical controls on biogeochemical processes in intertidal zones of beach aquifers. Water Resour. Res. 53, 9225–9244. doi: 10.1002/2017WR021110

CrossRef Full Text | Google Scholar

Hughes C. (1998). Characterisation of the hydrology of an estuarine wetland. J. Hydrol. 211, 34–49. doi: 10.1016/S0022-1694(98)00194-2

CrossRef Full Text | Google Scholar

Kim K., Heiss J., Michael H., Cai W., Laattoe T., Post V. (2017). Spatial patterns of groundwater biogeochemical reactivity in an intertidal beach aquifer. J. Geophys. Res.: Biogeosci. 122, 2548–2562. doi: 10.1002/2017JG003943

CrossRef Full Text | Google Scholar

Kreyns P., Geng X., Michael H. (2020). The influence of connected heterogeneity on groundwater flow and salinity distributions in coastal volcanic aquifers. J. Hydrol. 586, 124863. doi: 10.1016/j.jhydrol.2020.124863

CrossRef Full Text | Google Scholar

Kringel R., Rechenburg A., Kuitcha D., Fouepe A., Fomo M. (2016). Mass balance of nitrogen and potassium in urban groundwater in central Africa, Yaounde/Cameroon. Sci. Total Environ. 547, 382–395. doi: 10.1016/j.scitotenv.2015.12.090

PubMed Abstract | CrossRef Full Text | Google Scholar

Kuan W., Xin P., Jin G., Robinson C., Gibbes B., Li L. (2019). Combined effect of tides and varying inland groundwater input on flow and salinity distribution in unconfined coastal aquifers. Water Resour. Res. 55, 8864–8880. doi: 10.1029/2018WR024492

CrossRef Full Text | Google Scholar

Li J., Hua G., Liu S., Liu X., Huang Y., Shi Y. (2021). Effects of crab disturbance on nitrogen migration and transformation in a coastal tidal flat wetland. Environ. Sci. pollut. Res. 28, 52345–52356. doi: 10.1007/s11356-021-14393-5

CrossRef Full Text | Google Scholar

Liu Y., Jiao J., Luo X. (2016). Effects of inland water level oscillation on groundwater dynamics and land-sourced solute transport in a coastal aquifer. Coast. Eng. 114, 347–360. doi: 10.1016/j.coastaleng.2016.04.021

CrossRef Full Text | Google Scholar

Lu C., Werner A. (2013). Timescales of seawater intrusion and retreated. Adv. Water Resour. 59, 39–51. doi: 10.1016/j.advwatres.2013.05.005

CrossRef Full Text | Google Scholar

Lu J., Bai Z., Velthof G., Wu Z., Chadwick D., Ma L. (2019). Accumulation and leaching of nitrate in soils in wheat-maize production in China. Agric. Water Manage. 212, 407–415. doi: 10.1016/j.agwat.2018.08.039

CrossRef Full Text | Google Scholar

Meile C., Porubsky W., Walker R., Payne K. (2009). Natural attenuation of nitrogen loading from septic effluents: spatial and environmental controls. Water Resour. 44, 1399–1408. doi: 10.1016/j.watres.202009.11.019

CrossRef Full Text | Google Scholar

Michael H., Mulligan A., Harvey C. (2005). Seasonal oscillations in water exchange between aquifers and the coastal ocean. Nature 436, 1145–1148. doi: 10.1038/nature03935

PubMed Abstract | CrossRef Full Text | Google Scholar

Pan F., Xiao K., Guo Z., Li H. (2022). Effects of fiddler crab bioturbation on the geochemical migration and bioavailability of heavy metals in coastal wetlands. J. Hazard. Mater. 437, 129380. doi: 10.1016/j.jhazmat.2022.129380

PubMed Abstract | CrossRef Full Text | Google Scholar

Radfard M., Rahmatinia M., Tabatabaee H., Solimani H., Mahvi A., Azhdarpoor A. (2018). Data on health risk assessment to the nitrate in drinking water of rural areas in the Khashcity, Iran. Data Brief 21, 1918–1923. doi: 10.1016/j.dib.2018.11.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson C., Li L., Barry D. (2007a). Effect of tidal forcing on a subterranean estuary. Adv. Water Resour. 30, 851–865. doi: 10.1016/j.advwatres.2006.07.006

CrossRef Full Text | Google Scholar

Robinson C., Li L., Prommer H. (2007b). Tide-induced recirculation across the aquifer ocean interface. Water Resour. Res. 43. doi: 10.1029/2006WR005679

PubMed Abstract | CrossRef Full Text | Google Scholar

Robinson C., Xin P., Li L., Barry D. (2014). Groundwater flow and salt transport in a subterranean estuary driven by intensified wave conditions. Water Resour. Res. 50, 165–181. doi: 10.1002/wrcr.v50.1

CrossRef Full Text | Google Scholar

Rocha L. (2013). Seasonal enhancement of submarine groundwater discharge (SGD)-derived nitrate loading into the Ria Formosa coastal lagoon assessed by 1-D modeling of benthic NO₃^–profiles. Estuarine Coast. Shelf Sci. 132, 56–64. doi: 10.1016/j.ecss.2012.04.015

CrossRef Full Text | Google Scholar

Shen C., Zhang C., Kong J., Xin P., Lu C., Zhao Z., et al. (2019). Solute transport influenced by unstable flow in beach aquifers. Adv. Water Resour. 125, 68–81. doi: 10.1016/j.advwatres.2019.01.009

CrossRef Full Text | Google Scholar

Shen C., Zhang C., Xin P., Kong J., Li L. (2018). Salt dynamics in coastal marshes: Formation of hypersaline zones. Water Resour. Res. 54, 3259–3276. doi: 10.1029/2017WR022021

CrossRef Full Text | Google Scholar

Shuai P., Cardenas M., Knappett P., Bennett P., Neilson B. (2017). Denitrification in the banks of fluctuating rivers: the effects of river stage amplitude, sediment hydraulic conductivity and dispersivity, and ambient groundwater flow. Water Resour. Res. 53, 7951–7967. doi: 10.1002/2017WR020610

CrossRef Full Text | Google Scholar

Spiteri C., Slomp C., Charette M., Tuncay K., Meile C. (2008). Flow and nutrient dynamics in a subterranean estuary (Waquoit Bay, MA, USA): Field data and reactive transport modeling. Geochim. Cosmochim. Acta 72, 3398–3412. doi: 10.1016/j.gca.2008.04.027

CrossRef Full Text | Google Scholar

Sun Q., Zheng T., Zheng X., Chang Q., Walther M. (2019). Influence of a subsurface cutoff wall on nitrate contamination in an unconfined aquifer. J. Hydrol. 575, 234–243. doi: 10.1016/j.jhydrol.2019.05.030

CrossRef Full Text | Google Scholar

Sun Q., Zheng T., Zheng X., Walther M. (2021). Effects of physical barrier on seawater intrusion and nitrate accumulation in upstream aquifers. J. Contam. Hydrol. 243, 103913. doi: 10.1016/j.jconhyd.2021.103913

PubMed Abstract | CrossRef Full Text | Google Scholar

van Genuchten M. (1980). A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. America J. 44, 892–898. doi: 10.2136/sssaj1980.03615995004400050002x

CrossRef Full Text | Google Scholar

Williams M., Buda A., Elliott H., Hamlett J., Boyer E., Schmidt J. (2014). Groundwater flow path dynamics and nitrogen transport potential in the riparian zone of an agricultural headwater catchment. J. Hydrol. 511, 870–879. doi: 10.1016/j.jhydrol.2014.02.033

CrossRef Full Text | Google Scholar

Xiao K., Pan F., Santos I., Zheng Y., Zheng C., Chen N., et al. (2022). Crab bioturbation drives coupled iron-phosphate-sulfide cycling in mangrove and salt marsh soils. Geoderma 424, 115990. doi: 10.1016/j.geoderma.2022.115990

CrossRef Full Text | Google Scholar

Xiao K., Wilson A., Li H. (2019). Crab burrows as preferential flow conduits for groundwater flow and transport in salt marshes: A modeling study. Adv. Water Resour. 132, 103408. doi: 10.1016/j.advwatres.2019.103408

CrossRef Full Text | Google Scholar

Xie Y., Wang Y., Zhang J., Ye Y., Shen C., Zeng Y., et al. (2023). Seawater intrusion in fractured coastal aquifers influenced by tides: Laboratory and numerical investigations. J. Hydrol. 622, 129637. doi: 10.1016/j.jhydrol.2023.129637

CrossRef Full Text | Google Scholar

Xin P., Jin G., Li L., Barry D. (2009). Effects of crab burrows on pore water flows in salt marshes. Adv. Water Resour. 32, 439–449. doi: 10.1016/j.advwatres.2008.12.008

CrossRef Full Text | Google Scholar

Xin P., Robinson C., Li L., Barry D., Bakhtyar R. (2010). Effects of wave forcing on a subterranean estuary. Water Resour. Res. 46, W12505. doi: 10.1029/2010WR009632

CrossRef Full Text | Google Scholar

Xu X., Xin P., Zhou T., Xiao K. (2021). Effect of macropores on pore-water flow and soil conditions in salt marshes subject to evaporation and tides. Estuarine Coast. Shelf Sci. 261, 107558. doi: 10.1016/j.ecss.2021.107558

CrossRef Full Text | Google Scholar

Ying Z. (2021). Effect of macropores on water and salt exchange in coastal wetlands (China: Hohai University, Master's thesis). 11p.

Google Scholar

Zhang W., Tian Z., Zhang N., Li X. (1996). Nitrate pollution of groundwater in northern China. Agric. Ecosyst. Environ. 59, 223–231. doi: 10.1016/0167-8809(96)01052-3

CrossRef Full Text | Google Scholar

Zhang B., Zheng X., Zheng T., Xin J., Sui S., Zhang D. (2019). The influence of slope collapse on water exchange between a pit lake and a heterogeneous aquifer. Front. Environ. Sci. Eng. 13, 9. doi: 10.1007/s11783-019-1104-9

CrossRef Full Text | Google Scholar

Zheng T., Gao M., Chang Q., Zheng X., Walther M. (2022). Dynamic desalination of intruding seawater after construction of cut-off walls in a coastal unconfined aquifer. Front. Mar. Sci. 9. doi: 10.3389/fmars.2022.857807

PubMed Abstract | CrossRef Full Text | Google Scholar

Zheng T., Yuan F., Gao S., Zheng X., Liu T., Luo J. (2023). The impact of hydraulic conductivity anisotropy on the effectiveness of subsurface dam. J. Hydrol. 626, 130360. doi: 10.1016/j.jhydrol.2023.130360

CrossRef Full Text | Google Scholar

Zheng T., Zheng X., Chang Q., Zhan H. (2021). Timescale and effectiveness of residual saltwater desalinization behind subsurface dams in an unconfined aquifer. Water Resour. Res. 57, e2020WR028493. doi: 10.1029/2020WR028493

CrossRef Full Text | Google Scholar

Zheng T., Zheng X., Sun Q., Wang L., Walther M. (2020). Insights of variable permeability full-section wall for enhanced control of seawater intrusion and nitrate contamination in unconfined aquifers. J. Hydrol. 586, 124831. doi: 10.1016/j.jhydrol.2020.124831

CrossRef Full Text | Google Scholar