Mapping soil layers using electrical resistivity tomography and validation: Sandbox experiments

Journal of Hydrology 575 (2019) 523–536

Contents lists available at ScienceDirect

Journal of Hydrology journal homepage: www.elsevier.com/locate/jhydrol

Research papers

Mapping soil layers using electrical resistivity tomography and validation: Sandbox experiments

T

Dong Xua,d, Ronglin Sunb, Tian-Chyi Jim Yehc,d, Yu-Li Wangd, Moe Momayeze, Yonghong Haod, ⁎ Cheng-Haw Leef, Xiangyun Hua, a

Hubei Subsurface Multiscale Imaging Laboratory, Institute of Geophysics and Geomatics, China University of Geosciences, Wuhan 430074, China School of Environmental Studies, China University of Geosciences, Wuhan 430074, China c Tianjin Key Laboratory of Water Environment and Resources, Tianjin Normal University, Tianjin 300387, China d Department of Hydrology and Atmospheric Science, University of Arizona, Tucson, AZ, USA e Department of Mining and Geological Engineering, University of Arizona, AZ, USA f Department of Resources Engineering, National Cheng Kung University, Tainan, Taiwan, ROC b

A R T I C LE I N FO

A B S T R A C T

This manuscript was handled by Corrado Corradini, Editor-in-Chief, with the assistance of Renato Morbidelli, Associate Editor

Knowledge of the geologic structure at a field site is a useful piece of information for hydrologic modeling since it can serve as more site-specific prior information about the hydraulic parameter patterns at the site than generic spatial statistics. Widely accepted electrical resistivity tomography (ERT) survey for mapping subsurface anomalies could be a viable tool for acquiring this information. However, their ability to delineate geologic structures has not been thoroughly investigated. In this study, two-dimensional ERT numerical experiments were first conducted to study the effects of boundary conditions on the dipole-dipole and pole-pole array configurations. An ERT setup subsequently was implemented in a sandbox consisting of complex layers of different sands. A continuous copper wire was installed along the sides of the sandbox to impose potential boundaries. Using data collected with the pole-pole array in the sandbox under different degrees of drainage and the Successive Linear Estimator (SLE) algorithm, we show that ERT yields electrical conductivity estimates of complex layer structures with small uncertainties. In addition, using SLE with physically meaningful correlation scales as prior information can lead to an electrical conductivity field that is consistent with visually observed layer structures. The correlation scale concept also was demonstrated to provide guidance to the design of the electrode spacing in the surveys. Moreover, the estimated field was validated by predicting electrical potential fields from two independent ERT surveys using electrodes at different locations. Results of this study suggest that the combination of ERT and SLE is a viable geophysical survey tool for mapping geologic layer structures. Research Significance: This study develops a method to implement prescribed potential boundaries for enhancing the pole-pole ERT survey, illustrates the importance of correlation scales and develops an approach for validating ERT results.

Keywords: Electrical resistivity tomography Successive linear estimator Sandbox experiment Fine-resolution geophysics Uncertainty Validation

1. Introduction The increasing interest in predictions of groundwater flow and solute transport in the subsurface has prompted a large number of studies on an accurate characterization of subsurface heterogeneity. Although hydraulic tomography (HT) has been developed and proven as an effective technique to hydraulic heterogeneity in the subsurface over past years (Yeh and Liu, 2000; Bohling et al., 2002; Zhu and Yeh, 2005; Illman et al., 2007, 2009; Cardiff et al., 2009, 2012; Mao et al 2013), prior large-scale geologic conceptual models built upon borehole logs still play an important role when the number of pumping tests and

⁎

drawdown measuring points is sparse. Besides, these conceptual models improve HT analysis as well as hydrogeological forward modeling efforts (Tso et al. 2016; Zha et al., 2017; Zhao and Illman, 2017). Nonetheless, borehole information often is sparse and it provides discrete and local-scale descriptions of geology. To improve the large-scale geologic conceptual model, inexpensive electrical resistivity surveys (e.g., Hermans et al., 2015, 2017) may be the solution. Surface resistivity surveys have been widely used in the field. They measure the difference in voltage generated by injection of electric current between two electrodes implanted at the ground surface. Then, the electric potentials induced by many combinations of transmitting

Corresponding author. E-mail address: [email protected] (X. Hu).

https://doi.org/10.1016/j.jhydrol.2019.05.036 Received 14 March 2019; Received in revised form 8 May 2019; Accepted 10 May 2019 Available online 11 May 2019 0022-1694/ © 2019 Elsevier B.V. All rights reserved.

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

the electrical conductivity distribution from the observed potentials through the governing current and electric potential equations (Cockett et al., 2015; Gunther et al., 2006; Loke et al., 2013, 2014; Pidlisecky et al., 2007; Blaschek et al., 2008; Lelievre et al., 2009). For example, using a smoothness constraint (e.g., a discrete second difference matrix (Constable et al., 1987)), Occam’s algorithm aims to find the most featureless models that are consistent with the measurements (Gouveia and Scales, 1997). While such a smoothness constraint resolves the underdetermined nature of the inverse problem, the solution to the problem remains highly uncertain. That is, the smoothness constraint does not necessarily lead to an estimate which is the best-unbiased and has the minimum variance in comparison with the true field. To overcome this shortcoming, recent studies emphasized inversions based on regularization approach (e.g., covariance-based regularization) that incorporates prior information of the geologic structure and some direct measurements of parameters (Robinson et al., 2013; Johnson et al., 2007; Hermans et al., 2012, 2016; Bouchedda et al., 2017; Jordi et al., 2018). To incorporate prior information of the geologic structure and direct measurements of parameters, Yeh et al. (2002) introduced a geostatistics-based successive linear estimator (SLE) approach to the ERT inverse modeling. This iterative geostatistical technique incorporates the spatial statistics, point measurements of variables as prior information to derive the conditional mean estimates and residual covariance, which describes the uncertainty associated with the estimates. Liu and Yeh (2004) demonstrated the utility of the SLE method for monitoring 3-D moisture content distributions using numerical experiments. Solving mathematical ERT inverse problems also encounters the issues related to the unknown boundary conditions. In laboratory experiments, insulated boundary conditions are known because of isolation of the sandbox. However, a solution domain surrounded by insulated boundaries is a mathematically ill-posed forward problem. Constant potential boundaries or arbitrarily assigning a constant potential in the domain has often been used to avoid this mathematical issue for the laboratory experiment (Sentenac et al., 2015). Albeit precise interpretation of ERT as discussed above may be difficult, ERT remains as an attractive tool for mapping different-scale geologic structures in terms of different electrode spacing. This geologic structure information is beneficial to the characterization of geologic hydraulic properties since it serves as a site-specific prior geologic structure information for hydrologic inverse modeling efforts or hydraulic tomography as demonstrated by Tso et al. (2016), Zha et al. (2017). In spite of the fact that ERT has shown some promising results, such as Hermans and Irving (2017), its ability to map geologic layers at different scales remains to be further verified, in particular, in controlled laboratory-scale experiments. The purpose of this paper, therefore, is to demonstrate the ERT inversion using SLE for mapping soil layers and to validate its estimates using sandbox experiments. For this purpose, a design that facilitates constant electrical potential boundaries was implemented in the sandbox to allow fast/easy pole-pole array measurements. SLE was then used to interpret the observations to demonstrate the effects of electrical conductivity contrasts in the sandbox at different degrees of saturation. We then elucidate the importance of spatial correlation scales as well as measurement electrode spacing on ERT. The estimated electrical conductivity field was further validated by predicting the response during ERT surveys using different source electrode arrays while quantifying the uncertainty of parameter estimates. We believe these are new and important contributions to ERT surveys and analyses.

and receiving electrodes were measured to interpret the distribution of subsurface resistive or conductive zones (Yeh et al., 2002; Parsekian et al., 2015), which are then used to construct the geologic conceptual model. Similar to the surface surveys, Shima and Sakayama (1987) proposed “Electrical resistivity tomography (ERT)” based on crossborehole and borehole-to-surface resistivity surveys for reconstructing a resistivity or electrical conductivity cross-section. Sasaki (1992) and Yeh et al. (2002) showed the resolution of subsurface images from the tomography technique is much clear than that of the surface resistivity surveys. Notice the resolution is a measure of the level of detail that can be seen using an ERT survey. Over the past decades, ERT has increasingly been recognized as a useful geophysical method to investigate the subsurface resistivity anomalies (Lapenna et al., 2005; Miller et al., 2008). It also has been widely applied to environmental and hydrogeological studies (Binley et al., 2010; Binley et al., 2015; Hayley et al., 2009; Revil et al., 2012; Singha et al., 2015; Slater, 2007) because ERT surveys are sensitive to the conductivity contrast due to differences in porosity, moisture content (Xu et al., 2016), geologic layering (Loke et al., 2013), or presence and migration of contaminants (Deng et al., 2017; Mao et al., 2016) in the subsurface. Further, ERT is easy-to-implement and relatively inexpensive (Loke et al., 2013). Deciphering the data from an ERT survey, however, often involves large uncertainty. First, the volume of the geologic medium to be surveyed is large, and the density of ERT array is generally low. Additionally, the electrical potential response at a location to a transmitting current at a different location represents the cumulative response from all parts of the medium. While the tomographic survey aims to collect non-redundant responses using sparse measurement arrays to make the inverse problem better-defined (Yeh et al., 2015a,b), the uncertainty associated with the interpretation remains large. For hydrogeological applications, translating electrical resistivity to hydrologic properties (such as hydraulic conductivity, moisture content, and concentration) adds another uncertainty. First of all, the bulk electrical conductivity of geologic media is affected by many factors (such as soil texture, moisture content, porosity, temperature, mineralogy, and geochemistry, see Slater, 2007). Although some constitutive relationships (i.e., petrophysical model), such as Archie’s law (Archie, 1942) for the relationship between porosity and resistivity and log–log relationships between hydraulic conductivity and resistivity (Purvance and Andricevic, 2000a,b), are often used to relate the resistivity to these factors, these constitutive relationships are empirical, scale-dependent (Yeh et al., 2002; Slater, 2007). Moreover, the parameters (e.g., the cementation factor in Archie’s law) in the relationships are not uniform for a given site but spatially varying (Yeh et al., 2002) due to the local geologic and geochemical environment (Purvance and Andricevic, 2000a,b; Slater, 2007). To avoid the aforementioned uncertainty and assess the usefulness of ERT applications to hydrogeological problems, over the past few decades, many researchers have conducted controlled laboratory experiments. For example, Chen et al. (2016) used the ERT to estimate groundwater flow velocity and monitor the advection of a conductive plume in a sandbox using Archie’s law. Similarly, Pollock and Cirpka (2012) conducted a salt tracer experiment in a sandbox to determine the hydraulic conductivity distribution based on the fully coupled hydrogeophysical inversion (Pollock and Cirpka, 2008, 2010) of timelapse ERT and hydraulic head measurements. Slater et al. (2002) conducted a 3D cross-borehole ERT in a large, heterogeneous, experimental tank to study solute transport processes. They converted the time series of electrical images into estimated fluid conductivity breakthrough curves which were in good agreements with direct measures of fluid conductivity in the sampling ports. Other similar studies include Mao et al. (2015), Clement and Moreau (2016), Slater et al. (2000), and Binley et al. (1996). As for interpreting the results of ERT survey, the least squares optimization method with regularization has been widely used to estimate

2. Methodology 2.1. Governing equations The electric potential field in geologic media induced by a direct current transmission generally can be described by: 524

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

∇ ·(σ ∇u) + Iδ (x − xs ) δ (y − ys ) δ (z − z s ) = 0

The coefficient matrix ωT is determined by solving the following equation:

(1)

where σ = σ (x , y, z ) is a given electrical conductivity [S m−1] distribution in a Cartesian system of coordinates (x , y, z ) , u = u (x , y, z ) the electrostatic potential [V], I the source current [A] at a location (xs , ys , z s ) and δ the Dirac delta function. In order to solve Eq. (1), boundary conditions must be specified and they are

u|Γ1 = u∗

(r ) (r ) [Rdd + θ (r ) diag(Rdd )] ωf(r ) = Rdf(r )

is the auto-covariance of the potential measurement, and where Rdf(r ) is the cross-covariance between the conductivity and the potential measurements. The parameter θ (r ) is a dynamic stability multiplier, and (r ) ) is a stability matrix, which is a diagonal matrix, which is the diag(Rdd (r ) (r ) . The auto-covariance Rdd and same as the diagonal elements ofRdd (r ) cross-covariance Rdf are calculated using the first-order numerical approximation (Zhu and Yeh 2005; Sun et al., 2013):

(2)

or

σ ∇u·n|Γ2 = j

(3)

where u∗ is the electrical potential specified at the boundary Γ1, j denotes the prescribed current density per unit area (A m−2), and n is the unit vector normal to the boundary Γ2 (Yeh et al., 2002; Pidlisecky et al., 2007).

Rdf(r ) = Jd(r ) Rff(r ) , (r ) Rdd = Jd(r ) Rff(r ) Jd(r )T

where

Rff(r ) (N × N)

(6) is the auto-covariance of the conductivity.

Jd(r ) (m × N) is the sensitivity (Jacobian) matrix of potentials with respect to the element-wise parameter, evaluated using the conductivity estimated at the current iteration. An adjoint state variable approach is adopted (Sykes et al. 1985; Li and Yeh 1998; Mao et al. 2013) to compute the sensitivity matrix. For r ≥ 1 the covariance function is updated to represent the conditional (or residual) covariance of the conductivity according to

2.2. SLE for ERT tomography In this study, the inversion of ERT surveys was conducted using a stochastic estimation model – successive linear estimator (SLE), which has been successfully used in hydraulic tomography as well as ERT. For details, refer to Yeh and Liu (2000), Yeh et al. (2002); Liu and Yeh 2004; Zhu and Yeh (2005); Yeh et al. (2008); Xiang et al. (2009), and Zha et al. (2018). Due to spatial variability of the properties of geologic media, SLE first conceptualizes the natural logarithm of electrical conductivity (ln σ (x) ) as a spatial stochastic process (or random field), characterized by a joint probability distribution with given mean, variance, and correlation lengths. This conceptualization is a general description of the spatial varying conductivity of a geologic medium in the geostatistical framework (Yeh et al., 2015a,b) ——a major distinction between SLE and the covariance-based regularized inversion. That is, SLE assumes the conductivity properties at a site in terms of the average value, the deviation from the average, and the average dimensions (length, width, and thickness) of different conductivity clusters (i.e. geologic structures). Likewise, the electrical potential, u (x), is treated as a spatial stochastic process. These random fields are then expressed as ln σ (x) = F + f (x) and u(x) = 〈u〉 + d (x) , where F and 〈u〉 are the mean (i.e., F = 〈ln u (x )〉, the angle bracket denotes expected value); f (x) andd (x) are the perturbations with zero means (i.e., 〈f (x )〉 = 0 and 〈d (x )〉 = 0 ). Afterward, the estimated conductivity vector, given the measured electrical potentials is Fc (subscript c denotes conditional: given the measured potentials). Since the relationship between the conductivity and potential in inverse modeling is nonlinear, the conductivity vector is iteratively determined using the following stochastic linear estimator:

Fc(r + 1) = Fc(r ) + ωTf [D∗ (qk ) − G h (F (r) c , qk)]

(5)

(r ) Rdd

Rff(r + 1) = Rff(r ) − ωTf Rdf(r )

(7)

This update is to reflect the improvement (i.e., reduction of the uncertainty) of the estimate due to the inclusion of the measured potentials. Note that at iteration r = 0, Rff(0) (N × N) is the unconditional covariance of parameters, which could be characterized by an exponential spatial covariance function (or other covariance functions, see Zha et al., 2018):

Rff (P, P′) = ζ f2exp[− (x − x′)2 /λ x2 + (z − z′)2 /λ z2 ]

ζ f2

(8)

In Eq. (8), is the variance of ln σ ; λ x and λ z , are the correlation scales in x and z directions (i.e., the average length and width of geologic layers (Yeh et al., 2015b)), respectively. This covariance implies that the correlation relationship between theσ at the point P (x , z ) and P′ (x ′, z ′) decreases with an increasing separation distance between the two points, normalized by the correlation length. Two convergence criteria are used to terminate the iteration: (1) if a change in the spatial variance of the estimated conductivity field between the current and last iterations is smaller than a specified tolerance, implying that the SLE cannot improve the estimation any further. (2) If a change of simulated potentials between successive iterations is smaller than a given tolerance, indicating that the estimates will not significantly improve the potential field. Once one of the two criteria is met, the estimates are considered to be optimal, and the iteration procedure is terminated. With respect to resolving the ill-posed nature of ERT inversion, SLE takes a different approach to the traditional Occam (e.g., Hu et al., 2013), covariance-based regularized inversion (Hermans et al., 2016; Jordi et al., 2018) and Bayesian inversion of data (Gouveia and Scales, 1997). The Occam approach solves a regularized least-squares problem via a roughening operator (e.g., a discrete second-difference operator) to find the smoothest model. The covariance-based regularized approach yields solutions via using the covariance matrix and optimized damping factor as prior information (e.g., geostatistical regularization operator), while fixing the damping factor to unity would give the same solution as a true geostatistical solution (Tarantola, 2005). The Bayesian approach is philosophically similar to SLE. However, Bayesian inversion generally does not relate a priori model covariance matrix to the spatial statistics of the variability of geophysical properties. Furthermore, SLE is explicitly formulated from a covariance matrix, which represents the spatial characteristics of the variability of conductivity (e.g., mean, variance and correlation scale) of the heterogeneous

(4)

D∗ (qk ) is

a vector consisting of the meawhere r is the iteration index; sured electrical potentials, corresponding to the current transmission at the kth location (qk ); G h (·) denotes Eq. (1), which simulates the potentials induced by the current transmission, using the parameters, Fc(r ) , at iteration r. Note that r = 0, Fc(0) = F (unconditional mean value). The word “unconditional” implies that no measured potential has been used to modify them. Suppose we discretize an estimation domain into N elements. We have the current transmissions at p locations, and each transmission has l potential measurements. As a result, we have a total number of potential measurements m = p × l . With this survey, the dimension of Fc(r ) is an N × 1 vector; D∗ (qk ) andG h (·) are m × 1 vectors. The coefficient matrix, ωT , then has a dimension of N × m, and the superscript T denotes the transpose. This matrix (the weights) assigns the contribution of the difference between the observed and simulated potential at each measurement electrode during each transmission to the previously estimated conductivity value at each element. 525

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

selected reference electrodes (R1, R2, R3, R4, R5, and R6) and the other measuring electrodes. In Case 1, where the profile was surrounded by insulated boundary conditions, one of the reference electrodes was assigned a potential value so that a unique solution can be obtained for the governing equations. In Case 3, where pole-pole configurations were used, the reference electrodes of the transmitters and receivers were assumed to be at the constant potential boundaries. The current transmission locations were the locations of the six pairs of dipole source electrodes in Case 1 and Case 2. Therefore, Case 3 yielded 12 voltage-current datasets, which included a total of 1140 electrical potential measurements. The estimated conditional effective conductivity fields for the three cases were then derived from SLE with the true mean and covariance functions of the reference field as prior information. Estimates from Case 1 (dipole–dipole with insulated boundaries) are displayed in Fig. 1b, while the estimates from Case 2 (dipole–dipole with three prescribed potential boundaries) are illustrated in Fig. 1c. Estimates from Case 3 are shown in Fig. 1d. In comparison with the reference field and the estimated fields from the three cases, we notice that in spite of different configurations, using the constant potential boundaries, the estimated heterogeneous structures as well as the magnitudes of the estimates in Fig. 1c, and 1d are in agreements with the reference field. However, the configurations with fixed electrodes as a reference point in Case 1 yields a low-resolution image (Fig. 1b). Likely, this is attributed to the fact that at insulated boundaries, the potentials are unknown, and more unknowns are to be solved during the inverse modeling. The scatterplots with hexbin of the reference electrical conductivity field versus the estimated for the three cases are shown in Fig. 2a, 2b, and 2c. All the data are scattering around the 45-degree line, indicating that the estimates from SLE are unbiased. Notice that the case 2 and 3 have a higher coefficient of determination (R2): 0.77 and 0.79 respectively. Another statistics norm L2, which is the mean square error of the true field against the estimated field, also indicates that the surveys with the prescribed potential condition have smaller L2 norm (0.173 for dipole–dipole and 0.155 for pole-pole) than that with insulated boundaries (0.322 for case 1). As for norm L1, it has the same trend as norm L2 for the three cases (0.425, 0.312, and 0.293 for case 1, 2, and 3 respectively). Thus, the pole-pole array with three potential boundary sides yielded the best image (i.e., highest resolution). Contour plots of the residual variances (conditional variances), which represent the uncertainty of the estimate at each location, for the three cases are depicted in Fig. 3a, b, and c. We notice that the estimates using survey arrays with four insulated boundaries have relatively higher residual variances, especially around the boundaries, indicating that higher uncertainty of estimated electrical conductivity. These results are consistent with the estimated electrical conductivity fields and their scatterplots (Fig. 2). Additionally, contour plots of the error in the estimate for three cases were shown in Fig. 4a, b, and c. The error was defined by the square of the difference between the estimate and true field. The general behaviors between the residual variance and the error are the same: small differences near observation locations and vice versa. However, since the residual variance is an ensemble statistical concept, it does not necessarily correspond exactly the difference between the estimate and the true.

medium. Using this geostatistical representation of spatial variability, SLE then establishes the cross-correlation between the potential measurements and electrical conductivity everywhere in the domain via the sensitivity based on the governing equation. Afterward, SLE seeks the conditional effective electrical conductivity field (the most likely conductivity values given the measured potentials) at every location as opposed to the conductivity field, which best matches the observed potentials. The likely deviation (uncertainty) of the conditional effective conductivity from the true one is subsequently quantified by the diagonal terms of Rff in Eq. (7) (i.e., residual variance). In addition, the prior model covariance (the residual covariance matrix) is updated at each iteration (Eq. (7)) while it is not considered in the Bayesian approach (see Zha et al., 2018 for a detailed comparison). This update of the residual covariance is designed to reflect the reduced uncertainty of the estimates as more information about the heterogeneity is extracted from the voltage measurements during each iteration. However, this residual variance is an approximate estimate of the uncertainty under the assumptions of the SLE (i.e., multi-Gaussianity). 3. Numerical experiments In cross-hole resistivity tomography surveying, pole-pole and dipole–dipole arrays are widely used for 2D or 3D resistivity imaging, while other types of tri-electrode and four-electrode configurations (e.g., Schlumberger) are similar to a dipole–dipole array (Bing and Greenhalgh, 2000). As for the dipole–dipole configuration, the potential difference was measured between the two potential electrodes (measurement pairs) along with two opposite-polarity current electrodes (source pairs) at boreholes. However, the pole-pole array has two remote electrodes (i.e., one of the current electrodes and one of the potential electrodes). They are placed far away from the other electrodes regarded as zero-voltage reference points (Shima, 1992; Bing and Greenhalgh, 2000). Since the complete pole-pole data sets have maximum information and other configuration data sets can be derived from them, the pole-pole array theoretically should produce the best estimate of the study region under a low noise level condition (Bing and Greenhalgh, 2000). Prior to the sandbox experiments (unknown noise), numerical experiments (noise free) were conducted to investigate the effects of different ERT surveys (i.e., dipole–dipole or pole-pole configurations) under different boundary conditions of the laboratory sandbox experiments. For this purpose, the numerical experiments consider a two-dimensional, vertical profile (180 by 90 cm). This profile is discretized into 36 vertical and 36 horizontal elements (totally 1296 elements) with a dimension of 5 cm (length) ×2.5 cm (width). The electrical conductivity value for each element was created using a stochastic random field generator with a geometric mean of 0.01 S m−1 (-4.95 in log scale), a variance of 0.0001 S2 m−2 (0.7 in log scale), and an exponential correlation structure with a horizontal correlation scale of 90 cm and a vertical correlation scale of 5 cm. This heterogeneous conductivity field (Fig. 1a) is the reference field, in which six vertical lines of electrodes for the ERT survey were installed at a horizontal interval of 30 cm. Along each line, sixteen electrodes were placed at 5 cm interval (the solid circles in Fig. 1a). Three cases were examined in this numerical experiment. Case 1 represents the situation where dipole–dipole configures were used with insulated boundaries on the four sides of the profile; Case 2 is the same as Case 1 but the top, right-hand, and left-hand sides of the profile were replaced with a constant electrical potential boundary. Pole-Pole configurations were used in Case 3 in which the top, right-hand, and lefthand sides of the profile were assigned a constant electrical potential boundary. In Cases 1 and 2, where the dipole–dipole array was used, six different current transmitter pairs (A1-B1, A2-B2, A3-B3, A4-B4, A5-B5, and A6-B6 in Fig. 1a) were employed to yield six voltage-current datasets. Each set consists of 93 voltage difference measurements between the

4. Sandbox experiments 4.1. Setup To illustrate the effectiveness and merits of our new ERT design in a laboratory sandbox, we constructed a heterogeneous soil profile in a 180 cm long, 90 cm high, and 10 cm thick sandbox with 19 layers of sands of different grain sizes (Fig. 5a and b). Based on the numerical experiments, we installed 96 copper electrodes at 5 cm vertical spacing and 30 cm horizontal spacing on one side of the sandbox (black circles 526

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

(a) Reference Model

(c) Case 2: Dipole-dipole with prescribed potential boundaries

(S/m)

(S/m)

0.075

B4

B2 B5

Y (cm)

60 40

B6

R4 B1

B3

A3

R2 R5

R6 R3 A 6 A4

A1 A5

A2

20 0 0

0.065 0.055 0.045

R1

0.075

80

0.035

100

150

0.055 0.045

40

0.035 0.025

0.025

20

0.015

50

0.065

60

Y (cm)

80

0.005

0 0

0.015

50

X (cm)

100

0.005

150

X (cm)

(b) Case 1: Dipole-dipole with no-flux boundaries

(d) Case 3: Pole-pole with prescribed potential boundaries

(S/m)

(S/m)

0.075

80

0.065 0.055 0.045

40

0.035

0.065 0.055

60

Y (cm)

60

Y (cm)

0.075

80

0.045

40

0.035 0.025

0.025

20 0 0

0.015

50

100

150

0.005

X (cm)

20 0 0

0.015

50

100

150

0.005

X (cm)

Fig. 1. (a) The reference conductivity field. Black solid rectangles indicate the locations of six dipole sources which marked by A1-B1, A2-B2, A3-B3, A4-B4, A5-B5, and A6-B6. The red ellipses (R1, R2, R3, R4, R5, and R6) are the locations of associated reference electrode at each transmission. White circles show the distribution of 96 electrodes. Estimated conditional effective conductivity for cases: (b) dipole-dipole array with insulated boundary conditions at four sides; (c) dipole-dipole array with a prescribed potential boundary condition at left, right and top three sides; (d) pole-pole array with a prescribed potential boundary condition at left, right and top three sides.

in Fig. 5c) after the sandbox was packed with sands. The grain size, saturated hydraulic conductivity, and porosity for five materials were listed in Table 1. Since thickness-averaged electrical potentials were our primary interest, the length of the electrode was made long enough to penetrate 90% of the thickness of the sandbox, avoiding influences of the water on the other side of the sandbox wall. In addition, to decrease the contact resistance between the electrodes and the sands, we enlarged the electrode contact surface area as much as possible. Hence, copper electrodes were designed as 9 cm in length with 0.18 cm in diameter. The front part (approximately 8 cm) of the electrode was in direct contact with the sands, while the rest (about 1 cm) was wrapped with non-conductive PVC heat-shrink tube and fixed on the glass to avoid possible impacts of water between the sand and this side of the wall. As illustrated in Fig. 5a and c, both sides and top of the sandbox are compartments for water reservoirs (transparent boundaries in Fig. 5a or areas bounded by the solid black lines and dotted lines in Fig. 5c). Before the experiments, the sandbox was repeatedly filled with tap water by increasing the water level in the reservoirs on the two sides to completely fill the upper reservoir and then drained to ensure that the sand is fully compacted and settled. The water was drained from the lower left and lower right outlets of the two side reservoirs.

sandbox experiments. To create a constant potential boundary along the three sides of the sandbox, we inserted a continuous copper wire into the sands around the upper, left, and right sides of the sandbox (red line in Fig. 5c). This wire is connected to the negative pole of the current source as well as the negative pole of the receiver, as such this wire serves as a constant zero electrical potential boundaries. The positive pole of the source then introduced the direct current to one of the 96 electrodes, and the potentials at the other 95 electrodes are monitored (i.e., one snapshot). A given current was transmitted at selected ten different locations (the black solid rectangles and dashed ellipses in Fig. 5c, A1, A2, A3, A4, A5, A6, A7, and A8 for calibration; A9, and A10 for validation in Section 4.6) individually to scan the sandbox under drainage and fully saturated conditions. All the datasets acquired for above the configurations were taken by AGI SuperSting 8-channel resistivity meter. The measurement time was about 0.4 sec, which including current injection period and data recording time, for each individual. Following each measurement, the current source was shut down for 0.4 sec to reduce induced polarization effect. Generally, 10 snapshots for pole-pole array took less than 20 min.

4.2. ERT data acquisition

The aforementioned data acquisition method was applied to the sandbox under a fully saturated condition (T = 0), and unsaturated conditions with different drainage times (T = 0.5, 1.0, 1.5, 2.0 and 2.5 h). This drainage experiment was conducted based on the fact that

4.3. Electrical potential fields at different drainage times

Based on the results of our numerical experiments, we employed a pole-pole approach with constant potential boundaries during the 527

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

Fig. 2. The hexbin scatterplots of the reference field versus the estimated field for the three cases are shown as (a), (b), and (c), respectively.

after drainage, fine-textured sands will retain more water, and the coarse ones will retain less due to the difference in capillary forces as documented in the results of the field study by Ye et al. (2005). In effect, stochastic analysis of unsaturated flow by Yeh et al. (1985b) also suggested that variances of local-scale pressure head, moisture content, and unsaturated hydraulic conductivity increase as the average moisture content of a medium decreases. We thus expect that after drainage, the electric conductivity of the fine-textured sands should be higher than those of the coarse-textured sands, the contrasts in electric conductivity of different layers should increase, and the general shape of different layers could be detected. In Fig. 6, the contour maps show the electric potential fields after current injection at a location (0.75, 0.7) for the fully saturated condition (T = 0), and the unsaturated condition at T = 0.5, 1.0, 1.5, 2.0 and 2.5 h. All the potentials [in mV] were scaled to 1 mA current injection. For the saturated condition, the contour map shows concentric circles around the injection point due to the nearly homogeneous resistivity distributions. After 0.5 h of drainage, the potential distributions became elliptic. They stretched in horizontal directions at later times, suggesting anisotropic behaviors reflecting layering effects. Similar behaviors were observed when the sources were at other locations. These observations seem to corroborate the moisture-dependent anisotropy theory for unsaturated hydraulic conductivity (Yeh et al., 1985b). The theory stated that water in stratified soils tend to move horizontally as the soils become less saturated. That is, the electrical potential field behaves in a manner similar to the pressure head or moisture content field in flow through porous media, even though the underlining physics is different. Nevertheless, these potential maps did not clearly reveal the layer structure of the sands in the sandbox as some might expect. On the other hand, the smoothness of the map is indicative of the diffusive nature of the electrical potential: the change in potentials due to variation in resistivity is small. The smoothness of the potential field is similar to the drawdown distribution induced by a

pumping test in an aquifer: aquifer heterogeneity does not affect drawdown distribution significantly. The results are consistent with the results of the stochastic analysis of the effects of heterogeneity on the drawdown of groundwater flow (Gelhar, 1993). 4.4. SLE inversion Results. The raw data from the ERT surveys were then inputted to SLE algorithm to estimate the electric conductivity fields of the sandbox at the fully saturated condition, 0.5, 1.0, 1.5, 2.0, and 2.5 h after drainage. For these inverse modeling, the prior mean conductivity was 0.05 S m−1 (−3.2 in log scale), and its variance was 0.0018 S2 m−2 (0.5 in log scale); they were estimates based on apparent resistivity ρa [in ohm m] calculated by differential potentials (Δv ) and injected current (I ):

ρa = k

Δv I

(9)

where k is the geometric factor (see Loke et al., 2013) depending on the configurations. A value of 5 cm was assigned for the vertical correlation scale, and 90 cm for the horizontal correlation scale. These correlation scales were estimated based on the visual examination of the average thickness and length of the layers in the sandbox. The estimated conductivity field under a fully saturated condition is illustrated in Fig. 7a. It does not reveal the layer structure, likely due to the dominance of the electrical conductivity of the tap water and small contrast in the porosity of the different sands (around 0.37–0.39, see table 1) according to Archie’s law. Stochastic theories of groundwater flow (Gelhar, 1993) and unsaturated flow (Yeh et al., 1985a,b,c) may suggest a solution to this issue. These theories stated the head variance will increase in a heterogeneous media if the hydraulic gradient increases even if the heterogeneity remains the same. Thus, a large injected current should increase the potential variation reflecting heterogeneity such that it can be detected clearly by the sensor. Under this 528

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

(a) Case1: Dipole-dipole with no-flux boundaries

Var ln 80

Y (cm)

60

40

20

0

0

20

40

60

80

100

X (cm)

120

140

160

180

0.40 0.36 0.32 0.28 0.24 0.20 0.16 0.12 0.08 0.04 0.00

(b) Case2: Dipole-dipole with prescribed potential boundaries

Var ln 80

Y (cm)

60

40

20

0

0

20

40

60

80

100

X (cm)

120

140

160

180

0.40 0.36 0.32 0.28 0.24 0.20 0.16 0.12 0.08 0.04 0.00

(c) Case3: Pole-pole with prescribed potential boundaries

Var ln 80

Y (cm)

60

40

20

0

0

20

40

60

80

100

X (cm)

120

140

160

180

0.40 0.36 0.32 0.28 0.24 0.20 0.16 0.12 0.08 0.04 0.00

Fig. 3. Residual variances (conditional variances) of estimated lnσ for the three cases. White circles show the distribution of 96 electrodes.

conductive the medium is. The variance of the lnσ field is 0.3 at the saturation, rises rapidly and stabilizes at 2.0 with increasing of the drainage time, consistent with our earlier speculation. To further substantiate our speculation that increasing the resistivity contrast in the sandbox, ERT will be able to delineate the layering structure. The estimated conductivity distribution at 2.5 h after drainage and its residual variance map are illustrated in Fig. 8a and b, respectively. The black lines in Fig. 8a and b are the outlines of the 19 layers based on our visualization. Examining Fig. 8a, we notice that the layer structures from ERT estimate generally agree with the layer structure visually observed from the sandbox. The coarse sand layers have low electrical conductivity values and the fine sands have high

large input voltage, ERT may be able to detect the layer structure in the sandbox. Nevertheless, verification of this hypothesis is out of the scope of this study. As discussed earlier, the drainage of the sandbox could increase the resistivity contrasts in the sandbox. To support this hypothesis, the estimated variances and the spatial mean of the estimated conductivity fields at the fully saturated condition, 0.5, 1.0, 1.5, 2.0, and 2.5 h after drainage are shown in Fig. 7b. The spatial mean (average over the entire sandbox) of the estimated conductivity field changed from 0.1 S/ m at the saturation to around 0.05 S/m at the end of the drainage period (see the right-hand side axis). This decreasing trend in spatial mean agrees with the fact that less water in the sandbox, the less

529

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

(a) Case1: Dipole-dipole with no-flux boundaries Err. [ - ] 80

Y (cm)

60

40

20

0

0

20

40

60

80

100

X (cm)

120

140

160

180

1.00 0.90 0.80 0.70 0.60 0.50 0.40 0.30 0.20 0.10 0.00

(b) Case2: Dipole-dipole with prescribed potential boundaries Err. [ - ] 80

Y (cm)

60

40

20

0

0

20

40

60

80

100

X (cm)

120

140

160

180

1.00 0.90 0.80 0.70 0.60 0.50 0.40 0.30 0.20 0.10 0.00

(c) Case3: Pole-pole with prescribed potential boundaries Err. [ - ] 80

Y (cm)

60

40

20

0

0

20

40

60

80

100

X (cm)

120

140

160

180

1.00 0.90 0.80 0.70 0.60 0.50 0.40 0.30 0.20 0.10 0.00

Fig. 4. The error, which is defined by the square of the difference between the estimate and true field, in the parameter estimation for the three cases. White circles show the distribution of 96 electrodes.

(Sun et al., 2013). Such influences are controlled by the horizontal and vertical correlation scales as well as the sensitivity of potential to change in electric conductivity. A plot of the estimated conductivity field, using a vertical correlation scale of 5 cm and a horizontal correlation scale of 10 cm (smaller than the average length of the physical layers) as prior information, is displayed in Fig. 8c. The estimated field exhibits several localized anomalously high conductivity blocks in comparison with the estimates in Fig. 8a, which was obtained by using a horizontal correlation scale of 90 cm. Besides, the layer structures near the boundaries of the sandbox are not visible. These high conductivity zones and blurred boundaries suggest that the covariance function and the sensitivity do not capture the true spatial relationship between the measured potential and adjacent conductivity anomalies. On the other hand, based on the

values. The conductivities of the regions at the bottom of the sandbox (i.e., blue layer 1 and orange layer 2) were high due to the accumulation of undrained water at the bottom. 4.5. Importance of prior information for SLE SLE was used to interpret the ERT survey data. This inverse (or estimation) problem is underdetermined or ill-defined (more parameters than potential measurements) and cannot be solved exactly. In order to overcome this difficulty, SLE uses the spatial covariance function, Rff (spatial correlation between parameters) as prior information to regularize the problem. The covariance function is an exponential function, which extends the influences of the sensitivity of potential measurements to conductivity at different parts of the profile 530

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

Fig. 5. (a) The photograph of the sandbox, (b) a sketch of different grain size distribution, (c) The sketch of setup for pole-pole ERT survey with the continuous copper wire (red lines). Eight black rectangles are the locations of current injection for the calibration and two dashed ellipses for the validation.

i.e., root mean square error) between the potential predictions and observations (Linde et al., 2015). As expected, the RMS for the estimates in Fig. 8c has much larger value (about 23.07 mV) than that of the estimates in Fig. 8a (around 17.65 mV).

experiments not shown here, we notice that using values greater than 90 cm as the horizontal correlation scale did not further improve the estimates in comparison with those in Fig. 8a. Another way to assess the different estimated models is by comparing the data misfits (or RMS, 531

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

anomalies of our interests (Yeh et al., 2015b), and this was verified in Yeh and Liu (2000) for hydraulic tomography. Specifically, the optimal distance of the monitoring network at least should be smaller than 0.5 times the correlation scale of heterogeneity of our interest. For instance, the estimated conductivity field using data collected from the vertical electrode spacing of 10 cm is presented in Fig. 8d. The estimates in Fig. 8d are not as detailed as those in Fig. 8a. This finding, therefore, could be of practical importance to the design of the ERT layout.

Table 1 The grain size, saturated hydraulic conductivity and porosity for five materials. Grain size (mm)

0.1–0.25

0.25–0.4

0.3–0.6

0.6–1

1–4

K (cm/s) Porosity

0.0180 0.3892

0.0788 0.377

0.1395 0.3712

0.3352 0.3737

0.8527 0.3788

Although the covariance function is analogous to the regularization term in other inverse methods such as Occam’s inversion or Bayesian approach, the correlation scales in the covariance function of SLE or other geostatistics-based algorithms (e.g., Kitanidis, 1995) were derived from visual examinations of the sandbox (Fig. 5a or b). Under field situations, examinations of outcrops, geologic map, or borehole logs could suggest some approximate correlation scales. Yeh et al. (2015b) pointed out that the correlation scales do not have to be exact but they have to be reasonably close to the average dimensions of the heterogeneity. They are soft constraints to the inverse modeling efforts; the effects will be overridden once measurements become abundant (i.e., many samples within the correlation scale). Note that regularization approaches using covariance-based solution likely reach the similar conclusion about the effect of the choice of prior covariance model (e.g., Hansen et al., 2006). The prior information determines the unconstrained part of model with little measurements and is honored only when it is consistent with the observed data. The correlation scales also play an important role in the design of electrode spacing during ERT surveys. As a rule of thumb, the receiver electrode spacing should be smaller than the dimension of the dominant

4.6. Validation of the estimated parameters The true electrical conductivity distributions of a field site, as well as the sandbox in this study, are unknown. One possible approach to quantitatively assess the accuracy of the estimates from ERT is to validate the estimates using independent surveys. Independent surveys should be conducted at new transmitter and receiver locations, which were never used in the previous inverse modeling effects. That is, had the conductivity field been completely characterized by the ERT survey, the estimated field would have predicted exact electric field induced by the independent surveys. Otherwise, the predicted electric field will deviate from the observed as reported in hydraulic tomography by Huang et al. (2011). In order to carry out this validation exercise, we chose to use SLE to estimate a new conductivity field of the sandbox at 2.5 h after drainage, using potential measurements at 18 electrodes (white circles in Fig. 9a). This is because the aforementioned ERT inversion, all 96 electrodes were used as either transmitters or receivers. To create the new estimated field, the current transmission locations

Fig. 6. Contour maps of potential fields after current injection at point A5 (0.75, 0.7) for (a) the fully saturated condition, and (b)-(f) the unsaturated condition with different drainage time. All the potentials [in mV] were scaled to 1 mA current injection. White lines show the six vertical lines of electrodes. 532

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

Fig. 7. (a) Estimated conductivity field under a fully saturated condition. Red rectangles indicate the position of transmissions and white circles are the receivers. (b) Spatial variances and spatial mean of the estimated conductivity fields at different times after drainage of the sandbox. The shedded bands indicate the range of measuring time for each snapshot.

electrodes were not used as either transmission or receiver electrodes in the creation of the new estimated field. The validation results – scatter plots of the predicted vs. the measured potentials – are displayed in Fig. 9c and d, respectively, for A9 and A10 transmitters. As expected, the predicted potentials versus the observed are generally scattering along the 45-degree line, indicating that SLE yields unbiased conditional effective conductivity field but at a lower resolution (i.e. larger scatterings in the scatter plots in Fig. 9c and d than those in Fig. 9b. These large scatterings confirm that these two validation tests are independent of the tests used in the model calibration since they carry non-redundant information (Wen et al., 2019). More importantly, the new estimated electrical conductivity field can predict the potential

and the prior information of SLE remain the same as those in Fig. 8a (eight red rectangles in Fig. 9a). This arrangement left us a total of 72 electrodes for validation. The new estimated field is illustrated in Fig. 9a, and the scatter plot of the calibrated vs. measured potentials is shown in Fig. 9b. As shown in the scatter plot, the calibration is as good as those based on dense measurement electrodes. Nonetheless, the newly estimated field (Fig. 9a) has a lower resolution than those in Fig. 8a and d, because the potential measurement interval in the vertical direction is 25 cm, much larger than 5 cm and 10 cm in Fig. 8a. This newly estimated field was then used to predict the two potential fields, measured at those discarded electrode rows, induced by the current transmissions at locations A9 and A10. These transmission

Fig. 8. (a) Estimated electrical conductivity field. (b) The corresponding residual variance. The black lines overlaid on the maps (a) and (b) indicate the layer structures observed from the sandbox. (c) Estimated electrical conductivity field obtained using a vertical correlation scale of 5 cm and a horizontal correlation scale of 10 cm. (d) Estimated conductivity fields using data collected from the vertical electrode spacing of 10 cm (while correlation scale x = 90 cm, y = 5 cm). White circles show the distribution of 96 electrodes and red rectangles for current source locations. 533

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

Fig. 9. (a) Estimated conductivity fields using data collected from the vertical electrode spacing of 25 cm. (b) Scatterplots of the calibrated electric potentials versus measured. (c) and (d) validation tests using (c) A9 and (c) A10, respectively.

fields due to stresses at different locations unbiasedly with some discrepancies. These discrepancies are the effects of unresolved conductivity heterogeneity resulting from the low density of the survey. On the other hand, the results also imply that these two independent tests could increase the resolution of the estimates if they were included in the inverse modeling—they carry new information about the heterogeneity. While a similar concept has been used in hydraulic tomography [e.g., Liu et al. (2007), Xiang et al. (2009) and others], this independent test approach has not been reported in the literature of ERT surveys. This is particularly of importance in hydrogeophysics since ERT results have been generally compared with the estimates from borehole resistivity surveys. Since ERT and borehole geophysics are at different scales, the comparison may be limited to large-scale features only. For this reason, we believe that this validation approach may be a new way to substantiate the ERT survey estimates in the field of hydrogeophysics.

the current source) as long electrodes near the target to electrically monitor the waste plumes in the vadose zone at the T tank farm at the Hanford nuclear site, which overcomes the smearing issue by the upper conductive layer and obtained better deeper subsurface information. Yang et al. (2017) present a 3D ERT modeling with long electrode sources. Their results showed that the arrays with long electrode sources are more sensitive to the subsurface anomalies than those with point electrode sources. Based on these studies, it is rational to speculate that steel wells (i.e., long electrodes) or steel sheets can be used as constant potential boundary conditions by connecting to the source of the ERT survey. These studies and our results may lead to new research on field applications of this approach, although it might be a limitation to extend the method to the field at this moment. Unlike most ERT inversion methods, which are based on Tikhonov regularization, the SLE approach inherently uses physically meaningful correlation scales as prior information. These correlation scales do not have to be exact, they can be guessed based on borehole logs. For example, the vertical correlation scale could be easily obtained from the average thickness of layers in boreholes while the horizontal scales can be estimated from the continuity of stratifications between boreholes. With this prior information, SLE yields the unbiased conditional effective electrical conductivity field and quantifies the uncertainties of the estimated field by the residual (conditional) variance. The locations of large differences between the true and estimated electrical conductivity fields, on the average, coincide the locations where the SLE predicts high conditional variance. Also, the correlation scale concept further suggests the electrode spacing be smaller or equal to the correlation scales of dominant geologic heterogeneity of our interest. This study also demonstrates that independent tests could be a viable means to validate the estimates of ERT surveys. This validation concept involves the issue of redundant, moderately, and highly nonredundant datasets, which are related to the change of current flow fields (Wen et al., 2019), not just a jackknifing resampling approach.

5. Discussions and conclusions Our study suggests that the ERT survey using a pole-pole array with constant potential boundaries is more effective than using the dipole–dipole array. To implement constant potential boundaries in our sandbox experiments, we embedded a continuous copper wire in the sands on three sides of the sandbox and connected it to the negative pole of the current source as well as the negative pole of the receiver. Such a pole-pole setup facilitated a fast measurement and better conductivity estimation and lowered its uncertainties. Further, it expedited convergence of the estimation because of the precise boundary conditions (i.e., zero electrical potential). Albeit, the above approach appears to be limited to the sandbox experiment, it has the potential for field applications. As a matter of the fact, Rucker et al. (2010, 2011) used steel-casing wells (connected to 534

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

Validation means verifying if the estimated parameters from an inverse model can yield unbiased predictions of processes under different stresses. While the concept has been used in hydraulic tomography, this independent test approach has not been reported in the literature of ERT surveys. We further show that ERT with appropriate measurement setups and SLE can yield realistic images of layers of different sands in the sandbox under unsaturated conditions. This is attributed to that unsaturated conditions increase the resistivity contrast (due to the difference in moisture content) among different sands such that the layer structures are visible to ERT survey. The contrast in resistivity was due to the differences in water retained in each layer after drainage, which is related to the pore-size distribution of the layer. Again, such a condition is not limited to our sandbox experiment. For example, Ye et al., (2005) used neutron probes to detect 3D moisture distribution at Hanford site and they found that the field site was under steady-state, unit gradient, drainage conditions, fine-textured layers remain wet while coarse-textured ones remain dry. Moisture distributions at the site were consistent with the geologic structure. Therefore, it is likely that ERT surveys could be used to map the spatial continuity of layers at the site.

2016.02.008. Cockett, R., Kang, S., Heagy, L.J., Pidlisecky, A., Oldenburg, D.W., 2015. SimPEG: An open source framework for simulation and gradient based parameter estimation in geophysical applications. Comput. Geosci. 85, 142–154. https://doi.org/10.1016/j. cageo.2015.09.015. Constable, S.C., Parker, R.L., Constable, C.G., 1987. Occam’s inversion: a practical algorithm for generating smooth models from electromagnetic sounding data. Geophysics 52 (3), 289–300. https://doi.org/10.1190/1.1442303. Deng, Y.P., Shi, X.Q., Xu, H.X., Sun, Y.Y., Wu, J.C., Revil, A., 2017. Quantitative assessment of electrical resistivity tomography for monitoring DNAPLs migration - comparison with high-resolution light transmission visualization in laboratory sandbox. J. Hydrol. 544, 254–266. https://doi.org/10.1016/j.jhydrol.2016.11.036. Gelhar, L.W., 1993. Stochastic Subsurface Hydrology. Prentice-Hall, Upper Saddle River, N. J. Gouveia, W.P., Scales, J.A., 1997. Resolution of seismic waveform inversion: Bayes versus Occam. Inverse Prob. 13 (2), 323. https://doi.org/10.1088/0266-5611/13/2/009. Gunther, T., Rucker, C., Spitzer, K., 2006. Three-dimensional modelling and inversion of dc resistivity data incorporating topography - II. Inversion. Geophys. J. Int. 166 (2), 506–517. https://doi.org/10.1111/j.1365-246X.2006.03011.x. Hansen, T.M., Journel, A.G., Tarantola, A., Mosegaard, K., 2006. Linear inverse Gaussian theory and geostatistics. Geophysics 71 (6), R101–R111. https://doi.org/10.1190/1. 2345195. Hayley, K., Bentley, L.R., Gharibi, M., 2009. Time-lapse electrical resistivity monitoring of salt-affected soil and groundwater. Water Resour. Res. 45 (7), W07425. https://doi. org/10.1029/2008WR007616. Hermans, T., Irving, J., 2017. Facies discrimination with electrical resistivity tomography using a probabilistic methodology: effect of sensitivity and regularisation. Near Surf. Geophys. 15, 13–25. https://doi.org/10.3997/1873-0604.2016047. Hermans, T., Kemna, A., Nguyen, F., 2016. Covariance-constrained difference inversion of time-lapse electrical resistivity tomography data. https://doi.org/10.1190/ GEO2015-0491.1. Hermans, T., Nguyen, F., Caers, J., 2015. Uncertainty in training image-based inversion of hydraulic head data constrained to ERT data: Workflow and case study. Water Resour. Res. 51, 5332–5352. https://doi.org/10.1002/2014WR016460. Hermans, T., Vandenbohede, A., Lebbe, L., Martin, R., Kemna, A., Beaujean, J., Nguyen, F., 2012. Imaging artificial salt water infiltration using electrical resistivity tomography constrained by geostatistical data. J. Hydrol. 438–439, 168–180. https://doi. org/10.1016/J.JHYDROL.2012.03.021. Hu, X., Peng, R., Wu, G., Wang, W., Huo, G., Han, B., 2013. Mineral exploration using CSAMT data: application to Longmen region metallogenic belt, Guangdong Province, China. Geophysics 78 (3), B111–B119. https://doi.org/10.1190/geo2012-0115.1. Huang, S.-Y., Wen, J.-C., Yeh, T.-C.J., Lu, W., Juan, H.-L., Tseng, C.-M., Lee, J.-H., Chang, K.-C., 2011. Robustness of joint interpretation of sequential pumping tests: numerical and field experiments. Water Resour. Res. 47, W10530. https://doi.org/10.1029/ 2011WR010698. Illman, W.A., Liu, X., Craig, A., 2007. Steady-state hydraulic tomography in a laboratory aquifer with deterministic heterogeneity: Multi-method and multiscale validation of hydraulic conductivity tomograms. J. Hydrol. 341 (3), 222–234. https://doi.org/10. 1016/j.jhydrol.2007.05.011. Illman, W.A., Liu, X., Takeuchi, S., Yeh, T.-C.J., Ando, K., Saegusa, H., 2009. Hydraulic tomography in fractured granite: Mizunami Underground Research site, Japan. Water Resour. Res. 45, W01406. https://doi.org/10.1029/2007WR006715. Johnson, T.C., Routh, P.S., Clemo, T., Barrash, W., Clement, W.P., 2007. Incorporating geostatistical constraints in nonlinear inversion problems. Water Resour. Res 43, 10422. https://doi.org/10.1029/2006WR005185. Jordi, C., Doetsch, J., Günther, T., Schmelzbach, C., Robertsson, J.O.A., 2018. Geostatistical regularization operators for geophysical inverse problems on irregular meshes. Geophys. J. Int. Geophys. J. Int 213, 1374–1386. https://doi.org/10.1093/ gji/ggy055. Kitanidis, P.K., 1995. Quasi-Linear Geostatistical Theory for Inversing. Water Resour. Res. 31, 2411–2419. https://doi.org/10.1029/95WR01945. Lapenna, V., Lorenzo, P., Perrone, A., Piscitelli, S., Rizzo, E., Sdao, F., 2005. 2D electrical resistivity imaging of some complex landslides in the Lucanian Apennine chain, southern Italy. Geophysics 70 (3), B11–B18. https://doi.org/10.1190/1.1926571. Lelievre, P., Oldenburg, D., Williams, N., 2009. Integrating geological and geophysical data through advanced constrained inversions, ASEG Extended Abstracts 1, 1–6. https://doi.org/10.1071/ASEG2009ab044. Li, B., Yeh, T.-C., 1998. Sensitivity and moment analyses of the head in variably saturated regimes. Adv. Water Resour. 21 (6), 477–485. https://doi.org/10.1016/S03091708(97)00011-0. Linde, N., Lochbühler, T., Dogan, M., Van Dam, R.L., 2015. Tomogram-based comparison of geostatistical models: Application to the Macrodispersion Experiment (MADE) site. J. Hydrol. 531, 543–556. https://doi.org/10.1016/j.jhydrol.2015.10.073. Liu, S., Yeh, T.C.J., 2004. An integrative approach for monitoring water movement in the vadose zone. Vadose Zone J. 3 (2), 681–692. https://doi.org/10.2136/vzj2004.0681. Liu, X., Illman, W., Craig, A., Zhu, J., Yeh, T.C., 2007. Laboratory sandbox validation of transient hydraulic tomography. Water Resour. Res. 43 (5), W05404. https://doi. org/10.1029/2006WR005144. Loke, M.H., Chambers, J.E., Rucker, D.F., Kuras, O., Wilkinson, P.B., 2013. Recent developments in the direct-current geoelectrical imaging method. J. Appl. Geophys. 95, 135–156. https://doi.org/10.1016/j.jappgeo.2013.02.017. Loke, M.H., Dahlin, T., Rucker, D.F., 2014. Smoothness-constrained time-lapse inversion of data from 3D resistivity surveys. Near Surface Geophys. 12 (1), 5–24. https://doi. org/10.3997/1873-0604.2013025. Mao, D.Q., Lu, L., Revil, A., Zuo, Y., Hinton, J., Ren, Z.J., 2016. Geophysical Monitoring of Hydrocarbon-Contaminated Soils Remediated with a Bioelectrochemical System.

Declaration of Competing Interest None. Acknowledgments This research was supported by the China Scholarship Council (Grant 201706410014), the National Natural Science Foundation of China (Grant 41630317, 41474055); the Fundamental Research Funds for the for the Central Universities, China University of Geosciences, China (Grant CUGCJ1707) and the National Key R&D Program of China (Grant 2016YFC0600201-7). References Archie, G.E., 1942. The electrical resistivity log as an aid in determining some reservoir characteristics. Trans AIME 146 (01), 54–62. https://doi.org/10.2118/942054-G. Binley, A., Cassiani, G., Deiana, R., 2010. Hydrogeophysics: opportunities and challenges. Bollettino Di Geofisica Teorica Ed Applicata 51 (4), 267–284. Binley, A., Henry-Poulter, S., Shaw, B., 1996. Examination of Solute Transport in an Undisturbed Soil Column Using Electrical Resistance Tomography. Water Resour. Res. 32 (4), 763–769. https://doi.org/10.1029/95WR02995. Binley, A., Hubbard, S., Huisman, J., Revil, A., Robinson, D., Singha, K., Slater, L.D., 2015. The emergence of hydrogeophysics for improved understanding of subsurface processes over multiple scales. Water Resour. Res. 51, 3837–3866. https://doi.org/ 10.1002/2015WR017016. Bing, Z., Greenhalgh, S.A., 2000. Cross-hole resistivity tomography using different electrode configurations. Geophys. Prospect. 48 (5), 887–912. https://doi.org/10.1046/j. 1365-2478.2000.00220.x. Blaschek, R., Hördt, A., Kemna, A., 2008. A new sensitivity-controlled focusing regularization scheme for the inversion of induced polarization data based on the minimum gradient support. https://doi.org/10.1190/1.2824820. Bohling, G.C., Zhan, X.Y., Butler, J.J., Zheng, L., 2002. Steady shape analysis of tomographic pumping tests for characterization of aquifer heterogeneities. Water Resour. Res. 38 (12), 17. https://doi.org/10.1029/2001wr001176. Bouchedda, A., Bernard, G., Gloaguen, E., 2017. Constrained electrical resistivity tomography Bayesian inversion using inverse Matérn covariance matrix. Geophysics 82, E129–E141. https://doi.org/10.1190/geo2015-0673.1. Cardiff, M., Barrash, W., Kitanidis, P., Malama, B., Revil, A., Straface, S., Rizzo, E., 2009. A Potential-Based Inversion of Unconfined Steady-State Hydraulic Tomography. Groundwater 47, 259–270. https://doi.org/10.1111/j.1745-6584.2008.00541.x. Cardiff, M., Barrash, W., Kitanidis, P.K., 2012. A field proof-of-concept of aquifer imaging using 3-D transient hydraulic tomography with modular, temporarily-emplaced equipment. Water Resour. Res. 48 (5), 213–223. https://doi.org/10.1029/ 2011WR011704. Chen, J.L., Chen, C.H., Kuo, C.L., Fen, C.S., Wu, C.C., 2016. Estimating groundwater velocity using apparent resistivity tomography: a sandbox experiment. IOP Conference Series. Earth Environ. Sci. 39, 012056. https://doi.org/10.1088/17551315/39/1/012056. Clement, R., Moreau, S., 2016. How should an electrical resistivity tomography laboratory test cell be designed? Numerical investigation of error on electrical resistivity measurement. J. Appl. Geophys. 127, 45–55. https://doi.org/10.1016/j.jappgeo.

535

Journal of Hydrology 575 (2019) 523–536

D. Xu, et al.

211–229. https://doi.org/10.1016/S0926-9851(02)00124-6. Slater, L., Binley, A.M., Daily, W., Johnson, R., 2000. Cross-hole electrical imaging of a controlled saline tracer injection. J. Appl. Geophys. 44 (2–3), 85–102. https://doi. org/10.1016/S0926-9851(00)00002-1. Sun, R., Yeh, T.C.J., Mao, D., Jin, M., Lu, W., Hao, Y., 2013. A temporal sampling strategy for hydraulic tomography analysis. Water Resour. Res. 49 (7), 3881–3896. https:// doi.org/10.1002/wrcr.20337. Sykes, J., Wilson, J., Andrews, R., 1985. Sensitivity analysis for steady-state groundwater flow using adjoint operators. Water Resour. Res. 21 (3), 359–371. https://doi.org/10. 1029/WR021i003p00359. Tarantola, A., 2005. Inverse problem theory and methods for model parameter estimation. Society for Industrial and Applied Mathematics, Philadelphia, PA. Wen, J.C., Chen, J.-L., Yeh, T.-C.J., Wang, Y.-L., Huang, S.-Y., Tian, Z., Yu, C.-Y., 2019. Redundant and non-redundant information for Model Calibration or Hydraulic Tomography, in press. Ground Water. https://doi.org/10.1111/gwat.12879. Xiang, J.W., Yeh, T.C.J., Lee, C.H., Hsu, K.C., Wen, J.C., 2009. A simultaneous successive linear estimator and a guide for hydraulic tomography analysis. Water Resour. Res. 45 (2), W02432. https://doi.org/10.1029/2008WR007180. Xu, D., Hu, X.-Y., Shan, C.-L., Li, R.-H., 2016. Landslide monitoring in southwestern China via time-lapse electrical resistivity tomography. Appl. Geophys. 13 (1), 1–12. https:// doi.org/10.1007/s11770-016-0543-3. Yang, J., Liu, Y.Y., Wu, X.P., 2017. 3-D DC resistivity modeling with arbitrary long electrode sources using finite element method on unstructured grids. Geophys. J. Int. 211 (2), 1162–1176. https://doi.org/10.1093/gji/ggx356. Ye, M., Khaleel, R., Yeh, T.-C.J., 2005. Stochastic analysis of moisture plume dynamics of a field injection experiment. Water Resour. Res. 41, W03013. https://doi.org/10. 1029/2004WR003735. Yeh, T.-C., Khaleel, R., Carroll, K.C., 2015a. Flow through heterogeneous geologic media. Cambridge University Press. Yeh, T.C.J., Lee, C.H., Hsu, K.C., Wen, J.C., 2008. Fusion of hydrologic and geophysical tomographic surveys. Geosci. J. 12 (2), 159–167. https://doi.org/10.1007/s12303008-0017-6. Yeh, T.C.J., Liu, S., 2000. Hydraulic tomography: Development of a new aquifer test method. Water Resour. Res. 36 (8), 2095–2105. https://doi.org/10.1029/ 2000WR900114. Yeh, T.C.J., Liu, S., Glass, R.J., Baker, K., Brainard, J.R., Alumbaugh, D., LaBrecque, D., 2002. A geostatistically based inverse model for electrical resistivity surveys and its applications to vadose zone hydrology. Water Resour. Res. 38 (12), 1278. https://doi. org/10.1029/2001WR001204. Yeh, T.C.J., Mao, D.Q., Zha, Y.Y., Wen, J.C., Wan, L., Hsu, K.C., Lee, C.H., 2015b. Uniqueness, scale, and resolution issues in groundwater model parameter identification. Water Sci. Eng. 8 (3), 175–194. https://doi.org/10.1016/j.wse.2015.08.002. Yeh, T.-C.J., Gelhar, L.W., Gutjahr, A.L., 1985a. Stochastic Analysis of Unsaturated Flow in Heterogeneous Soils: 1. Statistically Isotropic Media. Water Resour. Res. 21 (4), 447–456. https://doi.org/10.1029/WR021i004p00447. Yeh, T.-C.J., Gelhar, L.W., Gutjahr, A.L., 1985b. Stochastic Analysis of Unsaturated Flow in Heterogeneous Soils: 2. Statistically Anisotropic Media With Variable α. Water Resour. Res. 21 (4), 457–464. https://doi.org/10.1029/WR021i004p00457. Yeh, T.-C.J., Gelhar, L.W., Gutjahr, A.L., 1985c. Stochastic Analysis of Unsaturated Flow in Heterogeneous Soils: 3. Observations and Applications. Water Resour. Res. 21 (4), 465–471. https://doi.org/10.1029/WR021i004p00465. Zha, Y., Yeh, T.-C.J., Illman, W.A., Onoe, H., Mok, C.M.W., Wen, J.-C., Huang, S.-Y., Wang, W., 2017. Incorporating geologic information into hydraulic tomography: a general framework based on geostatistical approach. Water Resour. Res. 53 (4), 2850–2876. https://doi.org/10.1002/2016WR019185. Zha, Y., Yeh, T.C.J., Illman, W.A., Zeng, W., Zhang, Y., Sun, F., Shi, L., 2018. A ReducedOrder Successive Linear Estimator for Geostatistical Inversion and its Application in Hydraulic Tomography. Water Resour. Res. 54 (3), 1616–1632. https://doi.org/10. 1002/2017WR021884. Zhao, Z., Illman, W.A., 2017. On the importance of geological data for three-dimensional steady-state hydraulic tomography analysis at a highly heterogeneous aquifer-aquitard system. Journal of Hydrology 544, 640–657. https://doi.org/10.1016/j.jhydrol. 2016.12.004. Zhu, J., Yeh, T.-C.J., 2005. Characterization of aquifer heterogeneity using transient hydraulic tomography. Water Resour. Res. 41 (7), W07028. https://doi.org/10. 1029/2004WR003790.

Environ. Sci. Technol. 50 (15), 8205–8213. Mao, D., Revil, A., Hort, R.D., Munakata-Marr, J., Atekwana, E.A., Kulessa, B., 2015. Resistivity and self-potential tomography applied to groundwater remediation and contaminant plumes: sandbox and field experiments. J. Hydrol. 530, 1–14. https:// doi.org/10.1016/j.jhydrol.2015.09.031. Mao, D., Yeh, T.C.J., Wan, L., Lee, C.H., Hsu, K.C., Wen, J.C., Lu, W., 2013. Cross-correlation analysis and information content of observed heads during pumping in unconfined aquifers. Water Resour. Res. 49 (2), 713–731. https://doi.org/10.1002/ wrcr.20066. Michael Tso, C.-H., Zha, Y., Yeh, T.C.J., Wen, J.-C., 2016. The relative importance of head, flux, and prior information in hydraulic tomography analysis. Water Resour. Res. 52, 3–20. https://doi.org/10.1002/2015WR017191. Miller, C.R., Routh, P.S., Brosten, T.R., McNamara, J.P., 2008. Application of time-lapse ERT imaging to watershed characterization. Geophysics 73 (3), G7–G17. https://doi. org/10.1190/1.2907156. Parsekian, A.D., Singha, K., Minsley, B.J., Holbrook, W.S., Slater, L., 2015. Multiscale geophysical imaging of the critical zone. Rev. Geophysics 53 (1), 1–26. https://doi. org/10.1002/2014rg000465. Pidlisecky, A., Haber, E., Knight, R., 2007. RESINVM3D: a 3D resistivity inversion package. Geophysics 72 (2), H1–H10. https://doi.org/10.1190/1.2402499. Pollock, D., Cirpka, O.A., 2008. Temporal moments in geoelectrical monitoring of salt tracer experiments. Water Resour. Res. 44 (12), W12416. https://doi.org/10.1029/ 2008WR007014. Pollock, D., Cirpka, O.A., 2010. Fully coupled hydrogeophysical inversion of synthetic salt tracer experiments. Water Resour. Res. 46 (7), W07501. https://doi.org/10.1029/ 2009WR008575. Pollock, D., Cirpka, O.A., 2012. Fully coupled hydrogeophysical inversion of a laboratory salt tracer experiment monitored by electrical resistivity tomography. Water Resour. Res. 48 (1), W01505. https://doi.org/10.1029/2011WR010779. Purvance, D.T., Andricevic, R., 2000a. On the electrical-hydraulic conductivity correlation in aquifers. Water Resour. Res. 36 (10), 2905–2913. https://doi.org/10.1029/ 2000WR900165. Purvance, D.T., Andricevic, R., 2000b. Geoelectric characterization of the hydraulic conductivity field and its spatial structure at variable scales. Water Resour. Res. 36 (10), 2915–2924. https://doi.org/10.1029/2000WR900187. Revil, A., Karaoulis, M., Johnson, T., Kemna, A., 2012. Review: Some low-frequency electrical methods for subsurface characterization and monitoring in hydrogeology. Hydrogeol. J. 20 (4), 617–658. https://doi.org/10.1007/s10040-011-0819-x. Robinson, J., Johnson, T., Slater, L., 2013. Evaluation of known-boundary and resistivity constraints for improving cross-borehole DC electrical resistivity imaging of discrete fractures. Geophysics 78 (3), D115–D127. https://doi.org/10.1190/geo2012-0333.1. Rucker, D.F., Loke, M.H., Levitt, M.T., Noonan, G.E., 2010. Electrical-resistivity characterization of an industrial site using long electrodes. Geophysics 75 (4), WA95WA104. https://doi.org/10.1190/1.3464806. Rucker, D.F., Fink, J.B., Loke, M.H., 2011. Environmental monitoring of leaks using timelapsed long electrode electrical resistivity. J. Appl. Geophys. 74 (4), 242–254. https://doi.org/10.1016/j.jappgeo.2011.06.005. Sasaki, Y., 1992. Resolution of resistivity tomography inferred from numerical simulation. Geophys. Prospect. 40, 453–463. https://doi.org/10.1111/j.1365-2478.1992. tb00536.x. Sentenac, P., Hogson, T., Keenan, H., Kulessa, B., 2015. Small-scale monitoring of a bioremediation barrier using miniature electrical resistivity tomography. J. Appl. Geophys. 115, 24–31. https://doi.org/10.1016/j.jappgeo.2014.11.006. Shima, H., 1992. 2-D and 3-D resistivity image reconstruction using crosshole data. Geophysics 57 (10), 1270–1281. https://doi.org/10.1190/1.1443195. Shima, H., Sakayama, T., 1987. Resistivity tomography: an approach to 2-D resistivity inverse problems: SEG Technical Program Expanded Abstracts. Society of Exploration Geophysicists 59–61. https://doi.org/10.1190/1.1892038. Singha, K., Day-Lewis, F.D., Johnson, T., Slater, L.D., 2015. Advances in the interpretation of subsurface processes with time-lapse electrical imaging. Hydrol. Process. 29 (6), 1549–1576. https://doi.org/10.1002/hyp.10280. Slater, L., 2007. Near-surface electrical characterization of hydraulic conductivity: from petrophysical properties to aquifer geometries - a review. Surv. Geophys. 28 (2–3), 169–197. https://doi.org/10.1007/s10712-007-9022-y. Slater, L., Binley, A., Versteeg, R., Cassiani, G., Birken, R., Sandberg, S., 2002. A 3D ERT study of solute transport in a large experimental tank. J. Appl. Geophys. 49 (4),

536