1. Introduction

Beaches mark the terminal end of terrestrial groundwater flow paths and comprise approximately 1/3 of the world’s ice-free coastline ^[1]. High rates of groundwater-surface water exchange and chemical fluxes across the permeable intertidal sand surface fuels biogeochemical transformations of ocean- and land-derived solutes in groundwater discharging to the coastal ocean ^[2][3][4][5]. Thus, it is important to investigate beach groundwater systems in a framework that considers the role of flow rates and patterns on controlling the spatial and temporal variability of dissolved species concentrations and porewater reactivity. Insights into hydro-biogeochemical processes within beach aquifers obtained using multiple methods (

Figure 1

Figure 1.

 Schematic of a beach aquifer with an intertidal circulation cell. Several instrumentation and methodological examples are illustrated, including augering and vibracoring, sensor networks, multi-level porewater samplers, and slotted sand column experiments.

Figure 1). Fresh groundwater flows through the beach aquifer before discharging to the coastal ocean as fresh submarine groundwater discharge (SGD). Due to high nutrient concentrations from agricultural practices and the prevalence of private septic systems along coastlines ^[6], nutrient loads in fresh SGD can approach or exceed riverine fluxes ^[7]. As a result, fresh SGD can degrade water quality and adversely impact ecosystem functioning, species composition, and diversity in coastal and estuarine waters ^[8][9][10].

Freshwater-saltwater (FW-SW) mixing in the intertidal zone was first investigated in the landmark work by Lebbe ^[11], where electrical resistivity logs revealed an “upper saltwater lens” overlaying a freshwater zone in a wide beach in Belgium. Since Lebbe ^[11], the upper saltwater lens has also been termed the saltwater (seawater) cell, tide/wave-driven nearshore circulation, intertidal circulation cell, upper saline plume, subterranean estuary, and tidally driven recirculation ^{[12][13][14][15]}. Here, we adopt the term “intertidal circulation cell” to describe the region of elevated salinity beneath the beach face because it describes: 1) the cross-shore coastal area where the hydrologic and chemical processes occur, 2) the flow pattern within the intertidal zone, and 3) the general geometry of the region of elevated salinity. The intertidal circulation cell forms as a result of wave and tidal input of seawater across the upper sandy beach face, leading to a super-elevation of the beach water table above mean sea level near the high tide line (e.g., ^[16]) and mixing between infiltrated seawater and underlying fresh groundwater. The saltwater is forced downward and seaward in a circulating pattern by the terrestrial freshwater hydraulic gradient and the hydraulic gradient generated by the water table super-elevation. Farther seaward near the low tide line, a fresh groundwater discharge zone separates the intertidal circulation cell from the lower saltwater-freshwater interface (

Figure 1). Since the early works that described the physical characteristics of the intertidal flow system (e.g., ^{[11][12][17][18][19]}, geochemical and microbial properties of the intertidal circulation cell have gained attention, with a marked increase in the total number of publications and citations in the mid-2000s (

Figure 2

Figure 2.

 Number of peer-reviewed publications and citations on beach groundwater systems. The counts were obtained by conducting a search using the keywords “beach aquifer” and “subterranean estuary” using the Web of Science database. From 1991 to 2019, a total of 796 papers have been published totaling 19,988 citations.

Intertidal freshwater-seawater mixing zones host redox gradients that develop as a result of the contrasting geochemical signatures of the fresh and saline endmembers. The redox gradients support abiotic or microbially mediated biogeochemical transformation of solutes in groundwater prior to discharge ^[2][20][21]. A wide range of biogeochemical reactions have been observed in beaches, including aerobic respiration ^[4], dissolved organic carbon degradation ^[22][23], denitrification ^[24], sulfate reduction ^[21], iron oxidation-reduction ^[20][23][25], nitrification ^[26], ammonification ^[27][28], and annamox ^[29]. The formation and spatial distribution of reactive zones in beach aquifers are linked to the physical flow and mixing processes driven by the freshwater hydraulic gradient and ocean forcing ^[15][30][31].

While chemical fluxes from rivers to the ocean are relatively well-constrained, quantifying groundwater nutrient loads to coastal water bodies is challenging because fluid and chemical fluxes occur along the global coastline, are diffuse, and are spatially heterogeneous ^{[32][33][34][35][36]}. Further, hydrological shifts in flow and solute transport in beach aquifers, occurring from wave to seasonal time scales ^[18][37][38], lead to dynamic biogeochemical cycling of nutrients, heavy metals, and other elements in intertidal porewater ^{[39][40][41][42][43][44]}. As a result, large uncertainties remain regarding the magnitude and spatial variability of chemical transformations within beach aquifers and associated chemical fluxes to coastal ecosystems.

In the past two decades, a wealth of measurement approaches has been applied, modified, or developed to investigate the interplay between groundwater flow processes and chemical cycling in beaches. However, instrumenting and sampling beach groundwater systems remains a difficult task. Timescales of hydrologic variability (e.g., waves, tides) can be short (seconds to hours), making sampling difficult when aiming to resolve variability due to ocean forcing. Beach morphology cycles through accretion and erosional periods and can affect the distribution of the FW-SW mixing zone (e.g., ^[45][46]). Erosion of the beach face also imposes issues related to instrument exposure, including vandalism, biofouling, or destruction due to the action of waves, tides, currents, or water-borne ice. Owing to the lack of cohesiveness of sandy sediments, instrument deployment of can also be problematic, as sand collapses in uncased boreholes immediately after sediment is removed below water table. For this reason, most studies that use hand-auguring to install equipment are limited in spatial coverage to the upper few meters of the intertidal subsurface.

The delivery of labile marine organic carbon from infiltrating seawater is an important mechanism that supports reactivity in carbon-poor beaches While dissolved organic carbon (DOC) is typically assumed to be most reactive, particulate forms of carbon (POC; particulate organic carbon) are also capable of supporting reactions. Beach wrack, algal fragments, and buried marine carcasses (e.g., whale strandings) are viable forms of in-situ organic carbon that can support biogeochemical processes, either in particulate form or as a source of leached DOC ^{[41][47][48][49][50][51][52]}. This introduces an additional layer of complexity that affects the hydrogeochemistry of the beach. Particulate matter is transported through intertidal sediments at rates that differ from advective transport rates due to filtration and sorption processes occurring along flow paths. The filtration, sorption, and mobility of particulates sourced from fresh groundwater or seawater will depend on sediment characteristics, geologic heterogeneity, anisotropy, particulate charge and size, and the groundwater flow rate. Recent field studies by Kim et al. ^[41][42] showed that particulate organic carbon can be heterogeneously distributed within beach sediments and experience different transport dynamics from dissolved solutes. Results indicated that pools of POC can be intermittently used as an electron donor in response to their varying overlap with other solute reactants as groundwater flow patterns shift. Thus, it is important to characterize the physical and electrochemical properties of mobile organic particulates in fresh and saline groundwater.

2. Measuring the Physical Parameters

2.1. Approaches for Characterizing Beach Geology

Regional geologic assessments can help ensure that the selected along-shore position of the measurement transect is a representative cross-section. Previous work has shown that cross-cutting paleochannels and related alongshore variability in shallow confining units can control offshore (<500 m) subsurface salinity distributions ^[34][53] and freshwater discharge of nitrate, ammonium, and dissolved iron ^[35]. These studies suggest that beaches that overly paleochannels and low permeability caps may contain flow patterns that are not representative of a non-channelized coastline. Electrical resistivity tomography ^[54], time-domain electromagnetic surveys ^[55], seismic profiles ^[27], temperature surveys ^[56], and combinations of geophysical methods ^[57] of the nearshore area can provide valuable guidance during site selection to reduce the likelihood of instrumenting a site with large-scale geologic features that could otherwise lead to inaccurate interpretations of measured parameters. Therefore, it is important to identify and, if present, characterize the underlying major geologic features.

Regional geologic assessments can help ensure that the selected along-shore position of the measurement transect is a representative cross-section. Previous work has shown that cross-cutting paleochannels and related alongshore variability in shallow confining units can control offshore (<500 m) subsurface salinity distributions [34,53] and freshwater discharge of nitrate, ammonium, and dissolved iron [35]. These studies suggest that beaches that overly paleochannels and low permeability caps may contain flow patterns that are not representative of a non-channelized coastline. Electrical resistivity tomography [54], time-domain electromagnetic surveys [55], seismic profiles [27], temperature surveys [56], and combinations of geophysical methods [57] of the nearshore area can provide valuable guidance during site selection to reduce the likelihood of instrumenting a site with large-scale geologic features that could otherwise lead to inaccurate interpretations of measured parameters. Therefore, it is important to identify and, if present, characterize the underlying major geologic features.

Small-scale geologic features affect the dynamics of the intertidal circulation cell. Recent modeling studies have shown that spatial variability in hydraulic conductivity on the order of centimeters to meters is a major control of flow, transport, and reactivity in beach and riparian aquifers ^{[52][58][59][60][61]}, highlighting the importance of characterizing site-specific hydrogeologic properties (e.g., porosity, grainsize, hydraulic conductivity). Vibracores provide valuable information about grainsize distributions and potential field sampling strategies in mixed sediments ^[62][63]. It is important to select a core diameter that can to adequately core through pebbly layers without clogging, along with a core catcher to prevent sediment loss from the core during retrieval ^[64][65]. Core catchers can be riveted at the end of the barrel to further minimize sediment loss. Sands attenuate the vibration of the barrel and reduce penetration depths and core lengths to 2–3 m. Modifications to traditional vibracoring techniques can be applied to improve core penetration and recovery in sandy sediments. Cores of 10–12 m can be recovered using a series of short closely spaced cores, a water supply, and pumps ^[66].

Small-scale geologic features affect the dynamics of the intertidal circulation cell. Recent modeling studies have shown that spatial variability in hydraulic conductivity on the order of centimeters to meters is a major control of flow, transport, and reactivity in beach and riparian aquifers [52,58,59,60,61], highlighting the importance of characterizing site-specific hydrogeologic properties (e.g., porosity, grainsize, hydraulic conductivity). Vibracores provide valuable information about grainsize distributions and potential field sampling strategies in mixed sediments [62,63]. It is important to select a core diameter that can to adequately core through pebbly layers without clogging, along with a core catcher to prevent sediment loss from the core during retrieval [64,65]. Core catchers can be riveted at the end of the barrel to further minimize sediment loss. Sands attenuate the vibration of the barrel and reduce penetration depths and core lengths to 2–3 m. Modifications to traditional vibracoring techniques can be applied to improve core penetration and recovery in sandy sediments. Cores of 10–12 m can be recovered using a series of short closely spaced cores, a water supply, and pumps [66].

Augering through casing is another low-cost method of characterizing site hydrogeology. Two casing sizes may be telescoped to achieve greater depths. For example, on a sandy beach at Cape Shores, DE, USA a 12 cm outer diameter auger was used to auger through a 1.5 m length by 15 cm OD PVC casing [45]. A 3.0 m length by 10 cm OD casing was then inserted into the cased borehole and a 7 cm OD auger was used to auger through the inner casing to 3.0 m. As this method produces a cased borehole to the desired target depth, porewater sampling equipment can be installed and the casing extracted to allow the sediment to backfill around the sampling ports.

2.2. Beach Profiling

Beach topography can vary up to several meters over days to months (e.g., [67]) due to periodic and episodic forcings and should be considered to minimize the likelihood of equipment burial or damage from exposure to waves. Beaches used and managed for recreational purposes may also be intermittently replenished and often undergo scraping where large machinery is used to rake the surficial foreshore and backshore sediments to remove beach wrack, fishing lines, and woody debris (

Figure 3

, right). These anthropogenic modifications to the composition and distribution of beach sediments disrupt the natural beach profile, remove and rearrange particulate organic matter (e.g., beach wrack), and may homogenize surficial sediments. Knowledge of the existence and frequency of these events will aid interpretation of abrupt shifts in flow patterns or solute concentrations.

Figure 3.

 (

Left

): Seasonal beach profile evolution at Cape Shores, DE, USA and example sampling port locations. Topographic lows that form in the backshore (x = ~30 m) during winter can allow seawater to pool and serve as a continuous source of saltwater to the aquifer. (

Right

): Beaches used recreationally may be subject to beach management practices, such as raking, artificially altering beach topography.

Measurements and numerical modeling studies have demonstrated that the beach profile affects solute transport paths ^[68][69][70] and aquifer reactivity ^[31]. Given these topographic controls, it is important to monitor the sand surface elevation over the duration of the sampling campaign. Light Detection and Ranging (LiDAR) ^[71] and Interferometric Synthetic Aperture Radar (InSAR) ^[72] can provide highly accurate elevation measurements, however leveling using the Emery Method with survey poles is often sufficient to obtain a beach profile, as local topographic variations tend to be greater than measurement errors ^[73][74]. Seasonal topographic variations at a beach at Cape Shores, DE, USA required seasonal topographic surveys (

Measurements and numerical modeling studies have demonstrated that the beach profile affects solute transport paths [68,69,70] and aquifer reactivity [31]. Given these topographic controls, it is important to monitor the sand surface elevation over the duration of the sampling campaign. Light Detection and Ranging (LiDAR) [71] and Interferometric Synthetic Aperture Radar (InSAR) [72] can provide highly accurate elevation measurements, however leveling using the Emery Method with survey poles is often sufficient to obtain a beach profile, as local topographic variations tend to be greater than measurement errors [73,74]. Seasonal topographic variations at a beach at Cape Shores, DE, USA required seasonal topographic surveys (

Figure 3

, left), which explained that elevated salinity observed below the upper beach during winter was caused by seawater pooling in a topographic low near the high tide line [45]. As wells and porewater samplers are often installed manually and thus to limited (~3 m) depths, large (>3 m; [75]) erosive events provide a window of opportunity to install instruments to greater depths relative to mean sea level. The effective depth of the sample port locations will increase as the beach face is naturally replenished. In these instances, it is important to extend casings to a height that is above the future projected profile.

3. Biogeochemical Characterization of the Intertidal Zone

3.1. Distribution of Nutrients and Dissolved Constituents

For studies with a focus on oxygen consumption, the Mettler Toledo InLab OptiOx (Columbus, OH, USA), WTW ProiLine Oxi Portable Oxygen Meter (Weilheim, Germany), and the PyroScience FireStingO2 (Aachen, Germany) provide high-accuracy dissolved oxygen measurements that can be used in-line during pumping. Other sensors, such as the Aanderaa Oxygen Optode (Bergen, Norway) can be wired to dataloggers and installed in the subsurface for continuous dissolved oxygen monitoring (e.g., [43]).

Dissolved chemical constituents that require laboratory analyses can be sampled using the multi-level porewater samplers described above. A filtration system using Glass Fiber Filters (GFF; ~0.7 μm) and filter holders (filter membranes) can be attached directly to the peristaltic pump using luer-lock plastic nozzles. To preserve the redox state or to avoid gas stripping, inertial or bladder pumps may be used ^[76]. The diameter and pore size of the filter can be chosen depending on the sampled material. GFFs may be suitable for filtering larger particulate matter (sediments or algal fragments), but sampling for colloidal, bacterial, or genetic matter (0.1–0.4 µm) may require a finer mesh pore size or specialized filter systems (e.g., Sterivex filters, Darmstadt, Germany) to obtain concentrated DNA (DeoxyriboNucleic Acid) samples. GFFs must be combusted (450 °C, 4–5 h) prior to sampling nutrient and carbon species to prevent contamination. Sterile filters pre-loaded into disposable filter holders may also be used. Reusable polycarbonate filter holders can become brittle over time due to frequent acid washing, which can lead to damaged threads and leakage during filtration. Polypropylene filter holders tend to be more durable under frequent use. Filter holders are also readily available in stainless steel and silicone.

Dissolved chemical constituents that require laboratory analyses can be sampled using the multi-level porewater samplers described above. A filtration system using Glass Fiber Filters (GFF; ~0.7 μm) and filter holders (filter membranes) can be attached directly to the peristaltic pump using luer-lock plastic nozzles. To preserve the redox state or to avoid gas stripping, inertial or bladder pumps may be used [87]. The diameter and pore size of the filter can be chosen depending on the sampled material. GFFs may be suitable for filtering larger particulate matter (sediments or algal fragments), but sampling for colloidal, bacterial, or genetic matter (0.1–0.4 µm) may require a finer mesh pore size or specialized filter systems (e.g., Sterivex filters, Darmstadt, Germany) to obtain concentrated DNA (DeoxyriboNucleic Acid) samples. GFFs must be combusted (450 °C, 4–5 h) prior to sampling nutrient and carbon species to prevent contamination. Sterile filters pre-loaded into disposable filter holders may also be used. Reusable polycarbonate filter holders can become brittle over time due to frequent acid washing, which can lead to damaged threads and leakage during filtration. Polypropylene filter holders tend to be more durable under frequent use. Filter holders are also readily available in stainless steel and silicone.

In eutrophic systems where nutrient loads and ocean primary production are high, measurements of N speciation (NO

₃⁻

, NO

₂⁻

, NH

₄⁺

, organic N), P (PO

₄³⁻

), and C (organic and inorganic, including chlorophyll-

a

) can be used to identify nutrient sources and types of chemical transformations occurring in the aquifer. Concentration measurements can also be used to quantify chemical fluxes to surface water when paired with hydraulic or seepage meter measurements [35]. Instruments such as the multi-channel SEAL Analyzer (Mequon, WI, USA) can be used to measure a suite of nutrient concentrations, including NO

₃⁻

, NO

₂⁻

, NH

₄⁺

, PO

₄³⁻

, and dissolved Si in the laboratory. Other constituents such as chlorophyll, Fe(II), total Fe, and H

₂

S can be measured using spectrophotometric methods (e.g., Hach portable spectrophotometer). While some constituents such as NO

₃⁻

, NH

₄⁺

, chlorophyll, and phycocyanin can be measured in situ using sondes (e.g., In-Situ Aqua TROLL Multiparameter Sonde, Fort Collins, CO, USA), sondes can require large sample volumes and are large in size relative to Conductivity-Temperature-Depth (CTD) sensors.

In carbon-poor beach systems, the source and reactivity of organic carbon (both dissolved and particulate; see 

Section 3.2 for particulate organic carbon) plays a dominant role in biogeochemical reactivity. While chlorophyll is indicative of young, labile carbon of marine origin, terrestrial organic carbon can be refractory due to a more complex molecular structure and longer transport paths ^[22][77]. The distribution of carbon types in the intertidal circulation cell can therefore be used to trace groundwater flow paths and sources of solutes being delivered to the beach, as well as the biogeochemical potential of the aquifer. Chlorophyll can be detected and quantified using fluorometers. Total carbon (TC), total organic carbon (TOC), and total inorganic carbon (TIC) can be quantified using TOC analyzers, such as the Shimadzu TOC Analyzer or Sievers TOC Analyzer. These and similar instruments measure TC by oxidizing the sample via combustion and measuring generated CO

 for particulate organic carbon) plays a dominant role in biogeochemical reactivity. While chlorophyll is indicative of young, labile carbon of marine origin, terrestrial organic carbon can be refractory due to a more complex molecular structure and longer transport paths [22,88]. The distribution of carbon types in the intertidal circulation cell can therefore be used to trace groundwater flow paths and sources of solutes being delivered to the beach, as well as the biogeochemical potential of the aquifer. Chlorophyll can be detected and quantified using fluorometers. Total carbon (TC), total organic carbon (TOC), and total inorganic carbon (TIC) can be quantified using TOC analyzers, such as the Shimadzu TOC Analyzer or Sievers TOC Analyzer. These and similar instruments measure TC by oxidizing the sample via combustion and measuring generated CO

₂

 with a gas detector or conductometer. TIC is measured by volatilizing inorganic carbon via acidification. TOC is quantified as the difference between the two measurements.

Characterization of dissolved inorganic carbon (DIC) chemistry and total alkalinity (T-Alk) can help elucidate reaction dynamics in systems with active carbonate dissolution or high rates of anaerobic reactivity. For example, aerobic respiration, denitrification, and sulfate reduction alter DIC:T-Alk ratios, so tracking changes to DIC:T-Alk spatially or over a range of salinities can complement other solute distribution datasets ^[78][79]. Porewater samples intended for DIC and T-Alk analysis should be fixed with HgCl

Characterization of dissolved inorganic carbon (DIC) chemistry and total alkalinity (T-Alk) can help elucidate reaction dynamics in systems with active carbonate dissolution or high rates of anaerobic reactivity. For example, aerobic respiration, denitrification, and sulfate reduction alter DIC:T-Alk ratios, so tracking changes to DIC:T-Alk spatially or over a range of salinities can complement other solute distribution datasets [89,90]. Porewater samples intended for DIC and T-Alk analysis should be fixed with HgCl

₂ immediately in the field and can be measured using DIC analyzers such as those offered by ApolloSci Tech ^[80]. We note that adding HgCl

 immediately in the field and can be measured using DIC analyzers such as those offered by ApolloSci Tech [91]. We note that adding HgCl

₂

 to high-sulfate samples can precipitate HgS and interfere with T-Alk measurements, thus high-salinity samples should be tested for sulfate concentrations prior to fixing with HgCl

₂

.