Methods in Capturing the Spatiotemporal Dynamics of Flow and Biogeochemical Reactivity in Sandy Beach Aquifers: A Review

Sandy beach aquifers are complex hydrological and biogeochemical systems where fresh groundwater and seawater mix. The extent of the intertidal mixing zone and the rates of circulating flows within beaches are a primary control on porewater chemistry and microbiology of the intertidal subsurface. Interplay between the hydrological and biogeochemical processes at these landsea transition zones moderate fluxes of chemicals, particulates, heavy metals, and biota across the aquifer-ocean interface, affecting coastal water quality and nutrient loads to marine ecosystems. Thus, it is important to characterize hydrological and biogeochemical processes in beach aquifers when estimating material fluxes to the ocean. This can be achieved through a suite of cross-disciplinary measurements of beach groundwater flow and chemistry. In this review, we present measurement approaches that have been developed and employed to characterize the physical (geology, topography, subsurface hydrology) and biogeochemical (solute and particulate distributions, reaction rates) properties of and processes occurring within sandy intertidal aquifers. As applied to beach systems, we discuss vibracoring, sample collection, laboratory experiments, variable-density considerations, instrument construction, and sensor technologies. We discuss advantages and limitations of typical hydrologic field sampling methods when used to investigate beach aquifers and provide a measurement framework for researchers seeking to sample and collect data from these systems.


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) can provide a stronger basis for quantifying groundwaterborne chemical fluxes to coastal surface water ecosystems.
Fresh groundwater delivers land-derived nutrients such as nitrate, phosphate, dissolved silica, and organic carbon (often refractory due to their complex structure; Seidel et al., 2015) to coastal aquifers ( 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 Fresh groundwater delivers land-derived nutrients such as nitrate, phosphate, dissolved silica, and organic carbon (often refractory due to their complex structure; Seidel et al., 2015) to coastal aquifers ( 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/wavedriven 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 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/wavedriven 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). 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 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 oxidationreduction [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.
In this review, we discuss sampling methods that have been developed and applied in literature to characterize the hydrogeology of intertidal sediments and the hydrological and chemical dynamics of beach aquifers. While the methods presented herein may not be applicable to all field sites, many of the techniques can be easily modified to obtain greater site-specific spatial or temporal sampling coverage. A sampling framework that incorporates a combination of the hydrological and chemical sampling techniques discussed in this review will greatly aid interpretation of geochemical signatures and will strengthen understanding of the relationships between the hydrological, biogeochemical, and microbial functioning of beach aquifers.

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 largescale 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].
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.

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.
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.

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. 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 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 Water 2021, 13, 782 6 of 16 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.

Porewater Sampling and Head Measurements
Multi-level samplers and piezometers are typically installed perpendicular to the shoreline from near the high tide line to tens of meters seaward of the low tide line to capture the cross-shore extent of the intertidal circulation cell, fresh discharge zone, and lower saltwater-freshwater interface (Figure 1, e.g., [13,17,68]). In settings with broad intertidal sand flats, the transect may need to be extended seaward of the beach face slope break > 100 m to capture fresh and saltwater fluxes across large ridge and tough features (e.g., [76]). To aid in selecting the seaward extent of the transect, shallow PushPoint Sippers (M.H.E Products, East Tawas, MI, USA) can be used to sample porewater for salinity in the upper 1.0 m of the sand flat or subtidal zone to determine the location of the fresh discharge zone and lower interface. PushPoint Sippers are inexpensive sampling systems that comprise a 0.3 or 0.6 cm OD 316 stainless steel pipe with~3 cm slotted screens and a handle that aids insertion into sediment. The strengths of PushPoint Sippers are the ease of use and ability to collect many shallow porewater samples quickly. Other sampling methods are warranted for deeper sampling depths.
There are primarily three types of equipment used for porewater sampling in beach sediments: (1) Stainless steel AMS Inc. (American Falls, ID, USA) drive-point gas samplers (temporary or permanent), (2) PVC and tubing multi-level samplers, and (3) nested piezometers. AMS gas samplers were initially designed to monitor gas vapors in the unsaturated zone and have been repurposed by coastal researchers to obtain porewater profiles [77]. The gas sampling kits have two options for sample collection: (1) a reusable AMS Retract-a-Tip for sampling at multiple depths in a single installation and (2) a permanent disposable tip that is installed at the target depth. Both systems utilize hollow 4130 alloy steel extension rods with internal tubing that extends from the screened interval (Retract-a-Tip or dedicated tip) to the top of the extension rods. The Retract-a-Tip consists of a stainless-steel sliding cover that is closed during installation and opens when the extension rod is lifted to expose a 6 cm screen. Multiple porewater samples can be collected at various target depths during installation by lifting the extension rods slightly to expose the screen before advancing to deeper sampling depths. The Retract-a-Tip is recovered with the removal of the extension rods. The dedicated tips rely on the same stainless extension rods, but are not threaded to the leading extension rod; the tips are held in place during installation by the weight of the extension rods. Once the desired depth is reached, the extension rods are removed and the dedicated tip as well as the attached internal tubing remains in the sediment. The Retract-a-Tip system is ideal for exploratory sampling because multiple depths can be sampled rapidly in a single installation. If porewater is to be monitored over time, then a separate Retract-a-Tip installation is required for each sampling campaign, and it is unlikely that the Retract-A-Tip screen will be inserted to same location in the aquifer across multiple surveys. The permanent dedicated tips are better suited for long-term monitoring at the expense of added labor and cost, as a separate installation and tip is needed for each port. The design of the AMS system has been modified for use with nylon and silicone parts to eliminate the presence of metal in the subsurface that may lead to heavy metal contamination [78].
Multi-level samplers constructed of PVC and tubing offer an inexpensive and metalfree approach to sample beach porewater. Heiss et al., [45] constructed permanent multilevel porewater samplers using sections of 7 mm OD polyethylene tubing that are attached to the outside of a 13 mm OD PVC pipe used as structural support. The ends of each section of tubing are screened over 2-3 cm and covered with fine nylon mesh. The polyethylene tubing can also be fed inside a larger diameter PVC pipe that has screened intervals that open to the bottom of each section of tubing [79]. The cross-shore width and depth of the multi-level sampling transect should be designed to capture the expected spatial extent of the targeted solute over the monitoring period, considering variability due to wave and tidal forcing, fresh groundwater input, and expansion and contraction of the circulation cell resulting from variability in the width of the intertidal zone due to beach profile changes. We have successfully maintained these systems for >9 years in a sandy beach and found it helpful to place the top of the samplers below the sand surface to prevent vandalism and damage from beach scraping.
Traditional nested piezometers are a third approach for sampling porewater [79]. However, nested piezometers can dampen rapid water table fluctuations and are not recommended for monitoring hydraulic responses due to wave forcing [79,80].
Motorized tools or manual methods can be used to install sampling equipment. Percussion drills, fence post drivers, or vibracores allow for installations to depths of 3-10 m [77,78,81], while manual methods such as push-point sampling, augering, or fence post hammering is generally limited to 1-3 m in sandy sediments [82,83].

Hydraulic Gradients and Ocean Water Levels
Measurement of water levels in piezometers installed in the beach can be used to determine the direction and magnitude of groundwater flow. Horizontal flow can be calculated from the difference in heads between two wells and the well separation distance. However, it is critical that the hydraulic head be converted to freshwater equivalent head and that the hydraulic gradient is evaluated from freshwater heads at the same elevation. This is because unlike uniform density systems, freshwater head may vary with depth in variable-density systems, even under hydrostatic conditions. The conversion from hydraulic head to freshwater equivalent head h f (L) is given by: where P i (M L −1 T −2 ) is the pressure at the piezometer screen, ρ f (M L −3 ) is the density of freshwater, g (L T −2 ) is gravitational acceleration, z i (L) is elevation head. Using Equation (1), horizontal flow q x (volume of water per unit cross-sectional area per unit time; L 3 L −2 T −1 ) in variable density systems can then be expressed as [84]: where k (L 2 ) is the intrinsic permeability, µ f is dynamic viscosity of fresh water (M L −1 T −1 ), µ is the dynamic viscosity of the porewater (M L −1 T −1 ), and K fx (L T −1 ) is the freshwater hydraulic conductivity in the horizontal direction. For most applications, it can be assumed that µ f /µ = 1 [85]. Simple conversion to freshwater heads is insufficient for calculating vertical flow because freshwater head is a function of depth in saltwater; two piezometers screened at different depths in a hydrostatic body of saltwater will produce different freshwater heads. To calculate vertical flow q z (L 3 L −2 T −1 ) the term (ρ − ρ f )/ρ f must be considered to account for buoyancy effects on vertical flow [84]: where ρ (M L −3 ) is the density of the pore water and K f z (L T −1 ) is the vertical hydraulic conductivity. Due to the added complexities of groundwater flow in variable-density systems, hydraulic head measurements in beach aquifers should be used cautiously to avoid misinterpretation of flow rates and directions. In addition to cased piezometers, a hydraulic potentiomanometer can be used to measure hydraulic gradients. The salinity of the porewater must be measured for conversion to freshwater head and buoyancy must be considered if salinity variations are suspected (Equation (3)). Hydraulic potentiomanometers can consist of any number of channels and have successfully been constructed and deployed with 18 channels connected to mini-piezometers in a sandy coastal aquifer [86].
In comparison to manual electric tape or potentiomanometer measurements, sensors are more appropriate when high-frequency measurements are required. A sample frequency of 1 2 the wave period is adequate to resolve wave-induced water level fluctuations. However, at these high sample frequencies (i.e., <~1 s), sensor options may be limited from manufacturers that specialize in hydrological monitoring. While few manufacturers (e.g., Solinst Inc., Georgetown, ON, USA) do produce pressure transducers that can record at intervals down to 1/8 s, manufacturers that cater to oceanographic applications produce sensors that are capable of sampling at very high frequencies (e.g., 32 Hz).

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 [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 3 − , NO 2 − , NH 4 + , organic N), P (PO 4 3− ), 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 3 − , NO 2 − , NH 4 + , PO 4 3− , and dissolved Si in the laboratory. Other constituents such as chlorophyll, Fe(II), total Fe, and H 2 S can be measured using spectrophotometric methods (e.g., Hach portable spectrophotometer). While some constituents such as NO 3 − , NH 4 + , 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,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 2 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 [89,90]. Porewater samples intended for DIC and T-Alk analysis should be fixed with HgCl 2 immediately in the field and can be measured using DIC analyzers such as those offered by ApolloSci Tech [91]. We note that adding HgCl 2 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 2 .

Quantifying Particulate Organic Carbon
Particulate matter is transported more slowly than the advective transport rate due to sorption and filtration processes, leading to challenges when sampling and tracking its transport relative to dissolved constituents. Laboratory flow-through column experiments using aquifer sediments or porous media with field-equivalent grainsizes can provide insight into retardation rates and adsorption behavior of particulates. Gast et al., [92] performed column experiments using fluorescent 1-µm plastic microspheres (Polyscience Fluoresbrite) to simulate bacterial transport in beach sands. The results demonstrated that the transport and adsorption behavior of the plastic microspheres were similar to enterococci. In the field, the microspheres were buried to a depth of~5 cm in the intertidal zone and were transported rapidly laterally and vertically by wave swash.
Slotted sand columns can be used in the field to investigate the influx of marine organic matter across the beach face. The columns can be constructed from slotted PVC and are pre-filled with homogenized and combusted beach sediments prior to installation ( Figure 1). Kim et al., [41] deployed 3 m long slotted sand columns in the intertidal zone adjacent to a porewater sampling transect. The columns were retrieved and sectioned once the sediment and porewater in the columns equilibrated with the surrounding aquifer (2 months). The retrieved sediment samples exhibited a fine layer of grey matter not present in control samples, indicating active delivery and transport of particulate matter. Microscopy indicated that the grey matter consisted of bacterial cells and reactive algal fragments of marine origin, and was subsequently confirmed to be highly reactive in incubation experiments.
Mobile portions of particulate matter can be sampled in the field using in-line filters (~0.7 µm). Particulate carbon (PC; including chlorophyll) and particulate nitrogen collected on the filters are measured using combustion-based instruments, such as the Costech CHN analyzer (Valencia, CA, USA). Particulate forms of carbon derived from buried wrack may leach dissolved organic carbon in situ, leading to heterogeneous distributions of DOC with concentrations above the seawater end-member concentration [41,52]. Additionally, particulate organic matter degradation to particulate N (PN) or particulate P (PP) can serve as a local source of nutrients to intertidal porewater (Kim et al., 2017). These studies suggest that, when coupled with information about flow patterns, spatially correlating the distribution of PC, PN, and PP may yield insight into particulate matter degradation and transport rates.
Loss-On-Ignition combustion methods can be used to quantify particulate organic carbon content in sediment samples. However, in high-quartz systems with low organic matter content and fine particulate carbon, weight differences between pre-and post-combusted samples may be minimal. Integration of in situ sampling, laboratory experiments, and porewater sampling as discussed in this section can be used to provide a clearer understanding of particulate transport, its interplay with the flow system, and its role in intertidal biogeochemical processes.

Porewater Reactivity
Reactant and chemical products measured in porewater represent the time-integrated result of on-going reactions. Direct quantification of reaction rates requires independent measurements. While reaction rates are difficult to measure in situ, flow-through cells [5] or incubation experiments can be used to quantify reaction rates in a laboratory setting [93,94].
In incubation experiments, replicate porewater samples can be collected in 12 mL Exetainers with gas-tight septa caps. Samples are filled carefully to prevent air bubbles. Samples are then submerged in water to minimize gas exchange during the incubation period. After incubating in the dark at 25 • C, each replicate sample is opened at a specified time interval across a period long enough to ensure oxygen depletion and denitrification (~2 weeks). Membrane Inlet Mass Spectrometer (MIMS; Bay Instrument, Easton, MD, USA) is used to track the consumption of O 2 (indicator of oxic respiration rate) and production of N 2 (indicator of denitrification rate).
Incubations can be conducted across a broad temperature range, but it is recommended that an incubation temperature be used that is similar to environmental conditions. This corresponds to an incubation temperature of 20-25 • C in temperate regions. Porewater temperature should be measured during sample collection, as groundwater temperature varies seasonally and spatially depending on flow path length and residence time. The incubated reaction rates can be converted to estimate in-situ reaction rates using the measured porewater temperature and the Arrhenius equation [95]: where k is the reaction rate, T is temperature in Kelvin (K), A is the pre-exponential factor, E (J/mol) is the activation energy, and R (J/mol-K) is the ideal gas constant.

Sediment Reactivity
Sediments from vibracores and slotted sand columns can be incubated to determine sediment reactivity. Kim et al. [41] sectioned two vibracores and five slotted sand columns from the FW-SW mixing zone. Sediment samples from each section were homogenized, drained, and weighed into five replicates. Samples were placed in sterilized borosilicate biological oxygen demand (BOD) bottles and filled to the rim with filtered (0.45 µm) seawater from the field site. These incubations yielded bulk (sediment and seawater) reaction rates. A control set of incubations were done simultaneously using only filtered seawater to quantify the seawater portion of reactivity. Aerobic respiration and denitrification rates in each incubation type were calculated from the change in O 2 and N 2 concentrations. Reaction rates from the control incubations (seawater only) were then subtracted from reaction rates from the incubations of sediment with filtered seawater to obtain sediment reaction rates. Vibracore incubations provide insight into in-situ sedimentary contributions to reactivity and sand column incubations are optimal for quantifying the reactivity of newer, more mobile, portions of particulate matter [41].

Long-Term Geochemical Monitoring
In-situ sensors or multi-sensor networks can complement porewater measurements by providing information into geochemical conditions between sampling events. Multi-level platinum redox sensors embedded into thin fiberglass rods (e.g., PaleoTerra, Amsterdam, Netherlands) can be installed to continuously monitor redox behavior in tidally influenced riparian aquifers (e.g., [96]) and in the intertidal mixing zone at high spatial (up to cm-scale) and temporal resolution. We installed eight redox sensor arrays, each consisting of four platinum redox sensors, across the intertidal zone on a tide-dominated beach (example array and dataset in Figure 4B,C). The sensors provided a broad yet spatially resolved depiction of geochemical conditions and showed that waves and tides can lead to extremely minor shifts in the boundary between water masses with contrasting chemical signatures. Sensors located at these jiving boundaries can exhibit large high-frequency fluctuations despite little porewater flow ( Figure 4C). Geochemical sensor deployments capable of sampling at high frequencies as demonstrated here highlight the importance of resolving porewater geochemistry across broad time scales to capture system responses to multiple ocean forcings.
Water 2021, 13, x FOR PEER REVIEW 12 of 16 lead to extremely minor shifts in the boundary between water masses with contrasting chemical signatures. Sensors located at these jiving boundaries can exhibit large high-frequency fluctuations despite little porewater flow ( Figure 4C). Geochemical sensor deployments capable of sampling at high frequencies as demonstrated here highlight the importance of resolving porewater geochemistry across broad time scales to capture system responses to multiple ocean forcings.

Conclusions
Beach aquifers are dynamic zones where mixing between seawater and fresh groundwater generates redox gradients that support biotic and abiotic chemical transformation of dissolved and particulate constituents in intertidal porewater. The extent, timing, and intensity of reactivity is a function of the degree of heterogeneity, beach slope, and the hydrologic, transport, and forcing regime. A suite of measurement approaches is required to capture the interplay between the geological, hydrological, and biogeochemical characteristics of the beach subsurface (Table 1). Hydrogeological characterization such as local and regional variations in hydraulic conductivity, grainsize, and beach topography are recommended during site selection and can be used to contextualize freshwater-saltwater mixing patterns and dissolved chemical distributions. Multiple porewater sampling techniques, including multi-level porewater samplers, repurposed gas vapor probe kits, and stainless-steel push-point samplers can be employed to obtain snapshots of the spatial distribution of chemical species. Chemical sampling and measurement methods vary greatly depending on the target parameter, but are generally limited in sample frequency because geochemical sampling is labor intensive. Sensor networks, such as distributed redox potential and dissolved oxygen sensors, capable of monitoring chemical parameters at high frequencies can greatly enhance interpretation of geochemical measurements. Future advances in the understanding of the hydrological and biogeochemical characteristics of beach aquifers are likely to benefit from more widespread use of sensor networks developed from low-cost and open-source microcontrollers, software, and sensing hardware. Ripe opportunities exist for continuous and long-term monitoring by integrating . While the two top sensors (blue and red) showed response to diurnal and spring-neap tidal fluctuations (grey), the deeper two sensors (orange and purple) were less affected by tides, oscillating slightly over a month due to spring-neap variability in tidal amplitude.

Conclusions
Beach aquifers are dynamic zones where mixing between seawater and fresh groundwater generates redox gradients that support biotic and abiotic chemical transformation of dissolved and particulate constituents in intertidal porewater. The extent, timing, and intensity of reactivity is a function of the degree of heterogeneity, beach slope, and the hydrologic, transport, and forcing regime. A suite of measurement approaches is required to capture the interplay between the geological, hydrological, and biogeochemical characteristics of the beach subsurface (Table 1). Hydrogeological characterization such as local and regional variations in hydraulic conductivity, grainsize, and beach topography are recommended during site selection and can be used to contextualize freshwater-saltwater mixing patterns and dissolved chemical distributions. Multiple porewater sampling techniques, including multi-level porewater samplers, repurposed gas vapor probe kits, and stainless-steel push-point samplers can be employed to obtain snapshots of the spatial distribution of chemical species. Chemical sampling and measurement methods vary greatly depending on the target parameter, but are generally limited in sample frequency because geochemical sampling is labor intensive. Sensor networks, such as distributed redox potential and dissolved oxygen sensors, capable of monitoring chemical parameters at high frequencies can greatly enhance interpretation of geochemical measurements. Future advances in the understanding of the hydrological and biogeochemical characteristics of beach aquifers are likely to benefit from more widespread use of sensor networks developed from low-cost and open-source microcontrollers, software, and sensing hardware. Ripe opportunities exist for continuous and long-term monitoring by integrating these technologies with wireless data communication systems.