Self-Potential as a Predictor of Seawater Intrusion in Coastal Groundwater Boreholes

Monitoring of self-potentials (SPs) in the Chalk of England has shown that a consistent electrical potential gradient exists within a coastal groundwater borehole previously affected by seawater intrusion (SI) and that this gradient is absent in boreholes further inland. Furthermore, a small but characteristic reduction in this gradient was observed several days prior to SI occurring. We present results from a combined hydrodynamic and electrodynamic model, which matches the observed phenomena for the ﬁ rst time and sheds light on the source mechanisms for the spatial and temporal distribution of SP. The model predictions are highly sensitive to the relative contribution of electrochemical exclusion and diffusion potentials, the exclusion ef ﬁ ciency , in different rock strata. Geoelectric heterogeneity, largely due to marls and hardgrounds with a relatively high exclusion ef ﬁ ciency, was the key factor in controlling the magnitude of the modeled SP gradient ahead of the saline front and its evolution prior to breakthrough. The model results suggest that, where suf ﬁ cient geoelectric heterogeneity exists, borehole SP may be used as an early warning mechanism for SI.


Introduction
Globally, groundwater provides the main source of water for human consumption and is critically important for agriculture in many countries (WWAP, 2014). Groundwater demand is particularly high in coastal areas, where population density is more than three times the global average (Small & Nicholls, 2003).
Many coastal aquifers are vulnerable to seawater intrusion (SI), with numerous incidences reported in every inhabited continent (Barlow & Reichard, 2010;Bocanegra et al., 2010;Custodio, 2010;Morgan & Werner, 2015;Shi & Jiao, 2014;Steyl & Dennis, 2010). SI risks are greatest when the water table is low relative to local sea levels, a situation exacerbated by increasing abstractions and climatic variability.
Traditional approaches for characterizing SI fall into three main categories (Werner et al., 2013): borehole hydrochemistry, monitoring of groundwater levels, and geophysical investigations. Monitoring of borehole water levels and hydrochemistry often fails to predict the timing of SI, particularly in heterogeneous aquifers, where seawater may be transported along a small number of preferential flow paths. Time-lapse resistivity and electromagnetic surveys have been used in various studies to investigate SI (e.g., Comte & Banton, 2007;Fitterman, 2014;McDonald et al., 1998), although these typically require a large footprint for the long-term installation of monitoring apparatus or repeated surveys during the predicted period of SI risk. This paper explores whether borehole measurements of self-potential (SP) may represent an alternative means of identifying a nearby saline front, as concentration gradients are known to generate SP, through the development of an electrochemical exclusion-diffusion potential V EED (e.g., Jackson, 2015;Jouniaux et al., 2009;Lanteri et al., 2009;Leinov & Jackson, 2014;Martínez-Pagán et al., 2010;Revil, 1999;Westermann-Clark & Christoforou, 1986). The relatively low cost of monitoring equipment, combined with a small footprint within a single borehole, represent significant advantages compared to traditional methods for identifying SI.
Electrochemically induced SPs have been used to track the position of an injected saline front in both field (Jougnot et al., 2015;Sandberg et al., 2002) and laboratory (Martínez-Pagán et al., 2010) experiments, with the latter suggesting that an SP signal occurs ahead of any increase in salinity. These experimental findings are supported by numerical modeling of SI in hydrocarbon reservoirs, where changes in SP occurred several tens of meters ahead of the saline front (Gulamali et al., 2011;. MacAllister et al. (2016) was the first study to demonstrate a link between tidal processes and borehole SP in a coastal aquifer. This was done by monitoring borehole SP and fluid electrical conductivity σ f in a coastal borehole in the UK Chalk subject to regular SI (the Saltdean Observation Borehole, OBH). MacAllister et al. (2018) then demonstrated that the tidal SP signal in the Saltdean OBH was dominated by the electrochemically induced component of SP across a remote saline front. MacAllister et al. (2018) also showed that the Saltdean OBH displays a consistent SP gradient ahead of the saline front and that this feature is absent in Chalk boreholes further inland. Moreover, a characteristic reduction in this gradient, or precursor, occurs several days prior to saline breakthrough in the borehole (MacAllister, 2016). Numerical modeling conducted by MacAllister (2016) and MacAllister et al. (2018) was unable to replicate the magnitude of the initial SP gradient and did not attempt to simulate the subsequent precursor. Consequently, the source mechanisms for these phenomena remain unexplained.
The aims of the present study are twofold. First, we aim to match the observed SP gradient within the borehole using a combined hydrodynamic and electrodynamic model and use this model to explain the key parameters that control this phenomenon. Second, we will use the combined numerical model to match the observed evolution of SP prior to breakthrough and investigate possible causes of the precursor, in order to understand the broader applicability of SP as a predictor of SI.

SP Source Mechanisms
Numerous mechanisms can generate SPs, including electrokinetic, electrochemical, thermoelectric, piezoelectric, and redox processes (e.g., Jouniaux et al., 2009). At conditions close to thermodynamic equilibrium (Revil & Linde, 2011), the contributions of each source mechanism i can be related to total current flow j via the following generalized equation (Jackson, 2015;Jackson, Gulamali, et al., 2012;Saunders et al., 2008;Sill, 1983): where σ T is the effective conductivity of the host material including matrix and pore space constituents (S/m), V is the electrical potential (V), U i is the potential of the source mechanism (e.g., pressure and concentration), and L i is the associated cross-coupling term. In the absence of SP source mechanisms, equation (1) is equivalent to Ohm's law. The cross-coupling term is often expressed in terms of a coupling coefficient C i (Jackson, Gulamali, et al., 2012): In coastal aquifers, V EED and V EK are likely to be the most significant components of SP (MacAllister et al., , 2018, because of respective variations in pressure induced by tidal fluctuations and in groundwater salinity, due to the presence of seawater at depth below comparatively fresh groundwater near the water table.

Electrokinetic Potential V EK
An excess of charge typically occurs at mineral-water interfaces, which is balanced by an adjoining layer of opposing charge (counterions) within the fluid (Hunter, 1981). This arrangement is often referred to as the electric double layer. The fluid layer closest to the mineral surface (the Stern layer) is characterized by strongly sorbed counterions; the layer furthest from the mineral surface is known as the diffuse layer (or Gouy-Chapman layer) and is characterized by a lower density of charge, which is mobile and can be transported by flow .
In chalk saturated with groundwater or seawater, the surface charge is typically negative (Jackson, Butler, & Vinogradov, 2012;MacAllister, 2016) and groundwater flow from nonhydrostatic pressure gradients transports positively charged ions along the diffuse layer, leading to an electrokinetic current (Jackson, 2015) ( Figure 1a). A conduction current arises to maintain overall electroneutrality, and V EK represents the electrical potential required to sustain this current (Jackson, Gulamali, et al., 2012). In the absence of other current sources, equations (1) and (2) give (Jackson, 2015) j where P n is nonhydrostatic pressure (Pa) and C EK is the electrokinetic coupling coefficient (V/Pa), which in the absence of current flow is defined by 10.1029/2018WR022972

Water Resources Research
The magnitude of C EK decreases with concentration up to the concentration of seawater; in highly saline groundwater, counterions are closer to the mineral surface on average and a lower charge density is transported under a given (nonhydrostatic) pressure gradient (Jaafar et al., 2009;Vinogradov et al., 2010).

Exclusion-Diffusion Potential V EED
In the absence of other current sources, equations (1) and (2) can be used to give V EED (e.g., Jackson, 2015;Revil & Linde, 2006), which comprises a diffusion (or liquid junction) and an exclusion (or membrane) component: where C f is the molar fluid concentration (M) and C EED is the electrochemical coupling coefficient (V). In the absence of current flow, this may be rewritten as (Leinov & Jackson, 2014;Revil, 1999) where T + is the macroscopic Hittorf number for cations (dimensionless), k B is the Boltzmann constant (J/K), T is temperature (K), and e is the charge on an electron (C).
In the case of a salinity gradient in an uncharged porous medium, where the solute is dominated by sodium and chloride ions, charge separation is induced by the greater mobility of chloride, giving rise to an electrochemical diffusion potential (V ED ; e.g., Jackson, 2015;Revil, 1999;Figure 1b). In this scenario, the SP signal is effectively a diffusion potential and T + is equivalent to the microscopic Hittorf number for sodium, t Na . Ignoring the weak concentration dependency of t Na , and assuming a constant temperature, gives the following expression (Jackson, 2015;Leinov & Jackson, 2014): The following relation can be used to represent the slight changes in t Na with ionic strength (Gulamali et al., 2011;Jackson, Gulamali, et al., 2012;Leinov & Jackson, 2014): If the salinity gradient lies within a negatively charged porous medium and the electrical double layer is thick relative to the pore-throat radius, chloride ions are excluded from the pore space, leading to an exclusion potential V EE (Figure 1c). In a perfect membrane, T + = 1, giving (Leinov & Jackson, 2014;Westermann-Clark & Christoforou, 1986): The relative importance of the exclusion and diffusion components is defined by the dimensionless exclusion efficiency η (

Site Characterization
A downhole array of SP electrodes was installed in May 2013 in the Saltdean OBH near Brighton on the south coast of England ( Figure 2). The Saltdean OBH, which lies in a dry valley approximately 1.8 km from the coast, was formerly used as an abstraction borehole, but was abandoned in 1936 because of repeated saline intrusion (Jones & Robins, 1999). Borehole logs from the area suggest that an adit with a diameter of 1.4 m intersects the Saltdean OBH at À2 m above Ordnance Datum (AOD) and extends 32 m to the northwest, connecting the Saltdean OBH to two other groundwater boreholes nearby (BGS, 2018). The Saltdean OBH now acts as a monitoring borehole for the Balsdean Pumping Station (PS) further inland, which provides drinking water to the eastern part of Brighton.

Geology
The site lies within the South Downs, an area dominated by chalk strata that form part of the wider Chalk Group. The Saltdean OBH extends to 60 m below ground level, intersecting a series of chalk, marl, and hardground layers within the Seaford and Lewes Nodular Chalk Members of the White Chalk Subgroup ( Figure 2). The strata dip gently toward the coast (Jones & Robins, 1999) at an angle of approximately 5°, based on the inferred depth of the Seaford-Lewes Nodular Chalk contact in the Saltdean OBH (AMEC, 2012) and at the PS (MWH, 2006).

Hydrogeology 3.2.1. Groundwater Flow
The majority of groundwater flow in the Chalk occurs within fractures (Allen et al., 1997;Bloomfield, 1996;MacDonald & Allen, 2001), where the hydraulic conductivity is typically several orders of magnitude higher than that of the surrounding matrix Jones & Robins, 1999). Fractures may exist as faults, joints, or bedding planes, although the latter are considered to be the most pervasive and the most significant conduits for flow (Bloomfield, 1996).
Numerous, laterally extensive marl seams have been logged within the Chalk, which act as barriers to flow due to their relatively low hydraulic conductivity (Gallagher et al., 2012;Molyneux, 2012;Zaidman et al., 1999). Hardgrounds are often considered flow conduits, due to their brittle nature and enhanced fracturing (Schurch & Buckley, 2002;Soley et al., 2012), although they may also be barriers to flow when fracturing is minimal (Jones & Robins, 1999).
In general, the hydraulic conductivity of the Chalk decreases with depth and the most significant flow horizons typically occur within 50 m of the water table Jones & Robins, 1999;Williams et al., 2006). Laterally, the Chalk is significantly more permeable in valleys than in interfluves (Jones & Robins, 1999;

Water Resources Research
MacDonald & Allen, 2001). Within the South Downs, the Seaford Chalk is considered to be the major waterbearing unit (Jones & Robins, 1999).
According to Jones and Robins (1999), logging of the Saltdean OBH has indicated a hydraulically significant fracture at its base, which constitutes the main conduit for SI. This was later supported by MacAllister et al. (2018), who showed through a series of σ f profiles that saline water enters near the base of the borehole and exits via the adit near the top of the water column.
Recharge estimates vary across the South Downs, although in the Brighton Block, where the Saltdean OBH is located, Jones and Robins (1999) report a mean value of around 480 mm/yr. The majority of recharge occurs during the winter months, with limited (Ireson et al., 2009;Wellings, 1984) or no recharge (Limbrick, 2002) between May and September. Despite the minimal flow occurring in the unsaturated zone during this period, the narrow pore-throat diameters in the Chalk matrix imply that it remains almost entirely saturated within 30 m of the water table (Price, 1987), although dewatering of the fractures occurs, typically accounting for less than 1% of the total rock volume Mathias, 2005;Price, 1987).

Aquifer Properties
The median storage coefficient in the South Downs Chalk, as reported by MacDonald and Allen (2001) is 0.0022 and the mean porosity of the White Chalk Subgroup in the south of England is 39% (Bloomfield et al., 1995).
Based on observed drawdown at the Balsdean PS, MWH (2006) inferred a hydraulic conductivity K of 900 m/day for the surrounding aquifer, significantly higher than the median value of 9 m/day for the South Downs as a whole (based on the median transmissivity of 880 m 2 /day reported by MacDonald and Allen (2001) and an assumed aquifer thickness of 50 m). We estimate here a hydraulic conductivity of 250 m/day, based on a mean lag τof 0.152 days between tidal peaks observed at the coast (in Newhaven, 7 km to the southeast; British Oceanographic Data Centre, 2015) and using (Jacob, 1950): where t 0 is the semidiurnal period of tidal fluctuation (0.518 days; MacAllister et al., 2016), S is storativity (0.0022), x is distance from the coast (1,800 m), and b is aquifer thickness (50 m).

Electrical Parameters
Rearranging equations (5)-(7), (9), and (10), assuming no current flow, gives the following expression for C EED in a fully saturated medium: MacAllister (2016) measured η in Seaford (0.01-0.12) and Lewes Nodular Chalk (0.02-0.06) cores saturated with local groundwater and seawater, although no direct measurements of η are available for marls and hardgrounds, which lie beneath the Saltdean study area. These values can be inferred based on observed changes in σ f and SP with depth in the Victoria Gardens borehole, 7 km west of Saltdean (see supporting information Figure S1).
Values of σ f in mS/cm are converted to total dissolved solids in mg/L (TDS), using the following relation (Walton, 1989): Measurements of TDS for seawater (34,113 mg/L) and groundwater (376 mg/L) are compared to calculated values of C f (0.673 for seawater; 0.00851 for groundwater) reported by MacAllister (2016) giving Cf ≈2:12Â10 À5 TDS: Applying the variations in C f and SP to equations (7) and (9) Christoforou, 1986), which suggest that lithologies with relatively narrow pore-throat diameters, such as shale (Nelson, 2009), marls, and hardgrounds (Fay-Gomord et al., 2016) typically have higher values of η. Indeed, the value of η for shale (0.24) reported by MacAllister (2016) is similar to the values derived for marls and hardgrounds. Vinogradov et al. (2010) showed that |ln(C EK )|varies linearly with ln(C f ) in fully saturated core samples below seawater salinity. Using measurements of C EK for local groundwater (À5.12 × 10 À6 V/Pa) and seawater (À1.76 × 10 À7 V/Pa) reported by MacAllister (2016), we obtain The electrical conductivity of the Chalk aquifer σ T can be estimated using Archie's law (Telford et al., 1990): where m is the cementation exponent, ϕ is porosity, andS w is saturation. MacAllister (2016) derived values for m of 2.1-2.6 from core samples of Seaford and Lewes Nodular Chalk.
The presence of clay minerals increases electrical conductivity; this can be accounted for using the empirical Waxman-Smits equation (Darling, 2005), originally developed for shaly sands: where ρ m is the mineral grain density (Fay-Gomord et al., 2016, report values of 2.7 g/cm for marl), CECis the cation exchange capacity (averaging 9.5 mEq/100 g for a typical marl; Cornell & Aksoyoglu, 1991) and where T C is temperature in degrees Celsius, which remains close to 11°C in the Saltdean OBH throughout the monitoring period.

Monitoring Apparatus and Data Processing
An array of 13 nonpolarizing Silvion Ag/AgCl WE300 electrodes was installed in the Saltdean OBH from May 2013 to February 2015 ( Figure 2). The shallowest of these was installed at À2.8 mAOD and used as a reference electrode throughout, with the remaining electrodes spaced at 2-m intervals below it. The deepest electrode was installed separately from the others at À26.8 mAOD, along with an AquaTroll 200 probe which recorded temperature T, σ f , and pressure P; several vertical profiles of SP and σ f were collected over the monitoring period using this traveling electrode. Two further AquaTroll probes were installed at À4.8 mAOD and À10.8 mAOD respectively to record vertical changes in T, σ f , and P. The monitoring apparatus was connected to a Campbell Scientific CR3000 data logger, which collected data at 5-min intervals. Further details on the equipment used and data collection methods are given by MacAllister et al. (2016).
The raw SP data include semidiurnal oscillations consistent with the M2 component of oceanic tides, which dominate the signal . A first-order Savitsky-Golay (SG) moving average filter, with a sampling window of 2.2 days, was applied to the data to assist in identifying longer-term trends in the SP data set. The SG filter has previously been used for SP analysis (e.g., MacAllister, 2016; Maineult et al., 2008) and has the advantages of maintaining the shape and amplitude of lower frequency oscillations, without introducing phase delay (Savitsky & Golay, 1964).
In addition to the long-term monitoring conducted at the Saltdean OBH, vertical SP profiles were collected from the Balsdean OBH around 1 km further inland ( Figure 2a) and from borehole PL10B, located at Trumpletts Farm in Berkshire, more than 60 km from the coast (MacAllister et al., 2018). All three boreholes intersect the Seaford and Lewes Nodular Chalk, allowing a comparison of SP profiles in a similar geological setting, at varying distances from the coast.

Summary of Field Data
The data recorded by the deepest T/σ f /P probe show that saline water entered the Saltdean OBH in late summer/early autumn of both 2013 and 2014 (Figure 3a). These events occurred following prolonged periods of low water levels in the borehole; intrusion ceased at the commencement of the winter recharge period, as shown by the higher water levels from November onward.
The potential gradient across the water column (based on vertical profiles taken by the traveling electrode) remained close to 0.2 mV/m throughout the monitoring period (Figure 3b). The respective SP gradients in the Balsdean OBH and Trumpletts borehole, which both lie further inland, are minimal by comparison (Figure 4). The σ f in the Saltdean OBH was approximately 0.87 mS/cm throughout the water column (~610 mg/L, from equation (13)).
SP across the water column decreases by 0.2-0.3 mV several days prior to intrusion ( Figure 5); this precursor is most apparent in the 2013 data. These precursors also appear in the 2014 data, although the occurrence of what appear to be two minor intrusion events prior to more obvious saline breakthrough on 12 September leads to a more complex pattern of both salinity and SP.
As well as the intrusion of saline water when water levels are low, a smaller increase in salinity occurs in conjunction with major recharge events, indicating the presence of elevated salinity above the water table ( Figure 6). This influx of more saline water is preceded by a sharp reduction in σ f of short duration. This may reflect distinct contributions from (i) rapid  bypass flow in the fractures, delivering relatively fresh recharge to the borehole, and (ii) piston flow, yielding more mineralized recharge, as a result of prolonged interactions between percolating rainwater and the rock matrix (see Jones & Robins, 1999).

Model Description
An electrodynamic model was produced to investigate the static SP gradient ( Figure 5) and variations in SP prior to intrusion (Figure 6) in the Saltdean OBH. The electrodynamic model relies on distributions of pressure and salinity (see equations (3)- (6)), which are provided by an accompanying hydrodynamic model of the coastal Chalk aquifer.

Hydrodynamic Model
Hydrodynamic simulations were conducted using SUTRA3D (Voss & Provost, 2002), one of the most widely used models for simulating density-driven flow and transport (Werner et al., 2013). Three model domains were used (Figure 7) to simulate pressure and salinity variations at the coast, while giving a detailed 3-D representation of the geology around the borehole, consistent with the stratigraphy and fracture zones shown in Figure 2.
The model was run initially with long-term average conditions until a steady state distribution of salinity was achieved. This steady state model provided the initial conditions for a transient model, which includes tidal oscillations in pressure at the coast and invokes saline intrusion by turning off recharge after 7 days of the simulation.
Parameterization of the model is consistent with the site characterization in sections 3.1 and 3.2. A detailed description of the parameters applied to the hydrodynamic model, along with time stepping and solver parameters, is given in supporting information Text S2.
Pressure and salinity distributions produced by the hydrodynamic model were used to calculate SP within the electrodynamic model described in the following section.

Electrodynamic Model 6.2.1. Initial Data Processing
The electrodynamic model was written in MATLAB, using the controlled volume finite difference method, and is based on the approach of Gulamali et al. (2011), Ijioma and Jackson (2014), and Jackson, Gulamali,   (5) for V EK and V EED , assuming no external current sources or sinks, such that at the model boundaries. Modeled SP is obtained by adding the contributions of V EK and V EED .
The model begins by mapping pressures and concentrations from the 2-D local model onto the logically rectangular mesh shown in Figure 8 and assuming no variation in the y direction. It then incorporates data for the region covered by the 3-D hydrodynamic model. The top 40 elements in Figure 8 lie above the extent of the 2-D local and 3-D hydrodynamic models; parameterization of this region is described in section 6.2.2.

Water Resources Research
Hydrostatic pressure P h was inferred by taking the nodal pressure values at the base of the fortieth element P 40 , calculating the arithmetic mean of fluid density ρ f between each element below them and invoking where g is the acceleration due to gravity (9.81 m/s 2 ), Δzis the change in elevation between elements, andρ f is linearly interpolated from the density of fresh water (1,000 kg/m 3 for 0 mg/L) and sea water (1,025 kg/m 3 for 35,000 mg/L). The nonhydrostatic pressure component P n was obtained by subtracting P h from the total pressure at the center of each element.
Spuriously low or high salinities produced by numerical oscillations in the SUTRA model were removed by limiting minimum and maximum concentrations to 609 mg/L and 35,000 mg/L, respectively. To mimic the rate of breakthrough shown in Figure 6 and avoid minor fluctuations in borehole salinity prior to this, all salinity values below 700 mg/L were set equal to 609 mg/L.
As the front within the 2-D local SUTRA model was more dispersed than that within the 3-D model (see supporting information Text S2), the advancing salinity within the upper fracture zone (shown in Figure 2) led to spurious elevated concentrations within the overlying strata at later time steps. To avoid this, the salinity distribution for the start of the transient model was maintained throughout for elements, which lie less than 2.5 m above the upper fracture zone.
Modeled porosity was 0.39 throughout the model domain, except in the borehole and adit, where a value of 1 was applied.

Unsaturated Zone
Where the water table lies above the top of the SUTRA local model, its elevation was inferred from equation (20) and values of P 40 , assuming hydrostatic conditions and a constant fluid density of 1,000.435 kg/m 3 (consistent with a fluid concentration of 609 mg/L). Saturation was reduced from 1 to 0.995 over the two elements immediately above the water table to reflect dewatering of fractures and concentration was increased linearly to 2,200 mg/L over a further two elements above this, reflecting the presence of elevated salinity in the unsaturated zone ( Figure 6). Constant values of saturation (0.995) and concentration (2,200 mg/L) were maintained above this. 6.2.3. Electrical Parameters TDS concentrations were converted to σ f and C f using equations (13) and (14).
The bulk electrical conductivity of the aquifer was simulated using equation (16), with a value of 2.5 for the cementation exponent m below the unnamed marl and a value of 1.5 above this, reflecting an increased degree of fracturing at shallow depths and, hence, greater connection between pores (Glover, 2009;Roubinet et al., 2018). For marl seams, equation (17) was applied, reflecting the contribution of surface conductance in clay minerals.

Water Resources Research
The values of the coupling coefficients C EED and C EK were obtained from equations (12) and (15); the values of η applied to equation (15) are shown in Table 1.

Model Results and Discussion
Steady state SP within the borehole is compared to the observed profile in Figure 9a. The modeled profile is dominated by V EED and the contribution of V EK is negligible by comparison, consistent with the previous study of MacAllister et al. (2018).
The sensitivity of the SP gradient to multiple input parameters is shown in Figure 9b. Increasing recharge by 25% produces a deeper saline front and reduces the SP gradient substantially. Removing variations in η has a more pronounced effect on SP and in the case of a constant value of η (0.04) throughout the model domain, the SP gradient ahead of the front disappears.
The removal of elevated salinity (and hence electrical conductivity) in the unsaturated zone also lowers the modeled gradient. The elevated electrical conductivity of clay minerals in the marl bands does not greatly affect the magnitude of the SP gradient, although it does appear to be responsible for the deviation in observed SP adjacent to the Shoreham Marl, just below À20 mAOD.
The results suggest that the observed SP gradient relies both on a nearby saline front and local variations in η. Given the extremely low permeability of the chalk matrix, the presence of a nearby saline front may reflect historic SI over a period of many years. While these historic events would suggest a high risk of future intrusion, it is difficult to infer the timing of future events from a single SP profile. To this end, the evolution of the SP gradient could provide a better indication of imminent SI.
Field observations of the precursor are compared to the results of the best fit model and sensitivity analyses in Figure 10. This includes (on the left-hand side) raw model outputs and field data with a short-period SG filter applied to remove high-frequency noise; and (on the right-hand side) field and model data with a 2.2-day SG filter applied to remove the effects of tidal oscillations. The results from the best fit model   (Figures 10a and 10b). The smaller magnitude of σ f peaks following breakthrough in the model highlights the propensity for numerical models to produce highly dispersed fronts (Konikow, 2011).
Sensitivity analysis in the model provides an important insight into the likely mechanisms that drive the precursor. The first analysis investigates the contribution of the inferred lower fracture zone in Figure 2 to the evolution of SP. This is achieved by treating the lower fracture zone as a region of low permeability, unfractured chalk, and modeling the initial salinity distribution in this layer throughout the transient model (while allowing the saline front to evolve freely throughout the rest of the model domain). As shown in Figures 10e and 10f, this completely removes any precursor to intrusion.
However, movement of salinity beneath the borehole is not sufficient by itself to produce a precursor. Applying η = 0.01 on both sides of the lower fracture zone also removes the modeled precursor ( Figures 10g and 10h). The precursor therefore reflects the movement of salinity through a heterogeneous geoelectric environment.
A precursor is also obtained by applying the variation in η values (from 0.01 below to 0.03 above) across the upper fracture zone (Figures 10i and 10j). In this scenario, a static front was maintained in the lower fracture zone, so movement of saline water within the upper fracture zone was the source of the modeled precursor. However, although the smoothed precursor is similar to that observed in the borehole, tidal oscillations are notably smaller in magnitude than those observed in the field and in the best fit model. Figure 11 shows the saline front near the borehole and changes in SP (relative to the steady state model) immediately prior to SI. In each case, the saline front within the upper fracture zone lies <1 m laterally from the base of the borehole, on its seaward side. Figures 11a and 11b, which represent the best fit model, show how the SP gradient in the borehole is affected ahead of the saline front within the upper fracture zone. When a static front is applied to the lower fracture zone only, changes to the vertical SP gradient occur directly above the saline front within the upper fracture zone, with minimal changes ahead of it (Figures 11c and 11d), thereby giving no clear precursor (as shown in Figure 10e).
When the best fit model is modified by applying η = 0.01 on either side of the lower fracture zone, changes in the vertical SP gradient lag behind the advancing saline front in the upper fracture zone (Figures 11e and 11f). Finally, by applying a change in η across the upper instead of lower fracture zone and implementing a static front in the lower fracture zone (Figures 11g and 11h), the changes in SP are similar to those seen in the field and in the best fit model.
The results from the transient model suggest that the precursor is driven by small spatial changes in η across one or more conduits for saline intrusion. As shown in Figures 11e and 11f, the SP source conduit may intersect the borehole, or it may lie immediately beneath it (Figures 11a and 11b).
In the latter scenario, the precursor does not directly reflect the movement of seawater that will enter the borehole. However, geophysical logging and hydraulic testing (e.g., Butler et al., 2009;Gallagher et al., 2012;Jones & Robins, 1999) suggests that there are typically numerous hydraulically significant fractures within the upper 100 m of the Chalk. SI in the Chalk (and other fractured coastal aquifers) may therefore be characterized by saline water filling progressively shallower fracture zones. Thus, the arrival of higher salinity groundwater immediately beneath the borehole would be a strong but indirect indicator of imminent breakthrough into the borehole itself.
At present, a precursor of this nature has been observed only at a single site, making it difficult to assess the broader applicability of SP monitoring to predict SI. Furthermore, there are very few data on η, so the prevalence and magnitude of small-scale variations in η around a typical borehole is unknown. Other parameters, such as clay content and permeability, are highly variable in the Chalk and most other common aquifer lithologies over distances of only a few millimeters. Consequently, there are grounds for optimism that variations in η are also common and so SP precursors to SI may indeed be widespread. Further data collection is required to identify the presence (or otherwise) of precursors in other boreholes.
Based on the sensitivity of the model and our experience of SP monitoring, it is important that nonpolarizing Ag/AgCl or Pb/PbCl electrodes (see Perrier et al., 1997) are used to collect long-term SP data sets, to reduce problems with electrode drift that may complicate subsequent analysis. MacAllister (2016) and MacAllister et al. (2016) provide additional detail on the requirements for monitoring SP in a coastal aquifer.

Conclusions
Data collected from a groundwater observation borehole in a coastal aquifer near the south coast of England show a strong and consistent SP gradient ahead of a saline front. There is a characteristic reduction in this gradient several days prior to SI in the borehole. If present in other boreholes, this precursor could act as a warning system for groundwater users and may allow sufficient time to mitigate the risk of SI occurring.
The combined hydrodynamic and electrodynamic model developed in this study closely matches the phenomena observed in the field, using plausible values for the numerous input parameters. Through sensitivity analysis, it is clear that the SP gradient requires a nearby saline front and local variations in the exclusion efficiency, a parameter which will require further measurements to understand its distribution and variability in different aquifer types. The precursor signal modeled here also relies on local variations in exclusion efficiency on either side of a fracture zone transmitting seawater inland. The fracture zone may either intersect the borehole, or be located immediately below it, although only the latter scenario replicates the observed magnitude of tidal SP variations. Further analysis should focus on tidal SP variations and how they evolve prior to SI.
The results are a promising first step in demonstrating the possible use of SP as a predictor of SI, although data from additional sites are required to demonstrate its widespread applicability.