# Induced Polarization as a Tool to Assess Mineral Deposits: A Review

André Revil, Pierre Vaudelet, Zhaoyang Su, Rujun Chen

<figure><img src="image_1.png" /></figure>

## To cite this version:

André Revil, Pierre Vaudelet, Zhaoyang Su, Rujun Chen. Induced Polarization as a Tool to Assess Mineral Deposits: A Review. Minerals, 2022, 12 (5), pp.571. 10.3390/min12050571. hal-04261971

HAL Id: hal-04261971

https://hal.science/hal-04261971

Submitted on 28 Jun 2024

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_1.png" /></figure>

HAL
open science

# Induced Polarization as a Tool to Assess Mineral Deposits: A Review

André Revil, Pierre Vaudelet, Zhaoyang Su, Rujun Chen

► To cite this version:

André Revil, Pierre Vaudelet, Zhaoyang Su, Rujun Chen. Induced Polarization as a Tool to Assess Mineral Deposits: A Review. Minerals, 2022, 12 (5), pp.571. 10.3390/min12050571. hal-04261971

HAL Id: hal-04261971

https://hal.science/hal-04261971v1

Submitted on 28 Jun 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Review

# Induced Polarization as a Tool to Assess Mineral Deposits:
A Review

André Revil<sup>1,*</sup> ◆, Pierre Vaudelet<sup>2,●</sup>, Zhaoyang Su<sup>3,4</sup> ◆ and Rujun Chen<sup>5</sup>

**Abstract:** Disseminated ores in porous or fractured media can be polarized under the application of an external low-frequency electrical field. This polarization is characterized by a dimensionless property that is called the "chargeability". Induced polarization is a nonintrusive geophysical sensing technique that be used in the field to image both the electrical conductivity and the chargeability of porous rocks together with a characteristic relaxation time. A petrophysical model of the induced polarization of metallic ores immersed in a porous conductive and polarizable material is reviewed, and its predictions are compared to a large dataset of experimental data. The model shows that the chargeability of the material is linearly dependent on the volume fraction of the ore and the chargeability of the background material, which can, in turn, be related to the conductivity of the pore water and the cation exchange capacity of the clay fraction. The relaxation time depends on the grain sizes of the ores and on the conductivity of the background material, which is close to the conductivity of the porous rock itself. Five applications of the induced-polarization method to ore and metallic bodies are discussed in order to show the usefulness of this technique. These applications include: (i) A sandbox experiment, in which cubes of pyrite are located in a specific area of the tank; (ii) The tomography of an iron slag at an archeological site in France; (iii) A study of partially frozen graphitic schists in the French Alps; (iv) The detection of a metallic tank through the tomography of the relaxation times; and (v) The detection and localization of a deep ore body that is associated with a tectonic fault. We also discuss the possibility of combining self-potential and induced-polarization tomography to better characterize ore bodies below the seafloor.

**Keywords:** induced polarization; complex conductivity; ore body; exploration

<figure><img src="image_3.png" /></figure>

check for
updates

**Citation:** Revil, A.; Vaudelet, P.; Su, Z.; Chen, R. Induced Polarization as a Tool to Assess Mineral Deposits: A Review. *Minerals* **2022**, *12*, 571.
https://doi.org/10.3390/min12050571

Academic Editors: Mario Zarroca,
Roberto Rodriguez and Isaac Corral

Received: 17 March 2022

Accepted: 26 April 2022

Published: 30 April 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_4.png" /></figure>

**Copyright:** © 2022 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article
distributed under the terms and
conditions of the Creative Commons
Attribution (CC BY) license (https://
creativecommons.org/licenses/by/
4.0/).

# 1. Introduction

Induced polarization is a geophysical method that can be used to image the low-frequency polarization properties of rocks from the surface of the Earth [1,2] and in boreholes [3,4]. By “low-frequency”, we refer to frequencies below ~10 kHz. The involved polarization processes are non-dielectric in nature and they refer to the ability of porous media to store reversibly electrical charges, such as capacitance [2,5,6]. Induced polarization has a long history in geophysics, and it was first used in the realm of mining geophysics to localize ore bodies [2,6–16]. This method was also used for the localization of buried slag heaps from ancient smelters in archeology [17–20], and for the detection of pyrite-rich sealing anticlines that serve as oil traps in sedimentary basins [21–24]. Other areas of interest concern: (i) The detection of partially pyrite-filled tectonic fractures in clay-rich

<footnote n="1">¹ Université Grenoble Alpes, Université Savoie Mont-Blanc, CNRS, UMR CNRS 5204, EDYTEM, 73370 Le Bourget du Lac, France</footnote>

<footnote n="2">² Naga Geophysics, Technolac, 73370 Le Bourget du Lac, France; pierre.vaudelet@naga-geophysics.com</footnote>

<footnote n="3">³ State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum (Beijing), Beijing 102249, China; 2019310419@student.cup.edu.cn</footnote>

<footnote n="4">⁴ Key Laboratory of Submarine Geosciences, Second Institute of Oceanography, MNR, Hangzhou 310012, China</footnote>

<footnote n="5">⁵ School of Geosicence and Info-Physics, Central South University, Changsha 410083, China;
chen.rujun@foxmail.com</footnote>

<footnote>* Correspondence: andre.revil@univ-smb.fr</footnote>
Review

# Induced Polarization as a Tool to Assess Mineral Deposits:
A Review

André Revil<sup>1,*</sup> and Pierre Vaudelet<sup>2,3,4</sup> and Zhaoyang Su<sup>3,4</sup> and Rujun Chen<sup>5</sup>

**Abstract:** Disseminated ores in porous or fractured media can be polarized under the application of an external low-frequency electrical field. This polarization is characterized by a dimensionless property that is called the "chargeability". Induced polarization is a nonintrusive geophysical sensing technique that be used in the field to image both the electrical conductivity and the chargeability of porous rocks together with a characteristic relaxation time. A petrophysical model of the induced polarization of metallic ores immersed in a porous conductive and polarizable material is reviewed, and its predictions are compared to a large dataset of experimental data. The model shows that the chargeability of the material is linearly dependent on the volume fraction of the ore and the chargeability of the background material, which can, in turn, be related to the conductivity of the pore water and the cation exchange capacity of the clay fraction. The relaxation time depends on the grain sizes of the ores and on the conductivity of the background material, which is close to the conductivity of the porous rock itself. Five applications of the induced-polarization method to ore and metallic bodies are discussed in order to show the usefulness of this technique. These applications include: (i) A sandbox experiment, in which cubes of pyrite are located in a specific area of the tank; (ii) The tomography of an iron slag at an archeological site in France; (iii) A study of partially frozen graphitic schists in the French Alps; (iv) The detection of a metallic tank through the tomography of the relaxation times; and (v) The detection and localization of a deep ore body that is associated with a tectonic fault. We also discuss the possibility of combining self-potential and induced-polarization tomography to better characterize ore bodies below the seafloor.

**Keywords:** induced polarization; complex conductivity; ore body; exploration

## 1. Introduction

Induced polarization is a geophysical method that can be used to image the low-frequency polarization properties of rocks from the surface of the Earth [1,2] and in boreholes [3,4]. By "low-frequency", we refer to frequencies below ~10 kHz. The involved polarization processes are non-dielectric in nature and they refer to the ability of porous media to store reversibly electrical charges, such as capacitance [2,5,6]. Induced polarization has a long history in geophysics, and it was first used in the realm of mining geophysics to localize ore bodies [2,6–16]. This method was also used for the localization of buried slag heaps from ancient smelters in archeology [17–20], and for the detection of pyrite-rich sealing anticlines that serve as oil traps in sedimentary basins [21–24]. Other areas of interest concern: (i) The detection of partially pyrite-filled tectonic fractures in clay-rich

<figure><img src="image_3.png" /></figure>

check for
updates

**Citation:** Revil, A.; Vaudelet, P.; Su, Z.; Chen, R. Induced Polarization as a Tool to Assess Mineral Deposits: A Review. *Minerals* **2022**, *12*, 571.
https://doi.org/10.3390/min12050571

Academic Editors: Mario Zarroca, Roberto Rodriguez and Isaac Corral

Received: 17 March 2022

Accepted: 26 April 2022

Published: 30 April 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_4.png" /></figure>

**Copyright:** © 2022 by the authors.
Licensee MDPI, Basel, Switzerland.
This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

<footnote n="1">¹ Université Grenoble Alpes, Université Savoie Mont-Blanc, CNRS, UMR CNRS 5204, EDYTEM, 73370 Le Bourget du Lac, France</footnote>

<footnote n="2">² Naga Geophysics, Technolac, 73370 Le Bourget du Lac, France; pierre.vaudelet@naga-geophysics.com</footnote>

<footnote n="3">³ State Key Laboratory of Petroleum Resources and Prospecting, China University of Petroleum (Beijing), Beijing 102249, China; 2019310419@student.cup.edu.cn</footnote>

<footnote n="4">⁴ Key Laboratory of Submarine Geosciences, Second Institute of Oceanography, MNR, Hangzhou 310012, China</footnote>

<footnote n="5">⁵ School of Geosicence and Info-Physics, Central South University, Changsha 410083, China; chen.rujun@foxmail.com</footnote>

<footnote>* Correspondence: andre.revil@univ-smb.fr</footnote>
environments that are useable for the long-term storage of nuclear wastes [25]; (ii) The de-
tection of potentially reuseable mineral resources that are associated with abandoned mine
dumps [26]; and (iii) The monitoring of the remediation of contaminant plumes [27–30].
Induced polarization was also developed for non-mining applications in the realm of
hydrogeophysics [31,32] and to investigate trees [33], including the water uptake by the
roots and the dynamic of saps.

Modeling the induced-polarization signature of disseminated metallic particles in porous media has been performed with parametric empirical models, such as the broadly used Cole–Cole model [11,34,35]. Semi-empirical models [36] and mechanistic approaches that are based on redox reactions at the interface between the pore water and the surface of the metallic particles [37,38] have also been developed, as well as models that are based on the polarization of the electrical double layer (or the formation of field-induced double layers) at the solid–water interface [39,40]. The physics that are involved in the analysis of Wong (1979) [37] are, however, still being debated, and especially the role of redox active species (for instance, the couple Fe²⁺/Fe³⁺), in the overall polarization of semiconductors, such as pyrite or magnetite. The model of Wong (1979) [37], as well as the models of Bücker et al. (2018) [39,40], are unable to explain a very important observation that is related to the dependence of the relaxation time on the conductivity of the background material, or on the salinity or the conductivity of the rock itself [41,42], at least for metallic grains with sizes larger than a few hundred micrometers. This observation implies that the polarization mechanism for metallic grains must be associated with the accumulation of the charge carriers inside the metallic grains [41,43], and not around them (see also [44] for a specific experiment that confirms this standpoint).

In the present paper, we review a simple physics-based model to account for the pres-
ence of disseminated metallic particles in a porous material, which can be polarizable [44–47].
Numerical simulations that are performed with the finite element method are used to vali-
date the theory and to obtain better insights into the underlying physical mechanisms at
play, and especially to decipher the observed dependence of the relaxation time on the
conductivity of the background material (i.e., the material around the metallic grains [41,43]).
This petrophysical model is compared to the experimental data from the literature.

Then, we discuss various approaches that are used to invert the electrical conductivity, the chargeability, and the relaxation time. We briefly discuss five applications of the method to the detection of ores and metallic bodies. Finally, we discuss how induced polarization could complement the self-potential technique that has recently been used to explore sulfide deposits on the seafloor by using deep-towed geophysical instruments.

## 2. Time-Domain and Frequency-Domain Induced Polarization

Induced polarization was first described in the seminal book of Conrad Schlumberger (1920) [2]. It refers to the reversible storage of electrical charges in the ground under the application of a primary (low-frequency) electrical field or electrical current. This reversible storage of electrical charges occurs at some polarization length scales of the porous medium that correspond to the grain and pore scales [37,39,40,48]. In turn, this reversible storage of electrical charges is responsible for a secondary electrical current or electrical field that is associated with the electro-diffusion of the accumulated charge carriers through the gradients in the electrochemical potential of the charge carriers. The secondary electrical field or current is therefore “induced” by the primary current or electrical field. With that being said, the term “induction” is usually avoided because it refers to another specific mechanism in electromagnetic studies.

This secondary electrical field can be remotely measured by using a network of elec-
trodes. The relationship between the output (measured) secondary electrical field and the
applied signal (e.g., the primary current or electrical field) can be linear or nonlinear. In this
paper, we consider that the Earth is a time-invariant linear filter and that it can therefore be
analyzed by using linear time-invariant system theory. Such an approach reveals the transfer
function of porous media (with or without metallic minerals) in the frequency domain and
environments that are useable for the long-term storage of nuclear wastes [25]; (ii) The de-
tection of potentially reuseable mineral resources that are associated with abandoned mine
dumps [26]; and (iii) The monitoring of the remediation of contaminant plumes [27–30].
Induced polarization was also developed for non-mining applications in the realm of
hydrogeophysics [31,32] and to investigate trees [33], including the water uptake by the
roots and the dynamic of saps.

Modeling the induced-polarization signature of disseminated metallic particles in porous media has been performed with parametric empirical models, such as the broadly used Cole–Cole model [11,34,35]. Semi-empirical models [36] and mechanistic approaches that are based on redox reactions at the interface between the pore water and the surface of the metallic particles [37,38] have also been developed, as well as models that are based on the polarization of the electrical double layer (or the formation of field-induced double layers) at the solid–water interface [39,40]. The physics that are involved in the analysis of Wong (1979) [37] are, however, still being debated, and especially the role of redox active species (for instance, the couple Fe²⁺/Fe³⁺), in the overall polarization of semiconductors, such as pyrite or magnetite. The model of Wong (1979) [37], as well as the models of Bücker et al. (2018) [39,40], are unable to explain a very important observation that is related to the dependence of the relaxation time on the conductivity of the background material, or on the salinity or the conductivity of the rock itself [41,42], at least for metallic grains with sizes larger than a few hundred micrometers. This observation implies that the polarization mechanism for metallic grains must be associated with the accumulation of the charge carriers inside the metallic grains [41,43], and not around them (see also [44] for a specific experiment that confirms this standpoint).

In the present paper, we review a simple physics-based model to account for the pres-
ence of disseminated metallic particles in a porous material, which can be polarizable [44–47].
Numerical simulations that are performed with the finite element method are used to vali-
date the theory and to obtain better insights into the underlying physical mechanisms at
play, and especially to decipher the observed dependence of the relaxation time on the
conductivity of the background material (i.e., the material around the metallic grains [41,43]).
This petrophysical model is compared to the experimental data from the literature.

Then, we discuss various approaches that are used to invert the electrical conductivity, the chargeability, and the relaxation time. We briefly discuss five applications of the method to the detection of ores and metallic bodies. Finally, we discuss how induced polarization could complement the self-potential technique that has recently been used to explore sulfide deposits on the seafloor by using deep-towed geophysical instruments.

## 2. Time-Domain and Frequency-Domain Induced Polarization

Induced polarization was first described in the seminal book of Conrad Schlumberger (1920) [2]. It refers to the reversible storage of electrical charges in the ground under the application of a primary (low-frequency) electrical field or electrical current. This reversible storage of electrical charges occurs at some polarization length scales of the porous medium that correspond to the grain and pore scales [37,39,40,48]. In turn, this reversible storage of electrical charges is responsible for a secondary electrical current or electrical field that is associated with the electro-diffusion of the accumulated charge carriers through the gradients in the electrochemical potential of the charge carriers. The secondary electrical field or current is therefore “induced” by the primary current or electrical field. With that being said, the term “induction” is usually avoided because it refers to another specific mechanism in electromagnetic studies.

This secondary electrical field can be remotely measured by using a network of elec-
trodes. The relationship between the output (measured) secondary electrical field and the
applied signal (e.g., the primary current or electrical field) can be linear or nonlinear. In this
paper, we consider that the Earth is a time-invariant linear filter and that it can therefore be
analyzed by using linear time-invariant system theory. Such an approach reveals the transfer
function of porous media (with or without metallic minerals) in the frequency domain and
their impulse responses in the time domain [49]. Causality is a necessary constraint on the transfer function of the Earth for induced polarization. As an example, an RC circuit (R denotes the resistance and C the capacitance) falls into this category. For nonlinear induced-polarization effects, the interested readers are directed to specific publications, such as Hall and Olhoeft (1986) [49] and Chu and Bazant (2006) [50]. With that being said, the nonlinear polarization mechanisms have not been explored very much in geophysics.

The procedure that is used to measure induced polarization in the laboratory or in the field is similar to one of the classical electrical methods [2]. An electrical current is injected into the ground or through a core sample through a bipole of (current) electrodes (A and B). We observe the associated voltage or electrical potential difference at a pair of receiving (voltage) electrodes (M and N). Figures 1 and 2 illustrate the concept of time-domain induced polarization and chargeability. The existence of various definitions of the chargeability makes the literature, at first, a bit cryptic to read for newbies in the field. Usually, the measurements are performed in the time domain in the field, or in the frequency domain in the laboratory (Figure 3, see [51]). With that being said, it is worth noticing that more and more field investigations are based on frequency-domain induced polarization measurements [28,52–55]. In our viewpoint, it does not really matter if the measurements are performed in the time domain or the frequency domain. The date can even be acquired in the time domain and can be inverted in the frequency domain, or vice versa.

<figure><img src="image_5.png" /><figcaption>**Figure 1.** Time-domain induced polarization. (a) The electrode array includes two current electrodes (A and B) and two voltage electrodes (M and N) connected to a high-impedance (>100 MOhm) voltmeter. (b) Injected current ($T_{on}$ is the period). (c) Recorded difference of potential between M and N. The buildup of the voltage and the decay of the secondary voltage are symmetrical. We distinguish an initial voltage buildup ($\psi_{\infty}$) (instantaneous potential), which defines the instantaneous conductivity ($\sigma_{\infty}$), and a steady-state or direct-current (DC) voltage buildup ($\psi_0$), which defines the DC conductivity ($\sigma_0$). (d) The definition of three chargeabilities are given: a total chargeability ($M$) (dimensionless), and a partial chargeability ($M_{12}$) (in s), which is defined for a given time window during the decay of the secondary voltage. $M = M_{12}$ when the times are small with respect to the relaxation time ($\tau$). The time from 0 to $t_0$ is called the “dead-time”. Finally, an “integral chargeability” can be obtained by integrating the decay voltage function.</figcaption></figure>
their impulse responses in the time domain [49]. Causality is a necessary constraint on the transfer function of the Earth for induced polarization. As an example, an RC circuit (R denotes the resistance and C the capacitance) falls into this category. For nonlinear induced-polarization effects, the interested readers are directed to specific publications, such as Hall and Olhoeft (1986) [49] and Chu and Bazant (2006) [50]. With that being said, the nonlinear polarization mechanisms have not been explored very much in geophysics.

The procedure that is used to measure induced polarization in the laboratory or in the field is similar to one of the classical electrical methods [2]. An electrical current is injected into the ground or through a core sample through a bipole of (current) electrodes (A and B). We observe the associated voltage or electrical potential difference at a pair of receiving (voltage) electrodes (M and N). Figures 1 and 2 illustrate the concept of time-domain induced polarization and chargeability. The existence of various definitions of the chargeability makes the literature, at first, a bit cryptic to read for newbies in the field. Usually, the measurements are performed in the time domain in the field, or in the frequency domain in the laboratory (Figure 3, see [51]). With that being said, it is worth noticing that more and more field investigations are based on frequency-domain induced polarization measurements [28,52–55]. In our viewpoint, it does not really matter if the measurements are performed in the time domain or the frequency domain. The date can even be acquired in the time domain and can be inverted in the frequency domain, or vice versa.

<figure><img src="image_5.png" /><figcaption>**Figure 1.** Time-domain induced polarization. (a) The electrode array includes two current electrodes (A and B) and two voltage electrodes (M and N) connected to a high-impedance (>100 MOhm) voltmeter. (b) Injected current ($T_{on}$ is the period). (c) Recorded difference of potential between M and N. The buildup of the voltage and the decay of the secondary voltage are symmetrical. We distinguish an initial voltage buildup ($\psi_{\infty}$) (instantaneous potential), which defines the instantaneous conductivity ($\sigma_{\infty}$), and a steady-state or direct-current (DC) voltage buildup ($\psi_{0}$), which defines the DC conductivity ($\sigma_{0}$). (d) The definition of three chargeabilities are given: a total chargeability ($M$) (dimensionless), and a partial chargeability ($M_{12}$) (in s), which is defined for a given time window during the decay of the secondary voltage. $M = M_{12}$ when the times are small with respect to the relaxation time ($\tau$). The time from 0 to $t_{0}$ is called the “dead-time”. Finally, an “integral chargeability” can be obtained by integrating the decay voltage function.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Laboratory experiment showing the secondary voltage decay for three periods of the injected current. (a) Difference in potential recorded between electrodes M and N. The current is injected between A and B for a period ($T_{on}$). (b) Time-window-integrated chargeability for the apparent chargeability ($M_a$) (here expressed in mV/V). This decay illustrates the decay of the secondary voltage over time, which is associated with the charge carriers returning to their equilibrium positions by electro-diffusion.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 3.** Spectral induced polarization (SIP) involves a set of close harmonic signals for which the amplitudes and phases are reported on a log scale (for instance: 3 points/decade). The phase is generally considered as negative. The current and the voltage are written as *I* and *ψ*, respectively. Note that the phase is related to the time lag ($\Delta t$) by $\text{Phase}(\text{rad}) = \Delta t \frac{2\pi}{T}$.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>**Figure 2.** Laboratory experiment showing the secondary voltage decay for three periods of the injected current. (a) Difference in potential recorded between electrodes M and N. The current is injected between A and B for a period ($T_{on}$). (b) Time-window-integrated chargeability for the apparent chargeability ($M_a$) (here expressed in mV/V). This decay illustrates the decay of the secondary voltage over time, which is associated with the charge carriers returning to their equilibrium positions by electro-diffusion.</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>**Figure 3.** Spectral induced polarization (SIP) involves a set of close harmonic signals for which the amplitudes and phases are reported on a log scale (for instance: 3 points/decade). The phase is generally considered as negative. The current and the voltage are written as *I* and *ψ*, respectively. Note that the phase is related to the time lag ($\Delta t$) by $\text{Phase}(\text{rad}) = \Delta t \frac{2\pi}{T}$.</figcaption></figure>
In time-domain induced polarization, the process of polarization is characterized by a chargeability, while, in the frequency domain, the induced polarization is characterized by complex conductivity spectra and a phase. The phase and the chargeability can be related to each other in a manner that is similar to the normalized chargeability (product of the chargeability by the conductivity) and the quadrature conductivity (the phase times and the conductivity also correspond to the imaginary part of the complex conductivity (e.g., van Voorhis et al., 1973 [56]; Revil et al., 2017b [57]), in the context of the constant phase model). The transfer function corresponds to the filter that is associated with the ground (or with a core sample). As far as induced polarization is concerned, the most famous empirical or parametric "filters" are known as the Debye, the Cole-Cole, or the constant-phase model transfer functions [58]. We will come back to some of these filters in the present paper.

**3. Model in Absence of Background Polarization**

**3.1. Implication of the Maxwell–Clausius–Mossotti Equation**

We consider a mixture of metallic grains in a background, such as is shown in Figure 4. The grains are considered to be spheres, and the material is assumed to be isotropic at all scales. The metallic grains and the clay particles polarize with distinct polarization mechanisms and polarization levels (Figure 5; see Marshall and Madden 1959 [5] for a discussion on this topic). The typical spectra for such mixtures (e.g., pyrite or slag with sand) are shown in Figures 6 and 7.

<figure><img src="image_8.png" /><figcaption>**Figure 4.** Sketch of the complex background and the mixture between the background and the metallic particles (with metallic grains, for instance, pyrite). The background material can be partially saturated (with air as a nonwetting immiscible fluid) at a partial saturation. It can comprise various types of grains. Polarizable bacteria can be present in the pore space. The quantities $\sigma^*$ and $\sigma_b^*$ denote the complex conductivity of the mixture and the complex conductivity of the background material around the metallic particles, respectively.</figcaption></figure>
In time-domain induced polarization, the process of polarization is characterized by a chargeability, while, in the frequency domain, the induced polarization is characterized by complex conductivity spectra and a phase. The phase and the chargeability can be related to each other in a manner that is similar to the normalized chargeability (product of the chargeability by the conductivity) and the quadrature conductivity (the phase times and the conductivity also correspond to the imaginary part of the complex conductivity (e.g., van Voorhis et al., 1973 [56]; Revil et al., 2017b [57]), in the context of the constant phase model). The transfer function corresponds to the filter that is associated with the ground (or with a core sample). As far as induced polarization is concerned, the most famous empirical or parametric "filters" are known as the Debye, the Cole-Cole, or the constant-phase model transfer functions [58]. We will come back to some of these filters in the present paper.

**3. Model in Absence of Background Polarization**

**3.1. Implication of the Maxwell–Clausius–Mossotti Equation**

We consider a mixture of metallic grains in a background, such as is shown in Figure 4. The grains are considered to be spheres, and the material is assumed to be isotropic at all scales. The metallic grains and the clay particles polarize with distinct polarization mechanisms and polarization levels (Figure 5; see Marshall and Madden 1959 [5] for a discussion on this topic). The typical spectra for such mixtures (e.g., pyrite or slag with sand) are shown in Figures 6 and 7.

<figure><img src="image_8.png" /></figure>

**Figure 4.** Sketch of the complex background and the mixture between the background and the metallic particles (with metallic grains, for instance, pyrite). The background material can be partially saturated (with air as a nonwetting immiscible fluid) at a partial saturation. It can comprise various types of grains. Polarizable bacteria can be present in the pore space. The quantities $\sigma^*$ and $\sigma_b^*$ denote the complex conductivity of the mixture and the complex conductivity of the background material around the metallic particles, respectively.
<figure><img src="image_9.png" /><figcaption>**Figure 5.** Induced polarization of metallic and nonmetallic grains. (a) Nonmetallic grains such as clays are characterized by an electrical double layer that coats their surface. Under the influence of an external electrical field (E₀), the polarization of this electrical double layer provides a dipole moment to the grain. The main mechanism of polarization is the polarization of the Stern layer (SL), which is the inner part of the electrical double layer. (b) Polarization of a pyrite cube. The electro-diffusion of the charge carriers inside the grain (electrons and p-holes) polarizes the grain, which provides a dipole moment to the grain. A second polarization component is associated with field-induced diffuse layers (DL) in the electrolyte at the solid-liquid interface. The left-hand-side figures define the instantaneous (high-frequency) conductivity for the nonmetallic and metallic grains, while the right-hand-side figures define the DC (direct current, low-frequency) conductivity that corresponds to a long application of the electrical field.</figcaption></figure>

We first consider that the background material (in which the metallic particles are embedded) is a pure conductor and does not polarize. This assumption will be relaxed later on. We write the conductivity of the background as $\sigma_b$ (pure real number) instead of $\sigma_b^*$ (Figure 4). The quantity $\sigma_b$ depends on the porosity, the cation exchange capacity, and the pore water conductivity; an expression of the conductivity of the background material will be provided later (see Section 4.2 below and Winsauer and McCardell, 1953 [59]; Waxman and Smits, 1968 [60]; Vinegar and Waxman, 1984 [61]). We start here with the well-established form of the Maxwell–Clausius–Mossotti equation [37,62,63] to analyze the complex conductivity of the mixture ($\sigma^*$):
<figure><img src="image_9.png" /><figcaption>**Figure 5.** Induced polarization of metallic and nonmetallic grains. (a) Nonmetallic grains such as clays are characterized by an electrical double layer that coats their surface. Under the influence of an external electrical field (E₀), the polarization of this electrical double layer provides a dipole moment to the grain. The main mechanism of polarization is the polarization of the Stern layer (SL), which is the inner part of the electrical double layer. (b) Polarization of a pyrite cube. The electro-diffusion of the charge carriers inside the grain (electrons and p-holes) polarizes the grain, which provides a dipole moment to the grain. A second polarization component is associated with field-induced diffuse layers (DL) in the electrolyte at the solid-liquid interface. The left-hand-side figures define the instantaneous (high-frequency) conductivity for the nonmetallic and metallic grains, while the right-hand-side figures define the DC (direct current, low-frequency) conductivity that corresponds to a long application of the electrical field.</figcaption></figure>

We first consider that the background material (in which the metallic particles are embedded) is a pure conductor and does not polarize. This assumption will be relaxed later on. We write the conductivity of the background as σb (pure real number) instead of σb* (Figure 4). The quantity σb depends on the porosity, the cation exchange capacity, and the pore water conductivity; an expression of the conductivity of the background material will be provided later (see Section 4.2 below and Winsauer and McCardell, 1953 [59]; Waxman and Smits, 1968 [60]; Vinegar and Waxman, 1984 [61]). We start here with the well-established form of the Maxwell–Clausius–Mossotti equation [37,62,63] to analyze the complex conductivity of the mixture (σ*):
$$ \frac{\sigma^* - \sigma_b}{\sigma^* + 2\sigma_b} = \varphi_m \left( \frac{\sigma_m^* - \sigma_b}{\sigma_m^* + 2\sigma_b} \right) \qquad (1) $$

where $\sigma^*$ denotes the complex conductivity of the mixture (background with metallic particles); $\varphi_m$ denotes the volume fraction of the metallic particle; and $\sigma_m^*$ is the complex conductivity of the metallic particles [46,64]. The in-phase and the quadrature conductivities of the composite material correspond to the real and imaginary parts of $\sigma^*$ (i.e., $\sigma' \equiv \text{Re}(\sigma^*)$ and $\sigma'' \equiv \text{Im}(\sigma^*)$). Equation (1) is only valid for small concentrations of metallic particles (typically, $\varphi_m \le 0.16$ (e.g., Wong, 1979 [37]), but we will see that, actually, the model developed below works up to $\varphi_m \sim 0.22$).

<figure><img src="image_10.png" /><figcaption>**Figure 6.** Complex resistivity spectrum (amplitude and phase) for a core sample with a volume fraction of 6.3% of disseminated pyrite (data from Mahan et al., 1986 [65]). At the peak frequency, the phase reaches a high value of 100 mrad, which is much higher than the phases of rocks and soils in the absence of metallic particles. The chargeability (M) describes the relative amplitude of the frequency dependence of the resistivity with respect to the DC resistivity ($\rho_0$). The lines are guides for the eyes. Modified from [64].</figcaption></figure>
$$ \frac{\sigma^* - \sigma_b}{\sigma^* + 2\sigma_b} = \varphi_m \left( \frac{\sigma_m^* - \sigma_b}{\sigma_m^* + 2\sigma_b} \right) \qquad (1) $$

where $\sigma^*$ denotes the complex conductivity of the mixture (background with metallic particles); $\varphi_m$ denotes the volume fraction of the metallic particle; and $\sigma_m^*$ is the complex conductivity of the metallic particles [46,64]. The in-phase and the quadrature conductivities of the composite material correspond to the real and imaginary parts of $\sigma^*$ (i.e., $\sigma' \equiv \text{Re}(\sigma^*)$ and $\sigma'' \equiv \text{Im}(\sigma^*)$). Equation (1) is only valid for small concentrations of metallic particles (typically, $\varphi_m \le 0.16$ (e.g., Wong, 1979 [37]), but we will see that, actually, the model developed below works up to $\varphi_m \sim 0.22$).

<figure><img src="image_10.png" /><figcaption>**Figure 6.** Complex resistivity spectrum (amplitude and phase) for a core sample with a volume fraction of 6.3% of disseminated pyrite (data from Mahan et al., 1986 [65]). At the peak frequency, the phase reaches a high value of 100 mrad, which is much higher than the phases of rocks and soils in the absence of metallic particles. The chargeability (M) describes the relative amplitude of the frequency dependence of the resistivity with respect to the DC resistivity ($\rho_0$). The lines are guides for the eyes. Modified from [64].</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>**Figure 7.** Complex resistivity spectrum (amplitude and phase) for a sample from a slag heap [17]. The particle-size distribution of the metallic particles is broad, with a mean radius around 2 cm, which explains the occurrence of a low-frequency peak in the phase. The chargeability ($M = 0.66$ in the present case) is much higher than the chargeability of rocks and soils in the absence of metallic particles (usually smaller than 0.10). The lines are guides for the eyes. Modified from [64].</figcaption></figure>

In Equation (1), the quantity $f(i\omega) = (\sigma_m^* - \sigma_b) / (\sigma_m^* + 2\sigma_b)$ denotes the (complex-valued) “reflection coefficient” [62], which is connected to the dipole moment of the metallic grains. If the conductivity of the metallic particles is very large with respect to the background ($|\sigma_m^*| \gg \sigma_b$), we have $f(i\omega) = 1$. For insulating particles ($|\sigma_m^*| = 0$), we obtain $f(i\omega) = -1/2$.

In semiconductors such as pyrite and magnetite, the charge carriers inside the metallic particle are electrons and positive holes (Figures 5 and 8). At low frequencies and in the absence of redox active species in the pore water, these carriers are blocked by the surface of the metallic grain (Figure 5). This blocking effect (perfect capacitance) is responsible for the accumulation of the charge carriers at the boundary of the metallic particle in the applied electrical field. In our approach, this intraparticle mechanism is responsible for the observed polarization. The case where these charge carriers can be exchanged with the surrounding electrolyte (leaky capacitance) will not be discussed in the present paper (see [37]).
<figure><img src="image_11.png" /><figcaption>**Figure 7.** Complex resistivity spectrum (amplitude and phase) for a sample from a slag heap [17]. The particle-size distribution of the metallic particles is broad, with a mean radius around 2 cm, which explains the occurrence of a low-frequency peak in the phase. The chargeability ($M = 0.66$ in the present case) is much higher than the chargeability of rocks and soils in the absence of metallic particles (usually smaller than 0.10). The lines are guides for the eyes. Modified from [64].</figcaption></figure>

In Equation (1), the quantity $f(i\omega) = (\sigma_m^* - \sigma_b) / (\sigma_m^* + 2\sigma_b)$ denotes the (complex-valued) “reflection coefficient” [62], which is connected to the dipole moment of the metallic grains. If the conductivity of the metallic particles is very large with respect to the background ($|\sigma_m^*| \gg \sigma_b$), we have $f(i\omega) = 1$. For insulating particles ($|\sigma_m^*| = 0$), we obtain $f(i\omega) = -1/2$.

In semiconductors such as pyrite and magnetite, the charge carriers inside the metallic particle are electrons and positive holes (Figures 5 and 8). At low frequencies and in the absence of redox active species in the pore water, these carriers are blocked by the surface of the metallic grain (Figure 5). This blocking effect (perfect capacitance) is responsible for the accumulation of the charge carriers at the boundary of the metallic particle in the applied electrical field. In our approach, this intraparticle mechanism is responsible for the observed polarization. The case where these charge carriers can be exchanged with the surrounding electrolyte (leaky capacitance) will not be discussed in the present paper (see [37]).
<figure><img src="image_12.png" /><figcaption>**Figure 8.** Conduction band, valence band, and the energy gap in an insulator, a semiconductor, and a metal. The distinction between insulators and semiconductors resides in the value of the energy-band gap ($E_g$). The energy band denotes the region of discrete energy levels with very close spacings. Electrons and positive holes in the energy gap can be considered as low-density particle ensembles, such as in an ideal electrolyte. These energy levels can be treated as ions in an electrolyte, or as particles in a box (Fermi gas model). In an intrinsic semiconductor (no impurities), the concentrations of the charge carriers depend on the energy gap from the chemical potential ($\mu$), which itself corresponds to the Fermi level (in the middle between the conduction ($E_c$) and valence ($E_v$) energy bands). Modified from [64].</figcaption></figure>

If the applied (external) electrical field is removed, the n- and p-charge carriers will go back to their statistical equilibrium positions through electro-diffusion (the driving force for the motion is an electrochemical potential, see Section 3.2 below). We can conclude that the metallic grains would appear as insulators at low frequencies. Therefore, we have $f(i\omega) = -1/2$ at low-frequencies (i.e., for the direct current (DC) condition). Indeed, once a grain is fully polarized, all of its charge carriers are frozen as long as the electrical field remains constant. In some sense, we can say that polarization is detrimental to conduction. At high frequencies, the conductivity of the metallic grain is pretty large (much greater than the background material [66]). Therefore, at high frequencies, $f(i\omega) = 1$.
<figure><img src="image_12.png" /><figcaption>**Figure 8.** Conduction band, valence band, and the energy gap in an insulator, a semiconductor, and a metal. The distinction between insulators and semiconductors resides in the value of the energy-band gap ($E_g$). The energy band denotes the region of discrete energy levels with very close spacings. Electrons and positive holes in the energy gap can be considered as low-density particle ensembles, such as in an ideal electrolyte. These energy levels can be treated as ions in an electrolyte, or as particles in a box (Fermi gas model). In an intrinsic semiconductor (no impurities), the concentrations of the charge carriers depend on the energy gap from the chemical potential ($\mu$), which itself corresponds to the Fermi level (in the middle between the conduction ($E_c$) and valence ($E_v$) energy bands). Modified from [64].</figcaption></figure>

If the applied (external) electrical field is removed, the n- and p-charge carriers will go back to their statistical equilibrium positions through electro-diffusion (the driving force for the motion is an electrochemical potential, see Section 3.2 below). We can conclude that the metallic grains would appear as insulators at low frequencies. Therefore, we have $f(i\omega) = -1/2$ at low-frequencies (i.e., for the direct current (DC) condition). Indeed, once a grain is fully polarized, all of its charge carriers are frozen as long as the electrical field remains constant. In some sense, we can say that polarization is detrimental to conduction. At high frequencies, the conductivity of the metallic grain is pretty large (much greater than the background material [66]). Therefore, at high frequencies, $f(i\omega) = 1$.
The transition between the low-frequency and high-frequency domains is controlled by a relaxation time. From a dimensional analysis, this relaxation time depends on the radius of the grains and on the diffusion coefficient of the charge carriers in the metallic particles. If the grains have different sizes, this transition is controlled by a (probability) distribution of the relaxation times that can be related to a (probability) distribution of the polarization length scales [67]. Finding an expression of the relaxation time is therefore a key objective of any model of the polarization of porous media.

Equation (1) provides an explicit expression of the complex conductivity of the mixture:

$$ \sigma^* = \sigma_b \left( \frac{1 + 2\varphi_m f(i\omega)}{1 - \varphi_m f(i\omega)} \right) \qquad (2) $$

For small concentrations of metallic particles, we can approximate Equation (3) by:

$$ \sigma^* = \sigma_b [1 + 3\varphi_m f(i\omega) + \dots] \qquad (3) $$

by using a first-order Taylor development in $\varphi_m$.

If the grain-size distribution of the metallic particles is close to a log-normal distribution, the Cole–Cole model is known to provide an accurate representation of the induced polarization spectra of the disseminated ores [35]. Figures 6 and 7 are examples of spectra that could be fitted with a Cole–Cole model [35]. We therefore write the complex conductivity of the mixture by using a Cole–Cole model as [46]:

$$ \sigma^* = \sigma_\infty + \frac{\sigma_0 - \sigma_\infty}{1 + (i\omega\tau_0)^c} \qquad (4) $$

$$ \sigma^* = \sigma_\infty \left( 1 - \frac{M}{1 + (i\omega\tau_0)^c} \right) \qquad (5) $$

$$ M = \frac{\sigma_\infty - \sigma_0}{\sigma_\infty} = \frac{\rho_0 - \rho_\infty}{\rho_0} \qquad (6) $$

where M denotes the (total) chargeability of the material; $\rho_0 = 1/\sigma_0$ (DC resistivity) and $\rho_\infty = 1/\sigma_\infty$ (the M is called the “polarizability” in the Russian literature); and where $\tau_0$ denotes the relaxation time that is associated with the depolarization of a fully polarized grain. The chargeability describes the ratio of the polarization to the conduction. Hereinafter, we will also use the normalized chargeability (i.e., the product of the chargeability by the electrical conductivity ($M_n = M \sigma_\infty = \sigma_\infty - \sigma_0$), which better characterizes the polarization process alone). Our goal is to now provide a physical meaning for the parameters that are involved in Equation (5) or Equation (6) (i.e., $M$, $\sigma_\infty$, $c$, and $\tau_0$). The exponent $c$, which is called the Cole–Cole exponent, can be related to the broadness of the relaxation-time distribution [34,68,69].

With the Cole–Cole model, the reflection coefficient is given by $f(i\omega) = 1 - (3/2)/[1 + (i\omega\tau_0)^c]$. This can be compared with models on the basis of the polarization of the diffuse layer (see Bücker et al., 2018a [39], their Appendix C; and Wong, 1979 [37]): $f(i\omega) \approx 1 - (3/2)/[1 + i\omega\tau_1]$(diffuse layer polarization) with $\tau_1 = a\lambda_D/2D_{(\pm)}$ and $f(i\omega) \approx 1 - (3/2)/[1 + (i\omega\tau_2)^{1/2}]$ (volume-diffusion relaxation), with $\tau_2 = a^2/4D_{(\pm)}$, where $D_{(\pm)}$ is the diffusion coefficient of the ions in the pore water (assuming a binary symmetric of 1:1 electrolyte with the same diffusion coefficient for the cations and anions: $\sim 10^{-9} \text{ m}^2 \text{ s}^{-1}$); and $\lambda_D$ is the Debye length (in m) that is associated with the thickness of the diffuse layer around the grain.

Our goal now is to determine an expression of the relaxation time ($\tau_0$). From the dimensional analysis, the main relaxation time (for the complex conductivity, Equation (5)) is written as: $\tau_0 = a^2/D_a$ [70], where $D_a$ is the (apparent) diffusion coefficient of the charge carriers that are responsible for the polarization of the conductive particles. In Equation (5), $\sigma_0$ and $\sigma_\infty$ denote, respectively, the DC conductivity ($\omega \ll 1/\tau_0$) and the high-frequency conductivity of the porous material ($\omega \gg 1/\tau_0$). Olhoeft (1985) [6] determined that the
The transition between the low-frequency and high-frequency domains is controlled by a relaxation time. From a dimensional analysis, this relaxation time depends on the radius of the grains and on the diffusion coefficient of the charge carriers in the metallic particles. If the grains have different sizes, this transition is controlled by a (probability) distribution of the relaxation times that can be related to a (probability) distribution of the polarization length scales [67]. Finding an expression of the relaxation time is therefore a key objective of any model of the polarization of porous media.

Equation (1) provides an explicit expression of the complex conductivity of the mixture:

$$ \sigma^* = \sigma_b \left( \frac{1 + 2\varphi_m f(i\omega)}{1 - \varphi_m f(i\omega)} \right) \qquad (2) $$

For small concentrations of metallic particles, we can approximate Equation (3) by:

$$ \sigma^* = \sigma_b [1 + 3\varphi_m f(i\omega) + \dots] \qquad (3) $$

by using a first-order Taylor development in $\varphi_m$.

If the grain-size distribution of the metallic particles is close to a log-normal distribution, the Cole–Cole model is known to provide an accurate representation of the induced polarization spectra of the disseminated ores [35]. Figures 6 and 7 are examples of spectra that could be fitted with a Cole–Cole model [35]. We therefore write the complex conductivity of the mixture by using a Cole–Cole model as [46]:

$$ \sigma^* = \sigma_\infty + \frac{\sigma_0 - \sigma_\infty}{1 + (i\omega\tau_0)^c} \qquad (4) $$

$$ \sigma^* = \sigma_\infty \left( 1 - \frac{M}{1 + (i\omega\tau_0)^c} \right) \qquad (5) $$

$$ M = \frac{\sigma_\infty - \sigma_0}{\sigma_\infty} = \frac{\rho_0 - \rho_\infty}{\rho_0} \qquad (6) $$

where M denotes the (total) chargeability of the material; $\rho_0 = 1/\sigma_0$ (DC resistivity) and $\rho_\infty = 1/\sigma_\infty$ (the M is called the "polarizability" in the Russian literature); and where $\tau_0$ denotes the relaxation time that is associated with the depolarization of a fully polarized grain. The chargeability describes the ratio of the polarization to the conduction. Hereinafter, we will also use the normalized chargeability (i.e., the product of the chargeability by the electrical conductivity ($M_n = M \sigma_\infty = \sigma_\infty - \sigma_0$), which better characterizes the polarization process alone). Our goal is to now provide a physical meaning for the parameters that are involved in Equation (5) or Equation (6) (i.e., $M$, $\sigma_\infty$, $c$, and $\tau_0$). The exponent $c$, which is called the Cole–Cole exponent, can be related to the broadness of the relaxation-time distribution [34,68,69].

With the Cole–Cole model, the reflection coefficient is given by $f(i\omega) = 1 - (3/2)/[1 + (i\omega\tau_0)^c]$. This can be compared with models on the basis of the polarization of the diffuse layer (see Bücker et al., 2018a [39], their Appendix C; and Wong, 1979 [37]): $f(i\omega) \approx 1 - (3/2)/[1 + i\omega\tau_1]$(diffuse layer polarization) with $\tau_1 = a\lambda_D/2D_{(\pm)}$ and $f(i\omega) \approx 1 - (3/2)/[1 + (i\omega\tau_2)^{1/2}]$ (volume-diffusion relaxation), with $\tau_2 = a^2/4D_{(\pm)}$, where $D_{(\pm)}$ is the diffusion coefficient of the ions in the pore water (assuming a binary symmetric of 1:1 electrolyte with the same diffusion coefficient for the cations and anions: $\sim 10^{-9} \text{ m}^2 \text{ s}^{-1}$); and $\lambda_D$ is the Debye length (in m) that is associated with the thickness of the diffuse layer around the grain.

Our goal now is to determine an expression of the relaxation time ($\tau_0$). From the dimensional analysis, the main relaxation time (for the complex conductivity, Equation (5)) is written as: $\tau_0 = a^2/D_a$ [70], where $D_a$ is the (apparent) diffusion coefficient of the charge carriers that are responsible for the polarization of the conductive particles. In Equation (5), $\sigma_0$ and $\sigma_\infty$ denote, respectively, the DC conductivity ($\omega \ll 1/\tau_0$) and the high-frequency conductivity of the porous material ($\omega \gg 1/\tau_0$). Olhoeft (1985) [6] determined that the
relationship between the relaxation time of the pyrite-bearing porous media and the radius (a) of the metallic particles is $\tau_0 \approx 3 \times 10^5 a^2$, and that, therefore, using $\tau_0 = a^2/D_a$ leads to $D_a = 3.0 \times 10^{-5} \text{ m}^2 \text{ s}^{-1}$.

Tarasov and Titov (2013) [69] show that the Cole–Cole relaxation time ($\tau_0$) is related to the time constant that corresponds to the peak of the phase ($\tau$) by $\tau_0 = \tau(1 - M)^{1/2c} < \tau$. If a Cole–Cole model is used for the resistivity (see, for instance, the Cole–Cole resistivity model that is used in [11]), the relaxation time in the resistivity is given by $\tau_\rho = \tau_0(1 - M)^{-1/c} > \tau_0$ [18]. Note that, for small chargeabilities, $\tau_\rho \approx \tau_0(1 + M/c + ...)$, and, therefore, this time constant increases linearly with the chargeability for materials that are characterized by the same radius for the metallic particles [11]. An explicit derivation of the relaxation time is found below in Section 3.2.

We are now going to provide an explicit formula for this chargeability [46,64]. At a very high frequency ($f(i\omega) \to 1$), the metallic particles appear infinitively conductive [37] since the mobility of the charge carriers in the metallic particles is high and the electromigration of the charge carriers is not limited by the boundary of the metallic particle. To the contrary, at a very low frequency ($f(i\omega) \to -1/2$), the metallic particles appear insulating [37]. Indeed, at very low frequencies, all of the available charge carriers have reached the boundary of the metallic grain. The grain is fully polarized, and all of these charges are lost for the conduction process. From these considerations and by using Equation (4), we obtain:

$$ \sigma_{\infty} = \sigma_{b}(1 + 3\varphi_{m}) \qquad (7) $$

$$ \sigma_0 = \sigma_b \left( 1 - \frac{3}{2} \varphi_m \right) \qquad (8) $$

Equations (8) and (9) can also be derived by using the differential effective medium theory, as is further discussed below. From Equations (7)–(9), the chargeability is approximately given by:

$$ M \approx \frac{9}{2} \varphi_m \qquad (9) $$

If we start directly with the complete Maxwell–Clausius–Mossotti equation, without any approximation, we obtain (after few algebraic manipulations) the following (more general) expression for the chargeability for spherical metallic particles immersed in a conductive background [37,71]:

$$ M = \frac{9\varphi_m}{2 + 5\varphi_m + 2\varphi_m^2} \qquad (10) $$

Since $M \le 1$, Equation (10) is strictly valid only for $\varphi_m \le 2/9 \approx 0.22$ (22% vol. content of the mixture).

## 3.2. The Relaxation Time

The spectral induced polarization measurements from the literature have demonstrated the existence of a relationship between the relaxation time ($\tau_0$) and the diameter of the metallic grains (a) [11,41,46,64,70,72]. In order to go one step further in the understanding of the relaxation time (or the distribution of the relaxation times), numerical simulations can be performed to gain a good physical understanding of the underlying physics [41,43]. The local equations to solve correspond to the Poisson–Nernst–Planck equations that describe both the diffusion and the electromigration of the charge carriers [5,73], and both inside and outside of the metallic particles. We consider the case of a metallic particle that is immersed in a background of 1:1 electrolyte (e.g., NaCl, KCl). In the metallic particle, the charge carriers are holes and electrons (Figure 8), and in the pore fluid, the charge carriers are anions and cations.
relationship between the relaxation time of the pyrite-bearing porous media and the radius (a) of the metallic particles is $\tau_0 \approx 3 \times 10^5 a^2$, and that, therefore, using $\tau_0 = a^2/D_a$ leads to $D_a = 3.0 \times 10^{-5} \text{ m}^2 \text{ s}^{-1}$.

Tarasov and Titov (2013) [69] show that the Cole–Cole relaxation time ($\tau_0$) is related to the time constant that corresponds to the peak of the phase ($\tau$) by $\tau_0 = \tau(1 - M)^{1/2c} < \tau$. If a Cole–Cole model is used for the resistivity (see, for instance, the Cole–Cole resistivity model that is used in [11]), the relaxation time in the resistivity is given by $\tau_\rho = \tau_0(1 - M)^{-1/c} > \tau_0$ [18]. Note that, for small chargeabilities, $\tau_\rho \approx \tau_0(1 + M/c + ...)$, and, therefore, this time constant increases linearly with the chargeability for materials that are characterized by the same radius for the metallic particles [11]. An explicit derivation of the relaxation time is found below in Section 3.2.

We are now going to provide an explicit formula for this chargeability [46,64]. At a very high frequency ($f(i\omega) \to 1$), the metallic particles appear infinitively conductive [37] since the mobility of the charge carriers in the metallic particles is high and the electromigration of the charge carriers is not limited by the boundary of the metallic particle. To the contrary, at a very low frequency ($f(i\omega) \to -1/2$), the metallic particles appear insulating [37]. Indeed, at very low frequencies, all of the available charge carriers have reached the boundary of the metallic grain. The grain is fully polarized, and all of these charges are lost for the conduction process. From these considerations and by using Equation (4), we obtain:

$$ \sigma_\infty = \sigma_b(1 + 3\varphi_m) \qquad (7) $$

$$ \sigma_0 = \sigma_b \left( 1 - \frac{3}{2} \varphi_m \right) \qquad (8) $$

Equations (8) and (9) can also be derived by using the differential effective medium theory, as is further discussed below. From Equations (7)–(9), the chargeability is approximately given by:

$$ M \approx \frac{9}{2} \varphi_m \qquad (9) $$

If we start directly with the complete Maxwell–Clausius–Mossotti equation, without any approximation, we obtain (after few algebraic manipulations) the following (more general) expression for the chargeability for spherical metallic particles immersed in a conductive background [37,71]:

$$ M = \frac{9\varphi_m}{2 + 5\varphi_m + 2\varphi_m^2} \qquad (10) $$

Since $M \le 1$, Equation (10) is strictly valid only for $\varphi_m \le 2/9 \approx 0.22$ (22% vol. content of the mixture).

## 3.2. The Relaxation Time

The spectral induced polarization measurements from the literature have demonstrated the existence of a relationship between the relaxation time ($\tau_0$) and the diameter of the metallic grains (a) [11,41,46,64,70,72]. In order to go one step further in the understanding of the relaxation time (or the distribution of the relaxation times), numerical simulations can be performed to gain a good physical understanding of the underlying physics [41,43]. The local equations to solve correspond to the Poisson–Nernst–Planck equations that describe both the diffusion and the electromigration of the charge carriers [5,73], and both inside and outside of the metallic particles. We consider the case of a metallic particle that is immersed in a background of 1:1 electrolyte (e.g., NaCl, KCl). In the metallic particle, the charge carriers are holes and electrons (Figure 8), and in the pore fluid, the charge carriers are anions and cations.
We consider an applied harmonic electrical field ($E_a = E_0 \exp(i\omega t)$) with $E_0$ as the amplitude. The concentrations and the electrical potential are given by the sum of the values at equilibrium plus perturbations [39,40]:

$$C_{(\pm)}(\mathbf{r}, t) = C_{(\pm)}^{0} + \delta C_{(\pm)}(\mathbf{r}, t)e^{i\omega t} \quad (11)$$

$$\psi_e(\mathbf{r}, t) = \delta\psi_e(\mathbf{r}, t)e^{i\omega t} \quad (12)$$

$$C_{(+)}^{0} = C_{(-)}^{0} \quad (13)$$

where $\psi_e$ (V) describes the electrical potential that is external to the metallic particle; $C_{(\pm)}^{0}$ (m$^{-3}$) describe the equilibrium concentrations of the ions in the electrolyte in the absence of an external applied field; and $\delta C_{(\pm)}$ and $\delta\psi_e$ represent the perturbations of the concentration and the potential around the equilibrium, respectively. In the electrolyte, the driving force that controls the movement of the ionic charge carriers is the gradient in their electrochemical potential:

$$\nabla \tilde{\mu}_{(\pm)} = -(\pm e)\mathbf{E} + k_b T \nabla \ln C_{(\pm)} \quad (14)$$

where $e = 1.6 \times 10^{-19}$ C (elementary charge); $k_b = 1.381 \times 10^{-23}$ m$^2$ kg s$^{-2}$ K$^{-1}$ (Boltzmann constant); and $T$ is the absolute temperature (in K). The mobilities and diffusion coefficients are related to each other by the Nernst–Einstein relationship:

$$\beta_{(\pm)} = |q_{(\pm)}| b_{(\pm)} = \frac{|q_{(\pm)}| D_{(\pm)}}{k_b T} \quad (15)$$

where $q_{(\pm)} = \pm e$ denotes the charge of the cations and anions. The diffusion flux densities are given by the Nernst–Planck equations [74–76]:

$$\mathbf{J}_{(\pm)} = -b_{(\pm)} C_{(\pm)} \nabla \tilde{\mu}_{(\pm)} \quad (16)$$

$$\mathbf{J}_{(\pm)} = \left[ (\pm e) \frac{D_{(\pm)} C_{(\pm)}^{0}}{k_{b} T} \nabla \delta \psi_{e}(\mathbf{r}, t) - D_{(\pm)} \nabla \delta C_{(\pm)}(\mathbf{r}, t) \right] e^{i\omega t} \quad (17)$$

The continuity equation for the ionic species can be expressed as [37]:

$$\nabla \cdot \mathbf{J}_{(\pm)} = - \frac{\partial C_{(\pm)}(\mathbf{r}, t)}{\partial t} \quad (18)$$

Combining Equations (18) and (19) yields the following partial differential equation:

$$(\pm e) \frac{D_{(\pm)} C_{(\pm)}^{0}}{k_{b} T} \nabla^{2} \delta \psi_{e}(\mathbf{r}, t) - D_{(\pm)} \nabla^{2} \delta C_{(\pm)}(\mathbf{r}, t) = -i\omega \delta C_{(\pm)}(\mathbf{r}, t) \quad (19)$$

The electrical potential and the ionic concentrations are also related to each other by Gauss's law [77]:

$$\varepsilon_f \nabla \cdot \mathbf{E} = \rho \quad (20)$$

where $\rho$ (C m$^{-3}$) denotes the charge density, which is defined as:

$$\rho = e \left[ C_{(+)}(\mathbf{r}, t) - C_{(-)}(\mathbf{r}, t) \right] \quad (21)$$

Combining Equations (12), (13), (21), and (22) yields:

$$\nabla^2 \delta \psi_e(\mathbf{r}, t) = - \frac{Ne}{\varepsilon_f} \left[ \delta C_{(+)}(\mathbf{r}, t) - \delta C_{(-)}(\mathbf{r}, t) \right] \quad (22)$$
We consider an applied harmonic electrical field ($E_a = E_0 \exp(i\omega t)$) with $E_0$ as the amplitude. The concentrations and the electrical potential are given by the sum of the values at equilibrium plus perturbations [39,40]:

$$C_{(\pm)}(\mathbf{r}, t) = C_{(\pm)}^{0} + \delta C_{(\pm)}(\mathbf{r}, t)e^{i\omega t} \quad (11)$$

$$\psi_e(\mathbf{r}, t) = \delta\psi_e(\mathbf{r}, t)e^{i\omega t} \quad (12)$$

$$C_{(+)}^{0} = C_{(-)}^{0} \quad (13)$$

where $\psi_e$ (V) describes the electrical potential that is external to the metallic particle; $C_{(\pm)}^0$ (m⁻³) describe the equilibrium concentrations of the ions in the electrolyte in the absence of an external applied field; and $\delta C_{(\pm)}$ and $\delta\psi_e$ represent the perturbations of the concentration and the potential around the equilibrium, respectively. In the electrolyte, the driving force that controls the movement of the ionic charge carriers is the gradient in their electrochemical potential:

$$\nabla \tilde{\mu}_{(\pm)} = -(\pm e)\mathbf{E} + k_b T \nabla \ln C_{(\pm)} \quad (14)$$

where $e = 1.6 \times 10^{-19}$ C (elementary charge); $k_b = 1.381 \times 10^{-23}$ m² kg s⁻² K⁻¹ (Boltzmann constant); and T is the absolute temperature (in K). The mobilities and diffusion coefficients are related to each other by the Nernst–Einstein relationship:

$$\beta_{(\pm)} = |q_{(\pm)}| b_{(\pm)} = \frac{|q_{(\pm)}| D_{(\pm)}}{k_b T} \quad (15)$$

where $q_{(\pm)} = \pm e$ denotes the charge of the cations and anions. The diffusion flux densities are given by the Nernst–Planck equations [74–76]:

$$\mathbf{J}_{(\pm)} = -b_{(\pm)} C_{(\pm)} \nabla \tilde{\mu}_{(\pm)} \quad (16)$$

$$\mathbf{J}_{(\pm)} = \left[ (\pm e) \frac{D_{(\pm)} C_{(\pm)}^{0}}{k_b T} \nabla \delta\psi_e(\mathbf{r}, t) - D_{(\pm)} \nabla \delta C_{(\pm)}(\mathbf{r}, t) \right] e^{i\omega t} \quad (17)$$

The continuity equation for the ionic species can be expressed as [37]:

$$\nabla \cdot \mathbf{J}_{(\pm)} = - \frac{\partial C_{(\pm)}(\mathbf{r}, t)}{\partial t} \quad (18)$$

Combining Equations (18) and (19) yields the following partial differential equation:

$$(\pm e) \frac{D_{(\pm)} C_{(\pm)}^{0}}{k_b T} \nabla^2 \delta\psi_e(\mathbf{r}, t) - D_{(\pm)} \nabla^2 \delta C_{(\pm)}(\mathbf{r}, t) = -i\omega\delta C_{(\pm)}(\mathbf{r}, t) \quad (19)$$

The electrical potential and the ionic concentrations are also related to each other by Gauss's law [77]:

$$\varepsilon_f \nabla \cdot \mathbf{E} = \rho \quad (20)$$

where $\rho$ (C m⁻³) denotes the charge density, which is defined as:

$$\rho = e \left[ C_{(+)}(\mathbf{r}, t) - C_{(-)}(\mathbf{r}, t) \right] \quad (21)$$

Combining Equations (12), (13), (21), and (22) yields:

$$\nabla^2 \delta\psi_e(\mathbf{r}, t) = - \frac{Ne}{\varepsilon_f} \left[ \delta C_{(+)}(\mathbf{r}, t) - \delta C_{(-)}(\mathbf{r}, t) \right] \quad (22)$$
where N denotes the Avogadro number (6.02 × 10²³ mol⁻¹). The same type of equations can be written inside of the metallic particle for the two charge carriers: the electrons and the p-hole (especially for semiconductors), which also obeys the same partial differential equation that is derived from the Poisson–Nernst–Planck equations, the continuity equations, and Gauss's law [78,79]:

$$ (\pm e) \frac{D_{(\pm)}^{m} C_{(\pm)}^{m}}{k_{b} T} \nabla^{2} \delta \psi_{i}(\mathbf{r}, t) - D_{(\pm)}^{m} \nabla^{2} \delta C_{(\pm)}^{m}(\mathbf{r}, t) = -i \omega \delta C_{(\pm)}^{m}(\mathbf{r}, t) \quad (23) $$

$$ \nabla^2 \delta\psi_i(\mathbf{r}, t) = -\frac{Ne}{\varepsilon_S} \left[ \delta C_{(+)}^m(\mathbf{r}, t) - \delta C_{(-)}^m(\mathbf{r}, t) \right] \quad (24) $$

In these equations, $D_{(\pm)}^m$ denote the diffusion coefficients of the positive charge carriers (p-holes) and the negative charge carriers (electrons) in the semiconductor; $\varepsilon_S$ denotes the dielectric constant of the solid metallic grain; and $\delta\psi_i$ denotes the perturbed potential inside the metallic particle.

For the numerical simulations, we consider a pyrite grain (FeS₂; grain size is 2 µm; electrons and p-holes) that is immersed in a background electrolyte with KCl (K⁺ and Cl⁻), with $\sigma_S^\infty$(pyrite) = 1000 S m⁻¹, $D_{e-} = D_{hole} = 2.9 \times 10^{-5}$ m² s⁻¹, $D_{K+} = D_{Cl-} = 1.7 \times 10^{-9}$ m² s⁻¹, and the relative dielectric constant for the background is 80, the relative dielectric constant for the solid grains (pyrite) is 10.9, and with a temperature of 293 K (20 °C), and various salinities in the range from $10^{-3}$ to 0.5 Mol m⁻³. The partial differential equations are solved in time by using a finite difference approximation that is based on Euler's method, whereas the spatial discretization is performed by using the finite element method (FEM) (see [41,43] for details).

The applied boundary conditions that were used in the numerical simulations are as follows (Figure 9): (1) The Dirichlet boundary condition to the right and left sides of the domain (a constant current source); (2) The Neumann condition at the top and the bottom of the model (where we assume that no flux in the charges is allowed (the normal derivative of the potential is null). This corresponds to the insulation of the boundary conditions at these frontiers; (3) The blocking boundary condition, which was also imposed on the electrolyte/grain interface ($\partial J_{(\pm)} / \partial \hat{n} = 0$ (where $\hat{n}$ is the unit vector on the surface of the grains)), and therefore we neglect here the redox processes that would be responsible for the leaking capacitance that is associated with the grain. The charge concentration and the charge mobility in each region (grain and electrolyte) are defined by using a space function, the ions are defined in the electrolyte, whereas the electrons and holes are defined in the metallic grain; (4) The initial condition (for a time-dependent problem); before the current injection (at $t < 0$), we imposed the following conditions: $\psi = 0$ everywhere in the medium; $c_{cation} = c_{anion}$ and homogenous in the solution; and $c_{electron} = c_{hole}$ and homogenous in the grain. The mesh that was used for the 2D simulations is shown in Figure 9.

Figure 10 shows a 2D map of the potential distribution that resulted from the solution of the Poisson–Nernst–Planck equations at two different times (left): directly after the beginning of the injection (0.03 µs), and (right) after 1.5 µs of the current injection. The particle is in the center of the model, and so the grain is obtaining constant potential after a few nanoseconds of the current injection. During the current injection, the charge accumulates at the interface grain/electrolyte until a critical equilibrium state is reached (1.5 µs). This forms a charge condensation at the interface between the metallic grain and the background electrolyte, which forms a dipole at the scale of the grain.
where N denotes the Avogadro number (6.02 × 10²³ mol⁻¹). The same type of equations can be written inside of the metallic particle for the two charge carriers: the electrons and the p-hole (especially for semiconductors), which also obeys the same partial differential equation that is derived from the Poisson–Nernst–Planck equations, the continuity equations, and Gauss's law [78,79]:

$$ (\pm e) \frac{D_{(\pm)}^{m} C_{(\pm)}^{m}}{k_{b} T} \nabla^{2} \delta \psi_{i}(\mathbf{r}, t) - D_{(\pm)}^{m} \nabla^{2} \delta C_{(\pm)}^{m}(\mathbf{r}, t) = -i \omega \delta C_{(\pm)}^{m}(\mathbf{r}, t) \quad (23) $$

$$ \nabla^2 \delta\psi_i(\mathbf{r}, t) = -\frac{Ne}{\epsilon_S} \left[ \delta C_{(+)}^m(\mathbf{r}, t) - \delta C_{(-)}^m(\mathbf{r}, t) \right] \quad (24) $$

In these equations, $D_{(\pm)}^m$ denote the diffusion coefficients of the positive charge carriers (p-holes) and the negative charge carriers (electrons) in the semiconductor; $\epsilon_S$ denotes the dielectric constant of the solid metallic grain; and $\delta\psi_i$ denotes the perturbed potential inside the metallic particle.

For the numerical simulations, we consider a pyrite grain (FeS₂; grain size is 2 µm; electrons and p-holes) that is immersed in a background electrolyte with KCl (K⁺ and Cl⁻), with $\sigma_S^\infty$(pyrite) = 1000 S m⁻¹, $D_{e-} = D_{hole} = 2.9 \times 10^{-5}$ m² s⁻¹, $D_{K+} = D_{Cl-} = 1.7 \times 10^{-9}$ m² s⁻¹, and the relative dielectric constant for the background is 80, the relative dielectric constant for the solid grains (pyrite) is 10.9, and with a temperature of 293 K (20 °C), and various salinities in the range from $10^{-3}$ to 0.5 Mol m⁻³. The partial differential equations are solved in time by using a finite difference approximation that is based on Euler's method, whereas the spatial discretization is performed by using the finite element method (FEM) (see [41,43] for details).

The applied boundary conditions that were used in the numerical simulations are as follows (Figure 9): (1) The Dirichlet boundary condition to the right and left sides of the domain (a constant current source); (2) The Neumann condition at the top and the bottom of the model (where we assume that no flux in the charges is allowed (the normal derivative of the potential is null). This corresponds to the insulation of the boundary conditions at these frontiers; (3) The blocking boundary condition, which was also imposed on the electrolyte/grain interface ($\partial J_{(\pm)} / \partial \hat{n} = 0$ (where $\hat{n}$ is the unit vector on the surface of the grains)), and therefore we neglect here the redox processes that would be responsible for the leaking capacitance that is associated with the grain. The charge concentration and the charge mobility in each region (grain and electrolyte) are defined by using a space function, the ions are defined in the electrolyte, whereas the electrons and holes are defined in the metallic grain; (4) The initial condition (for a time-dependent problem); before the current injection (at $t < 0$), we imposed the following conditions: $\psi = 0$ everywhere in the medium; $c_{cation} = c_{anion}$ and homogenous in the solution; and $c_{electron} = c_{hole}$ and homogenous in the grain. The mesh that was used for the 2D simulations is shown in Figure 9.

Figure 10 shows a 2D map of the potential distribution that resulted from the solution of the Poisson–Nernst–Planck equations at two different times (left): directly after the beginning of the injection (0.03 µs), and (right) after 1.5 µs of the current injection. The particle is in the center of the model, and so the grain is obtaining constant potential after a few nanoseconds of the current injection. During the current injection, the charge accumulates at the interface grain/electrolyte until a critical equilibrium state is reached (1.5 µs). This forms a charge condensation at the interface between the metallic grain and the background electrolyte, which forms a dipole at the scale of the grain.
<figure><img src="image_13.png" /><figcaption>**Figure 9.** Finite element mesh for the 2D numerical simulations together with the boundary conditions used in the numerical model. The potential is constant inside the metallic particle and is imposed at the external boundaries that correspond to the current electrodes. The potential is computed around the metallic grain in the background material. Modified from [45].</figcaption></figure>

When the current is stopped, the accumulated charges relax to their initial (equilibrium) positions, both outside and inside the metallic particle, which generates a secondary potential decaying in time. Figure 10 shows maps of the secondary potential distribution and its fluctuations over time. By using such numerical simulations, Abdulsamad et al. (2020) [43] were able to reproduce the dependence of the relaxation time on the background conductivity, which was controlled here by the salinity (Figures 11 and 12). The vertical shift between the time relaxation that resulted from the simulation and that which resulted from the measurements was caused by the difference in the grain sizes. The observed dependence of the relaxation time on the salinity could be explained by considering the metallic grain as an electrical dipole that is formed during the current injection. This electrical dipole impacts the electric potential and, consequently, the charge distribution in its vicinity. Figure 10 shows the variation in the charge concentration away from the metallic particle/electrolyte interface at the end of the current injection time. We notice that, as long as the electrolyte concentration is low, the influence of the electric dipole is farther from the grain/electrolyte interface, which provides a higher relaxation time (larger distance for ions to go back to the equilibrium state). These results are consistent with the model that was discussed in the previous sections.
<figure><img src="image_13.png" /><figcaption>**Figure 9.** Finite element mesh for the 2D numerical simulations together with the boundary conditions used in the numerical model. The potential is constant inside the metallic particle and is imposed at the external boundaries that correspond to the current electrodes. The potential is computed around the metallic grain in the background material. Modified from [45].</figcaption></figure>

When the current is stopped, the accumulated charges relax to their initial (equilibrium) positions, both outside and inside the metallic particle, which generates a secondary potential decaying in time. Figure 10 shows maps of the secondary potential distribution and its fluctuations over time. By using such numerical simulations, Abdulsamad et al. (2020) [43] were able to reproduce the dependence of the relaxation time on the background conductivity, which was controlled here by the salinity (Figures 11 and 12). The vertical shift between the time relaxation that resulted from the simulation and that which resulted from the measurements was caused by the difference in the grain sizes. The observed dependence of the relaxation time on the salinity could be explained by considering the metallic grain as an electrical dipole that is formed during the current injection. This electrical dipole impacts the electric potential and, consequently, the charge distribution in its vicinity. Figure 10 shows the variation in the charge concentration away from the metallic particle/electrolyte interface at the end of the current injection time. We notice that, as long as the electrolyte concentration is low, the influence of the electric dipole is farther from the grain/electrolyte interface, which provides a higher relaxation time (larger distance for ions to go back to the equilibrium state). These results are consistent with the model that was discussed in the previous sections.
<figure><img src="image_14.png" /><figcaption>**Figure 10.** Secondary potential created after current injection during 150 ns. Electrolyte is KCl with concentration of 0.1 mol/L. The electrical dipole formed inside the grain appears clearly just after shutting down the current injection. The electrodes M and N correspond to the two voltage electrodes on which the potential difference is “measured”. The volume of the metal is determined by the surface fraction of the metallic body inside a square passing by M and N. Therefore, increasing the distance from these electrodes from the center of the particle is similar to a decrease in the volume fraction of the metallic particle.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>**Figure 10.** Secondary potential created after current injection during 150 ns. Electrolyte is KCl with concentration of 0.1 mol/L. The electrical dipole formed inside the grain appears clearly just after shutting down the current injection. The electrodes M and N correspond to the two voltage electrodes on which the potential difference is “measured”. The volume of the metal is determined by the surface fraction of the metallic body inside a square passing by M and N. Therefore, increasing the distance from these electrodes from the center of the particle is similar to a decrease in the volume fraction of the metallic particle.</figcaption></figure>
For small grains (<100 µm), the induced-polarization response is dominated by field-induced diffuse layers outside of the metallic grains, and the model that was developed by Bücker et al. (2018a,b) [39,40] can be used to understand the process, as discussed above. For large grains (>500 µm at low salinity), we can neglect the formation of field-induced diffuse layers outside of the grains. For a grain with a radius of 1.5 cm, the relaxation time is in the order of a few seconds. The numerical simulations show that the relaxation time decreases with the increase in the conductivity of the background electrolyte, and that it increases with the increase in the metallic particle conductivity itself (Figure 12). This dependence is written as:

$$ \tau_0 = \frac{a^2}{D_m} \left( 1 + \frac{\sigma_S^\infty}{\sigma_b} \right) \approx a^2 \frac{\sigma_S^\infty}{D_m \sigma_b} \qquad (25) $$

where $\sigma_S^\infty$ denotes the instantaneous conductivity of the metallic grain, and $\sigma_b$ denotes the background conductivity. Thus, the apparent diffusion coefficient is given by $D_a = D_m\sigma_b/\sigma_S^\infty$. Since the conductivity of the grains depends on the diffusion coefficient of the charge carriers ($D_m$), this means that the $D_a$ is independent of the diffusion coefficient of the charge carriers in the metallic particle; however, it depends on the concentration of the charge carriers in these metallic particles [66,80]. This is perhaps a paradox of the proposed polarization process.

<figure><img src="image_15.png" /><figcaption>**Figure 11.** Dependence of the relaxation time on the conductivity of the background $\sigma_b = \sigma_w$ where $\sigma_w$ denotes the conductivity of the KCl electrolyte that is used as the background in the numerical simulations. The vertical shift between the $\tau$ measured and the $\tau$ simulated corresponds to the distinct grain sizes in the physical and numerical experiments. These experiments show that the relaxation time is inversely proportional to the conductivity of the background.</figcaption></figure>
For small grains (<100 µm), the induced-polarization response is dominated by field-induced diffuse layers outside of the metallic grains, and the model that was developed by Bücker et al. (2018a,b) [39,40] can be used to understand the process, as discussed above. For large grains (>500 µm at low salinity), we can neglect the formation of field-induced diffuse layers outside of the grains. For a grain with a radius of 1.5 cm, the relaxation time is in the order of a few seconds. The numerical simulations show that the relaxation time decreases with the increase in the conductivity of the background electrolyte, and that it increases with the increase in the metallic particle conductivity itself (Figure 12). This dependence is written as:

$$ \tau_0 = \frac{a^2}{D_m} \left( 1 + \frac{\sigma_S^\infty}{\sigma_b} \right) \approx a^2 \frac{\sigma_S^\infty}{D_m \sigma_b} \qquad (25) $$

where $\sigma_S^\infty$ denotes the instantaneous conductivity of the metallic grain, and $\sigma_b$ denotes the background conductivity. Thus, the apparent diffusion coefficient is given by $D_a = D_m\sigma_b/\sigma_S^\infty$. Since the conductivity of the grains depends on the diffusion coefficient of the charge carriers ($D_m$), this means that the $D_a$ is independent of the diffusion coefficient of the charge carriers in the metallic particle; however, it depends on the concentration of the charge carriers in these metallic particles [66,80]. This is perhaps a paradox of the proposed polarization process.

<figure><img src="image_15.png" /><figcaption>**Figure 11.** Dependence of the relaxation time on the conductivity of the background $\sigma_b = \sigma_w$ where $\sigma_w$ denotes the conductivity of the KCl electrolyte that is used as the background in the numerical simulations. The vertical shift between the $\tau$ measured and the $\tau$ simulated corresponds to the distinct grain sizes in the physical and numerical experiments. These experiments show that the relaxation time is inversely proportional to the conductivity of the background.</figcaption></figure>
<figure><img src="image_16.png" /></figure>

**Figure 12.** Finite element simulations (filled symbols) of the polarization time constant (τ) associated with a single metallic grain (diameter: 1 cm). The background corresponds to an electrolyte of conductivity (σ<sub>w</sub>). (a) Relaxation time. We note σ<sub>S</sub><sup>∞</sup> as the instantaneous conductivity of the metallic grain. We observe the relaxation time scales as (1 + σ<sub>S</sub><sup>∞</sup> / σ<sub>b</sub>) (solid line) (i.e., as the ratio between the instantaneous conductivity of the metallic grain and the conductivity of the background (σ<sub>b</sub> = σ<sub>w</sub>)). This model explains very well the observed dependence of the relaxation time on the salinity of the pore water. (b) Chargeability versus volume fraction of metallic particles for the 2D simulations. Modified from Abdulsamad et al. (2020) [43].

The model that is developed in Section 3.1 above shows a linear relationship between the chargeability and the volume fraction of the metallic particles in the medium [41,43,46,64,70]. From the numerical simulations, and by using the potential difference that is picked up between two points (M and N) (Figure 9), the chargeability (M) is calculated as the ratio between the secondary voltage (in V), which is measured just after stopping the current, and the primary voltage (in V), which is measured at the end of the current injection [81]. We will come back to this relationship below in Section 6.3. Figure 12b shows the results of the numerical simulations for the chargeability versus the volume of the metallic particles in the medium. The linear relationship that is evidenced by the numerical simulations is consistent with both the theory and the experimental data (see Figure 13).
<figure><img src="image_16.png" /><figcaption>**Figure 12.** Finite element simulations (filled symbols) of the polarization time constant ($\tau$) associated with a single metallic grain (diameter: 1 cm). The background corresponds to an electrolyte of conductivity ($\sigma_w$). (a) Relaxation time. We note $\sigma_S^\infty$ as the instantaneous conductivity of the metallic grain. We observe the relaxation time scales as $(1 + \sigma_S^\infty / \sigma_b)$ (solid line) (i.e., as the ratio between the instantaneous conductivity of the metallic grain and the conductivity of the background ($\sigma_b = \sigma_w$)). This model explains very well the observed dependence of the relaxation time on the salinity of the pore water. (b) Chargeability versus volume fraction of metallic particles for the 2D simulations. Modified from Abdulsamad et al. (2020) [43].</figcaption></figure>

The model that is developed in Section 3.1 above shows a linear relationship between the chargeability and the volume fraction of the metallic particles in the medium [41,43,46,64,70]. From the numerical simulations, and by using the potential difference that is picked up between two points (M and N) (Figure 9), the chargeability (M) is calculated as the ratio between the secondary voltage (in V), which is measured just after stopping the current, and the primary voltage (in V), which is measured at the end of the current injection [81]. We will come back to this relationship below in Section 6.3. Figure 12b shows the results of the numerical simulations for the chargeability versus the volume of the metallic particles in the medium. The linear relationship that is evidenced by the numerical simulations is consistent with both the theory and the experimental data (see Figure 13).
<figure><img src="image_17.png" /><figcaption>**Figure 13.** Total chargeability versus pyrite content. Comparison between the prediction of the linear relationship (plain line) and the experimental data from Gurin et al. (2013) [70] (filled symbols except the triangle), Pelton et al. (1978) [11] (crosses), and Mahan et al. (1986) [65] (triangles). In the first case, the ore is a mixture of (by volume) pyrrhotite (30%), pyrite (30%), magnetite (30%), and chalcopyrite (10%). The sample with the pyrite grains that are from 1 to 2 mm in radius is from Phillips (2010) (pyrite, sample K01, sericite quartz, 10% vol. pyrite) [82]. In the second case, the metallic grains correspond to pyrite. The chargeability is expressed as a function of the volume content of the metallic particles ($\varphi_m$). This relationship is independent of the grain radius (a) of the metallic particles, as predicted by the model. For this dataset, the background porous material is characterized by a low chargeability, and so we are in the case that is described in Section 3. Modified from [64].</figcaption></figure>

# 4. Model in Presence of Background Polarization

## 4.1. Background Polarization and Chargeability

Our goal in this section is to explore the effect of the chargeability of the background on the overall chargeability of the material. For this purpose, a slightly more general approach than the Maxwell–Clausius–Mossotti equation is preferred, such as the solution of the differential effective medium (DEM) theory [83]. For metallic particles that are suspended in a polarizable background (Figure 4), we have below the iso-conductivity point [83,84]:

$$ \sigma^* = \sigma_b^*(1 - \varphi_m)^{\frac{m_m}{1-m_m}} \left( \frac{\sigma_m^* - \sigma_b^*}{\sigma_m^* - \sigma^*} \right)^{\frac{m_m}{1-m_m}} \quad (26) $$

where $\sigma_b^*$ denotes the complex conductivity of the background rock, and $m_m$ (dimensionless) denotes a grain-shape exponent of the metallic particles [83]. If the metallic particles are much more conductive with respect to the background material (and the conductivity of the mixture), Equation (27) reduces to:

$$ \lim_{|\sigma_m^*| \gg |\sigma_b^*|} \sigma^* \approx \sigma_b^*(1 - \varphi_m)^{\frac{m_m}{1-m_m}} \quad (27) $$
<figure><img src="image_17.png" /><figcaption>**Figure 13.** Total chargeability versus pyrite content. Comparison between the prediction of the linear relationship (plain line) and the experimental data from Gurin et al. (2013) [70] (filled symbols except the triangle), Pelton et al. (1978) [11] (crosses), and Mahan et al. (1986) [65] (triangles). In the first case, the ore is a mixture of (by volume) pyrrhotite (30%), pyrite (30%), magnetite (30%), and chalcopyrite (10%). The sample with the pyrite grains that are from 1 to 2 mm in radius is from Phillips (2010) (pyrite, sample K01, sericite quartz, 10% vol. pyrite) [82]. In the second case, the metallic grains correspond to pyrite. The chargeability is expressed as a function of the volume content of the metallic particles ($\varphi_m$). This relationship is independent of the grain radius (a) of the metallic particles, as predicted by the model. For this dataset, the background porous material is characterized by a low chargeability, and so we are in the case that is described in Section 3. Modified from [64].</figcaption></figure>

# 4. Model in Presence of Background Polarization

## 4.1. Background Polarization and Chargeability

Our goal in this section is to explore the effect of the chargeability of the background on the overall chargeability of the material. For this purpose, a slightly more general approach than the Maxwell–Clausius–Mossotti equation is preferred, such as the solution of the differential effective medium (DEM) theory [83]. For metallic particles that are suspended in a polarizable background (Figure 4), we have below the iso-conductivity point [83,84]:

$$ \sigma^* = \sigma_b^*(1 - \varphi_m)^{\frac{m_m}{1-m_m}} \left( \frac{\sigma_m^* - \sigma_b^*}{\sigma_m^* - \sigma^*} \right)^{\frac{m_m}{1-m_m}} \quad (26) $$

where $\sigma_b^*$ denotes the complex conductivity of the background rock, and $m_m$ (dimensionless) denotes a grain-shape exponent of the metallic particles [83]. If the metallic particles are much more conductive with respect to the background material (and the conductivity of the mixture), Equation (27) reduces to:

$$ \lim_{|\sigma_m^*| \gg |\sigma_b^*|} \sigma^* \approx \sigma_b^*(1 - \varphi_m)^{\frac{m_m}{1-m_m}} \quad (27) $$
If the number of the metallic particles is small (in the case of a dilute suspension of metallic particles that are immersed in a background material), we have:

$$ \lim_{|\sigma_m^*| \gg |\sigma_b^*|} \sigma^* \approx \sigma_b^* \left( 1 + \frac{m_m}{m_m - 1} \varphi_m + \dots \right) \quad (28) $$

This case corresponds to the high-frequency limit that is discussed above in Section 3, for which we obtain an expression of the instantaneous conductivity of the mixture:

$$ \sigma_{\infty} \approx \sigma_{b}^{\infty} \left( 1 + \frac{m_{m}}{m_{m}-1} \varphi_{m} + \dots \right) \quad (29) $$

In this case, there is an increase in the instantaneous conductivity with the volumetric concentration of metallic particles. The second case corresponds to the case that the metallic particles appear insulating (i.e., at low frequencies). In this case, all of the charge carriers of the metallic particle have accumulated and they do not contribute anymore to the conductivity of the metallic particle. We have:

$$ \sigma^* = \sigma_b^*(1 - \varphi_m)^{m_m} \left( \frac{1 - \sigma_m^*/\sigma_b^*}{1 - \sigma_m^*/\sigma^*} \right)^{m_m} \quad (30) $$

$$ \sigma^* \approx \sigma_b^*(1 - m_m \varphi_m + \dots) \quad (31) $$

where Equation (32) is obtained from Equation (31) by using a first-order Taylor expansion in $\varphi_m$. In this case, there is a decrease in the DC conductivity with the volumetric concentration of metallic particles. Equation (32) can be used to determine the low-frequency (DC) limit of the complex conductivity:

$$ \sigma_0 \approx \sigma_b^0 (1 - m_m \varphi_m + \dots) \quad (32) $$

By using Equations (30) and (33) and the definition of the chargeability (see Equation (7)) yields, after a few algebraic manipulations, we obtain the following expression for the chargeability:

$$ M \approx 1 - \frac{\sigma_b^0}{\sigma_b^\infty} \left[ 1 - \left( \frac{m_m^2}{m_m - 1} \right) \varphi_m \right] \quad (33) $$

Two extreme situations can be discussed. First, in the absence of metallic particles, the chargeability is given by the chargeability of the background ($M_b$); indeed:

$$ \lim_{\varphi_m=0} M = 1 - \frac{\sigma_b^0}{\sigma_b^\infty} = M_b \quad (34) $$

The second case is the case of spherical ($m_m = 3/2$) metallic grains in a nonpolarizable background. We have:

$$ \lim_{M_b=0} M = \left( \frac{m_m^2}{m_m - 1} \right) \varphi_m = \frac{9}{2} \varphi_m = M_m \quad (35) $$

where $M_m$ is the chargeability that is associated with the metallic particles only (see Section 3 above, Equations (10) and (11)). From Equations (34)–(36), the chargeability of the material is given by:

$$ M = 1 - (1 - M_b)(1 - M_m) \quad (36) $$

An expression of the background chargeability will be given below in Equation (48). For small background and metallic chargeabilities ($M_b, M_m$) $\ll 1$, Equation (37) simplifies to [46,64]:

$$ M \approx M_b + M_m \quad (37) $$
If the number of the metallic particles is small (in the case of a dilute suspension of metallic particles that are immersed in a background material), we have:

$$ \lim_{|\sigma_m^*| \gg |\sigma_b^*|} \sigma^* \approx \sigma_b^* \left( 1 + \frac{m_m}{m_m - 1} \varphi_m + \dots \right) \quad (28) $$

This case corresponds to the high-frequency limit that is discussed above in Section 3, for which we obtain an expression of the instantaneous conductivity of the mixture:

$$ \sigma_{\infty} \approx \sigma_b^{\infty} \left( 1 + \frac{m_m}{m_m - 1} \varphi_m + \dots \right) \quad (29) $$

In this case, there is an increase in the instantaneous conductivity with the volumetric concentration of metallic particles. The second case corresponds to the case that the metallic particles appear insulating (i.e., at low frequencies). In this case, all of the charge carriers of the metallic particle have accumulated and they do not contribute anymore to the conductivity of the metallic particle. We have:

$$ \sigma^* = \sigma_b^*(1 - \varphi_m)^{m_m} \left( \frac{1 - \sigma_m^*/\sigma_b^*}{1 - \sigma_m^*/\sigma^*} \right)^{m_m} \quad (30) $$

$$ \sigma^* \approx \sigma_b^*(1 - m_m \varphi_m + \dots) \quad (31) $$

where Equation (32) is obtained from Equation (31) by using a first-order Taylor expansion in $\varphi_m$. In this case, there is a decrease in the DC conductivity with the volumetric concentration of metallic particles. Equation (32) can be used to determine the low-frequency (DC) limit of the complex conductivity:

$$ \sigma_0 \approx \sigma_b^0 (1 - m_m \varphi_m + \dots) \quad (32) $$

By using Equations (30) and (33) and the definition of the chargeability (see Equation (7)) yields, after a few algebraic manipulations, we obtain the following expression for the chargeability:

$$ M \approx 1 - \frac{\sigma_b^0}{\sigma_b^\infty} \left[ 1 - \left( \frac{m_m^2}{m_m - 1} \right) \varphi_m \right] \quad (33) $$

Two extreme situations can be discussed. First, in the absence of metallic particles, the chargeability is given by the chargeability of the background ($M_b$); indeed:

$$ \lim_{\varphi_m=0} M = 1 - \frac{\sigma_b^0}{\sigma_b^\infty} = M_b \quad (34) $$

The second case is the case of spherical ($m_m = 3/2$) metallic grains in a nonpolarizable background. We have:

$$ \lim_{M_b=0} M = \left( \frac{m_m^2}{m_m - 1} \right) \varphi_m = \frac{9}{2} \varphi_m = M_m \quad (35) $$

where $M_m$ is the chargeability that is associated with the metallic particles only (see Section 3 above, Equations (10) and (11)). From Equations (34)–(36), the chargeability of the material is given by:

$$ M = 1 - (1 - M_b)(1 - M_m) \quad (36) $$

An expression of the background chargeability will be given below in Equation (48). For small background and metallic chargeabilities ($M_b, M_m$) $\ll 1$, Equation (37) simplifies to [46,64]:

$$ M \approx M_b + M_m \quad (37) $$
In other words, the chargeabilities are cumulative as long as they are small. Equations (37) and (38) are very important equations for understanding the chargeability of complex porous media with various sources of polarization. We now need to generalize the Cole–Cole model for such media. This is performed in the next section.

## 4.2. A Generalized Mixture Model

We assume, in this section, a porous material with K-polarization processes, and each of them is described by a Cole–Cole transfer function. Generally speaking, this could be the polarization of the nonmetallic grains, the polarization of the metallic grains, and even the polarization of bacteria in the pore space for mining wastes, for instance (Figure 4). A generalization of the complex conductivity of the Cole–Cole model for K-induced polarization distributions can be written as:

$$ \sigma^* = \sigma_{\infty} \left( 1 - \sum_{k=1}^{K} \frac{M_k}{1 + (i\omega\tau_k)^{c_k}} \right) \qquad (38) $$

$$ M_k = \frac{\sigma_k^\infty - \sigma_k^0}{\sigma_k^\infty} \qquad (39) $$

$$ \sigma_0 = \sigma_{\infty} \left( 1 - \sum_{k=1}^{K} M_k \right) \qquad (40) $$

where $M_k$, $c_k$, and $\tau_k$ denote the chargeability, the Cole–Cole exponent, and the relaxation time for the process (k) for a total number of K polarization processes, respectively. In the case of mineralized rocks, we consider two polarization mechanisms, as discussed in the previous sections. For two polarization processes (metallic and non-metallic grains), we have:

$$ \sigma^* = \sigma_{\infty} \left( 1 - \frac{M_1}{1 + (i\omega\tau_1)^{c_1}} - \frac{M_2}{1 + (i\omega\tau_2)^{c_2}} \right) \qquad (41) $$

$$ \sigma_{\infty} = \sigma_{1}^{\infty} + \sigma_{2}^{\infty} \qquad (42) $$

$$ \sigma_0 = \sigma_{\infty}(1 - M_1 - M_2) \qquad (43) $$

where $M_1$ and $M_2$ are the chargeabilities (dimensionless), $c_1$ and $c_2$ are the two Cole–Cole exponents (dimensionless), and $\tau_1$ and $\tau_2$ are the relaxation times (expressed in s). Indexes 1 and 2 refer to two dispersion processes, respectively (metallic and nonmetallic grains). With these relationships, we have, $\sigma_0 \neq \sigma_1^0 + \sigma_2^0$ and $\sigma_1^\infty \neq \sigma_2^0$. If the total chargeability is defined by $M = (\sigma_\infty - \sigma_0)/\sigma_\infty$, we have:

$$ M = M_1 + M_2 \qquad (44) $$

and since $0 \le M \le 1$, we have $0 \le M_1 + M_2 \le 1$. In our model, one of the two mechanisms of polarization (for instance, Mechanism 1) corresponds to the polarization of metallic particles, while the second mechanism (therefore, Mechanism 2) corresponds to the polarization of the electrical double layer of the insulating grains that form the skeleton of the background material (Figure 5).

By using a volume-averaging approach, Revil (2013) [85] obtained the following expressions for the high- and low-frequency conductivities of the background material:

$$ \sigma_2^\infty = \frac{1}{F} \sigma_w + \left( \frac{1}{F\phi} \right) \rho_g B \text{ CEC} \qquad (45) $$

$$ \sigma_2^0 = \frac{1}{F} \sigma_w + \left( \frac{1}{F\phi} \right) \rho_g (B - \lambda) \text{ CEC} \qquad (46) $$

Equation (46) can be used for the conductivity of the background material in Section 3. In these equations, $\sigma_w$ (in S m⁻¹) denotes the pore water conductivity (which depends
In other words, the chargeabilities are cumulative as long as they are small. Equations (37) and (38) are very important equations for understanding the chargeability of complex porous media with various sources of polarization. We now need to generalize the Cole–Cole model for such media. This is performed in the next section.

## 4.2. A Generalized Mixture Model

We assume, in this section, a porous material with K-polarization processes, and each of them is described by a Cole–Cole transfer function. Generally speaking, this could be the polarization of the nonmetallic grains, the polarization of the metallic grains, and even the polarization of bacteria in the pore space for mining wastes, for instance (Figure 4). A generalization of the complex conductivity of the Cole–Cole model for K-induced polarization distributions can be written as:

$$ \sigma^* = \sigma_{\infty} \left( 1 - \sum_{k=1}^{K} \frac{M_k}{1 + (i\omega\tau_k)^{c_k}} \right) \qquad (38) $$

$$ M_k = \frac{\sigma_k^\infty - \sigma_k^0}{\sigma_k^\infty} \qquad (39) $$

$$ \sigma_0 = \sigma_{\infty} \left( 1 - \sum_{k=1}^{K} M_k \right) \qquad (40) $$

where $M_k$, $c_k$, and $\tau_k$ denote the chargeability, the Cole–Cole exponent, and the relaxation time for the process (k) for a total number of K polarization processes, respectively. In the case of mineralized rocks, we consider two polarization mechanisms, as discussed in the previous sections. For two polarization processes (metallic and non-metallic grains), we have:

$$ \sigma^* = \sigma_{\infty} \left( 1 - \frac{M_1}{1 + (i\omega\tau_1)^{c_1}} - \frac{M_2}{1 + (i\omega\tau_2)^{c_2}} \right) \qquad (41) $$

$$ \sigma_{\infty} = \sigma_{1}^{\infty} + \sigma_{2}^{\infty} \qquad (42) $$

$$ \sigma_0 = \sigma_{\infty}(1 - M_1 - M_2) \qquad (43) $$

where $M_1$ and $M_2$ are the chargeabilities (dimensionless), $c_1$ and $c_2$ are the two Cole–Cole exponents (dimensionless), and $\tau_1$ and $\tau_2$ are the relaxation times (expressed in s). Indexes 1 and 2 refer to two dispersion processes, respectively (metallic and nonmetallic grains). With these relationships, we have, $\sigma_0 \neq \sigma_1^0 + \sigma_2^0$ and $\sigma_1^\infty \neq \sigma_2^0$. If the total chargeability is defined by $M = (\sigma_\infty - \sigma_0)/\sigma_\infty$, we have:

$$ M = M_1 + M_2 \qquad (44) $$

and since $0 \le M \le 1$, we have $0 \le M_1 + M_2 \le 1$. In our model, one of the two mechanisms of polarization (for instance, Mechanism 1) corresponds to the polarization of metallic particles, while the second mechanism (therefore, Mechanism 2) corresponds to the polarization of the electrical double layer of the insulating grains that form the skeleton of the background material (Figure 5).

By using a volume-averaging approach, Revil (2013) [85] obtained the following expressions for the high- and low-frequency conductivities of the background material:

$$ \sigma_2^\infty = \frac{1}{F} \sigma_w + \left( \frac{1}{F\phi} \right) \rho_g B \text{ CEC} \qquad (45) $$

$$ \sigma_2^0 = \frac{1}{F} \sigma_w + \left( \frac{1}{F\phi} \right) \rho_g (B - \lambda) \text{ CEC} \qquad (46) $$

Equation (46) can be used for the conductivity of the background material in Section 3. In these equations, $\sigma_w$ (in S m$^{-1}$) denotes the pore water conductivity (which depends
on both salinity and temperature); F (dimensionless) is the intrinsic formation factor that is related to the connected porosity ($\phi$) (dimensionless) by the first Archie's law ($F = \phi^{-m}$, where m (dimensionless) is called the "first Archie exponent" or the "porosity exponent" [86]); m is the cementation exponent (the porosity exponent or first Archie's exponent) that is associated with the grain shape of the nonmetallic grains; $\rho_g$ is the grain density (in kg m$^{-3}$, usually $\rho_g \sim 2650$ kg m$^{-3}$); and CEC denotes the cation exchange capacity (in C kg$^{-1}$, and often expressed in meq/100 g with 1 meq/100 g = 963.20 C kg$^{-1}$). The CEC is mostly sensitive to the weight fraction of the clay minerals and the clay type (e.g., chlorite, kaolinite, illite, smectite). In Equations (46) and (47), the quantity B (expressed in m$^2$ s$^{-1}$ V$^{-1}$) denotes the apparent mobility of the counterions for the surface conduction, and the quantity $\lambda$ (expressed in m$^2$ s$^{-1}$ V$^{-1}$) denotes the apparent mobility of the counterions for the polarization that is associated with the quadrature conductivity (see Revil et al., 2017a [87], and the references therein). These two quantities are related to the mobility of the counterions in the Stern and diffuse layers, and to the fraction of the counterions in the Stern layer.

A dimensionless number (R) is also introduced with $R = \lambda/B$ (see Revil et al., 2017a [87], for further explanations). This number is introduced because it corresponds to the maximum value of the chargeability in the absence of metallic particles. From previous experimental works [88], we have: $B(\text{Na}^+, 25 \,^{\circ}\text{C}) = 3.1 \pm 0.3 \times 10^{-9} \text{ m}^{-2} \text{ s}^{-1} \text{ V}^{-1}$ and $\lambda(\text{Na}^+, 25 \,^{\circ}\text{C}) = 3.0 \pm 0.7 \times 10^{-10} \text{ m}^{-2} \text{ s}^{-1} \text{ V}^{-1}$. These two quantities were determined by using large datasets of rock samples (including hundreds of core samples). It follows that the dimensionless number (R) is typically around $0.10 \pm 0.02$ (independent of the temperature and saturation [45]).

From the previous system of equations, we obtain the following expression of the chargeability of the background:

$$M_2 = \frac{\rho_g \lambda \text{CEC}}{\phi \sigma_w + \rho_g B \text{CEC}} \quad (47)$$

This equation explicitly shows the dependence of the chargeability of the background ($M_2$) on the pore-water conductivity (and hence, the salinity) and the CEC. Because the dependence of $\lambda$, $B$, and $\sigma_w$ on the temperature are roughly the same, the $M_2$ is essentially independent of the temperature. We also have the property $M_2 \le R \sim 0.09 \pm 0.01$ ($M_2 \approx R$ at low salinities, or for materials that are characterized by a very high CEC). The Cole–Cole time constant ($\tau_2$) is associated with a characteristic pore or grain size (in m) [67,89]. The Cole–Cole exponent ($c_2$) measures the broadness of the pore- or grain-size distribution.

For a mixture of nonmetallic grains, pore water, and metallic grains (Figure 4), we already obtained $M = (9/2)\phi_m + M_2$, where $\phi_m$ is the volume fraction of metallic particles. The instantaneous and steady-state conductivities of the mixture are related to the instantaneous and steady-state conductivities of the background material, according to [46]:

$$\sigma_{\infty} = \sigma_2^{\infty} (1 + 3\phi_m + \dots) \quad (48)$$

$$\sigma_0 = \sigma_2^0 \left( 1 - \frac{3}{2} \phi_m + \dots \right) \quad (49)$$

The previous model yields:

$$M_1 = \frac{9}{2} \phi_m \quad (50)$$

$$\sigma_1^\infty = 3\phi_m\sigma_2^\infty \quad (51)$$

$$\sigma_1^0 = \sigma_1^\infty (1 - M_1) \quad (52)$$

The last parameter to obtain is the relaxation time ($\tau_1$) (in s), which is obtained by [45] $\tau_1 \approx a^2 e^2 C_m / (k_b T \sigma_2^\infty)$, where a denotes the diameter of the metallic particle (in m), $C_m$ denotes the concentration of the charge carriers in the solid metallic particle (in m$^{-3}$), and
on both salinity and temperature); F (dimensionless) is the intrinsic formation factor that is related to the connected porosity ($\phi$) (dimensionless) by the first Archie's law ($F = \phi^{-m}$, where m (dimensionless) is called the "first Archie exponent" or the "porosity exponent" [86]); m is the cementation exponent (the porosity exponent or first Archie's exponent) that is associated with the grain shape of the nonmetallic grains; $\rho_g$ is the grain density (in kg m$^{-3}$, usually $\rho_g \sim 2650$ kg m$^{-3}$); and CEC denotes the cation exchange capacity (in C kg$^{-1}$, and often expressed in meq/100 g with 1 meq/100 g = 963.20 C kg$^{-1}$). The CEC is mostly sensitive to the weight fraction of the clay minerals and the clay type (e.g., chlorite, kaolinite, illite, smectite). In Equations (46) and (47), the quantity B (expressed in m$^2$ s$^{-1}$ V$^{-1}$) denotes the apparent mobility of the counterions for the surface conduction, and the quantity $\lambda$ (expressed in m$^2$ s$^{-1}$ V$^{-1}$) denotes the apparent mobility of the counterions for the polarization that is associated with the quadrature conductivity (see Revil et al., 2017a [87], and the references therein). These two quantities are related to the mobility of the counterions in the Stern and diffuse layers, and to the fraction of the counterions in the Stern layer.

A dimensionless number (R) is also introduced with $R = \lambda/B$ (see Revil et al., 2017a [87], for further explanations). This number is introduced because it corresponds to the maximum value of the chargeability in the absence of metallic particles. From previous experimental works [88], we have: $B(\text{Na}^+, 25 \,^{\circ}\text{C}) = 3.1 \pm 0.3 \times 10^{-9} \text{ m}^{-2} \text{ s}^{-1} \text{ V}^{-1}$ and $\lambda(\text{Na}^+, 25 \,^{\circ}\text{C}) = 3.0 \pm 0.7 \times 10^{-10} \text{ m}^{-2} \text{ s}^{-1} \text{ V}^{-1}$. These two quantities were determined by using large datasets of rock samples (including hundreds of core samples). It follows that the dimensionless number (R) is typically around $0.10 \pm 0.02$ (independent of the temperature and saturation [45]).

From the previous system of equations, we obtain the following expression of the chargeability of the background:

$$M_2 = \frac{\rho_g \lambda \text{CEC}}{\phi \sigma_w + \rho_g B \text{CEC}} \quad (47)$$

This equation explicitly shows the dependence of the chargeability of the background ($M_2$) on the pore-water conductivity (and hence, the salinity) and the CEC. Because the dependence of $\lambda$, $B$, and $\sigma_w$ on the temperature are roughly the same, the $M_2$ is essentially independent of the temperature. We also have the property $M_2 \le R \sim 0.09 \pm 0.01$ ($M_2 \approx R$ at low salinities, or for materials that are characterized by a very high CEC). The Cole–Cole time constant ($\tau_2$) is associated with a characteristic pore or grain size (in m) [67,89]. The Cole–Cole exponent ($c_2$) measures the broadness of the pore- or grain-size distribution.

For a mixture of nonmetallic grains, pore water, and metallic grains (Figure 4), we already obtained $M = (9/2)\phi_m + M_2$, where $\phi_m$ is the volume fraction of metallic particles. The instantaneous and steady-state conductivities of the mixture are related to the instantaneous and steady-state conductivities of the background material, according to [46]:

$$\sigma_{\infty} = \sigma_2^{\infty} (1 + 3\phi_m + \dots) \quad (48)$$

$$\sigma_0 = \sigma_2^0 \left( 1 - \frac{3}{2} \phi_m + \dots \right) \quad (49)$$

The previous model yields:

$$M_1 = \frac{9}{2} \phi_m \quad (50)$$

$$\sigma_1^\infty = 3\phi_m\sigma_2^\infty \quad (51)$$

$$\sigma_1^0 = \sigma_1^\infty (1 - M_1) \quad (52)$$

The last parameter to obtain is the relaxation time ($\tau_1$) (in s), which is obtained by [45] $\tau_1 \approx a^2 e^2 C_m / (k_b T \sigma_2^\infty)$, where a denotes the diameter of the metallic particle (in m), $C_m$ denotes the concentration of the charge carriers in the solid metallic particle (in m$^{-3}$), and
$e$ denotes the elementary charge ($1.6 \times 10^{-19}$ C). The relaxation ($\tau_2$) can be related to the grain sizes of the nonmetallic grains.

**5. Comparison to Experimental Data**

Our goal is to compare the model to experimental data by respecting the assumptions that were used to build our model.

**5.1. Relaxation Time**

Our first goal is to test the relationship $\tau_0 = a^2/D_a$. The relaxation time for the complex conductivity that is characterized by a low chargeability ($<0.10$) can be plotted as a function of the radius of the metallic grains (Figures 14 and 15). The values of the apparent diffusion coefficient ($D_a$) that is found in Figure 15 ($2.9 \times 10^{-5} \text{ m}^2 \text{ s}^{-1}$ at 25 °C, consistent with Olhoeft, 1983 [13]) and Figure 16 ($4.0 \times 10^{-6} \text{ m}^2 \text{ s}^{-1}$ at 25 °C) are very high and are incompatible with the diffusion coefficient of the ions (typically $\sim 1 \times 10^{-9} \text{ m}^2 \text{ s}^{-1}$ at 25 °C). The data that are shown in Figures 14 and 15 indicate clearly that the relaxation time is proportional to the square of the grain radius ($\tau_0 \sim a^2$). This is true as long as the conductivity of the background material is constant (otherwise, see Figure 12).

<figure><img src="image_18.png" /></figure>

**Figure 14.** Time constant ($\tau_0 = 1/(2\pi f_{\text{peak}})$) for the complex conductivity model versus the radius of the metallic particles. The quantity ($f_{\text{peak}}$) denotes the critical frequency that corresponds to the maximum of the phase. The diffusion coefficient ($D$) is an apparent diffusion coefficient for the charge carriers responsible for the polarization process. Experimental data: Pelton et al. (1978) [11] (pyrite); Gurin et al. (2013) [70] (a mixture of pyrrhotite (v.30%), pyrite (v.30%), magnetite (v.30%), and chalcopyrite (v.10%)); Abdel Aal et al. (2014) [72] (pyrite and magnetite); Ostrander and Zonge [90] (1978) (pyrite); Vanhala (1997) [91] (pyrite from his Figure 16); De Witt (1979) [71] (pyrite grains); Phillips (2010) [82] (pyrite, sample K01, sericite quartz, 10% vol. pyrite). The smallest relaxation times are more scattered, possibly because of the influence of Maxwell-Wagner polarization at high frequencies. Background conductivity: $\sigma_b \sim 10^{-2}-10^{-3} \text{ S m}^{-1}$. Modified from Revil et al. (2015c) [47].
$e$ denotes the elementary charge (1.6 × 10⁻¹⁹ C). The relaxation ($\tau_2$) can be related to the grain sizes of the nonmetallic grains.

**5. Comparison to Experimental Data**

Our goal is to compare the model to experimental data by respecting the assumptions that were used to build our model.

**5.1. Relaxation Time**

Our first goal is to test the relationship $\tau_0 = a^2/D_a$. The relaxation time for the complex conductivity that is characterized by a low chargeability (<0.10) can be plotted as a function of the radius of the metallic grains (Figures 14 and 15). The values of the apparent diffusion coefficient ($D_a$) that is found in Figure 15 (2.9 × 10⁻⁵ m² s⁻¹ at 25 °C, consistent with Olhoeft, 1983 [13]) and Figure 16 (4.0 × 10⁻⁶ m² s⁻¹ at 25 °C) are very high and are incompatible with the diffusion coefficient of the ions (typically ~1 × 10⁻⁹ m² s⁻¹ at 25 °C). The data that are shown in Figures 14 and 15 indicate clearly that the relaxation time is proportional to the square of the grain radius ($\tau_0 \sim a^2$). This is true as long as the conductivity of the background material is constant (otherwise, see Figure 12).

<figure><img src="image_18.png" /></figure>

**Figure 14.** Time constant ($\tau_0 = 1/(2\pi f_{peak})$) for the complex conductivity model versus the radius of the metallic particles. The quantity ($f_{peak}$) denotes the critical frequency that corresponds to the maximum of the phase. The diffusion coefficient ($D$) is an apparent diffusion coefficient for the charge carriers responsible for the polarization process. Experimental data: Pelton et al. (1978) [11] (pyrite); Gurin et al. (2013) [70] (a mixture of pyrrhotite (v.30%), pyrite (v.30%), magnetite (v.30%), and chalcopyrite (v.10%)); Abdel Aal et al. (2014) [72] (pyrite and magnetite); Ostrander and Zonge [90] (1978) (pyrite); Vanhala (1997) [91] (pyrite from his Figure 16); De Witt (1979) [71] (pyrite grains); Phillips (2010) [82] (pyrite, sample K01, sericite quartz, 10% vol. pyrite). The smallest relaxation times are more scattered, possibly because of the influence of Maxwell-Wagner polarization at high frequencies. Background conductivity: $\sigma_b \sim 10^{-2}-10^{-3} \text{ S m}^{-1}$. Modified from Revil et al. (2015c) [47].
<figure><img src="image_19.png" /></figure>

**Figure 15.** Time constant ($\tau = 1/(2 \pi f_{peak})$) for the complex conductivity model versus the radius of the metallic particles for chalcopyrite and galena. The quantity ($f_{peak}$) denotes the critical frequency that corresponds to the maximum of the phase. Experimental data for chalcopyrite and galena: Ostrander and Zonge (1978) [90] (chalcopyrite); Wong (1979) [37] (chalcopyrite); Collett (1959) [92] (3% galena); Grisseman (1971) [93] (6.3% chalcopyrite). Background conductivity: $\sigma_b \sim 10^{-2}-10^{-3} \text{ S m}^{-1}$. Modified from [47].

Figure 16 shows a set of data for slag heap. Because we expect slags to be metallic electronic conductors (rather than semiconductors), we expect to see a higher diffusion coefficient for slag in comparison to pyrite and magnetite. With that being said, slag heaps also contain some small amounts of magnetite since they are characterized by some of the non-negligible values of their magnetic susceptibility [17]. Additional data are required to see if the trend that is shown in Figure 16 is general or not for non-glassy slags.
<figure><img src="image_19.png" /></figure>

**Figure 15.** Time constant ($\tau = 1/(2 \pi f_{\text{peak}})$) for the complex conductivity model versus the radius of the metallic particles for chalcopyrite and galena. The quantity ($f_{\text{peak}}$) denotes the critical frequency that corresponds to the maximum of the phase. Experimental data for chalcopyrite and galena: Ostrander and Zonge (1978) [90] (chalcopyrite); Wong (1979) [37] (chalcopyrite); Collett (1959) [92] (3% galena); Grisseman (1971) [93] (6.3% chalcopyrite). Background conductivity: $\sigma_b \sim 10^{-2}-10^{-3} \text{ S m}^{-1}$. Modified from [47].

Figure 16 shows a set of data for slag heap. Because we expect slags to be metallic electronic conductors (rather than semiconductors), we expect to see a higher diffusion coefficient for slag in comparison to pyrite and magnetite. With that being said, slag heaps also contain some small amounts of magnetite since they are characterized by some of the non-negligible values of their magnetic susceptibility [17]. Additional data are required to see if the trend that is shown in Figure 16 is general or not for non-glassy slags.
<figure><img src="image_20.png" /></figure>

**Figure 16.** Time constant ($\tau = 1/(2\pi f_{\text{peak}})$) for the complex conductivity model versus the radius of the metallic particles for iron slag, and $f_{\text{peak}}$ is the critical frequency that corresponds to the maximum of the phase. Experimental data from Nordsiek and Weller (2008) [20], Grissemann et al. (2000) [94], and Florsch et al. (2011) [17]. In the case of the data from Nordsiek and Weller (2008) [20], we use the relationship $\tau = \tau_{\rho}(1 - M)^{1/2c}$ to convert the Cole–Cole time constant for the complex resistivity ($\tau_{\rho}$) into the time constant ($\tau$) that corresponds to the peak of the phase. The diffusion coefficient ($D$) is in the order of $1 \times 10^{-3} \text{ m}^2 \text{ s}^{-1}$, according to the grey band that corresponds to the model. The picture shows a block of slag (vertical bar: 10 cm) from the Castel-Minier site (Aulus-les-Bains, France, Courtesy from F. Téreygeol).

## 5.2. Chargeability and Volume Content of Metallic Particles

Equations (10) and (11) imply that the chargeability is a simple and direct indicator of the volumetric amount of the metallic particle in the pore space of the material (e.g., $M = (9/2)\phi_m$, where $\phi_m$ denotes the volumetric content of the metallic particles). This relationship is independent of the radius of the metallic particles, and it is independent of the physico-chemistry of the pore water. This relationship is not affected much by the true shape of the metallic particles, as long as the particles are spheroids with an eccentricity that is smaller than 10 and as long as they are randomly oriented [47]. Figure 13
<figure><img src="image_20.png" /></figure>

**Figure 16.** Time constant ($\tau = 1/(2\pi f_{\text{peak}})$) for the complex conductivity model versus the radius of the metallic particles for iron slag, and $f_{\text{peak}}$ is the critical frequency that corresponds to the maximum of the phase. Experimental data from Nordsiek and Weller (2008) [20], Grissemann et al. (2000) [94], and Florsch et al. (2011) [17]. In the case of the data from Nordsiek and Weller (2008) [20], we use the relationship $\tau = \tau_{\rho}(1 - M)^{1/2c}$ to convert the Cole–Cole time constant for the complex resistivity ($\tau_{\rho}$) into the time constant ($\tau$) that corresponds to the peak of the phase. The diffusion coefficient ($D$) is in the order of $1 \times 10^{-3} \text{ m}^2 \text{ s}^{-1}$, according to the grey band that corresponds to the model. The picture shows a block of slag (vertical bar: 10 cm) from the Castel-Minier site (Aulus-les-Bains, France, Courtesy from F. Téreygeol).

## 5.2. Chargeability and Volume Content of Metallic Particles

Equations (10) and (11) imply that the chargeability is a simple and direct indicator of the volumetric amount of the metallic particle in the pore space of the material (e.g., $M = (9/2)\phi_m$, where $\phi_m$ denotes the volumetric content of the metallic particles). This relationship is independent of the radius of the metallic particles, and it is independent of the physico-chemistry of the pore water. This relationship is not affected much by the true shape of the metallic particles, as long as the particles are spheroids with an eccentricity that is smaller than 10 and as long as they are randomly oriented [47]. Figure 13
shows a comparison between the prediction of our model (plain line) and a variety of experimental data.

Figure 13 shows very clearly that the relationship between the chargeability and the volume content of metallic particles is mostly linear and is independent of the grain radius (a) of the metallic particles (at least for concentrations of semiconductors < 22% vol.). For metallic particle concentrations higher than 22%, the chargeability reaches a plateau that is close to 1. Figures 17 and 18 show the effect of the content of metallic particles on the resistivity and chargeability, with fairly good agreement between the prediction of the proposed model and the laboratory data. Finally, the proposed model also predicts that the chargeability is independent of the pore-water conductivity, which is in agreement with the observations that were made by Clavier et al. (1976) [21].

<figure><img src="image_21.png" /></figure>

**Figure 17.** Influence of the volumetric content of pyrite on the amplitude of the resistivity. (a) Influence of the volumetric pyrite content on the resistivity spectra. The resistivity is normalized by the resistivity of the background porous material. The plain lines are just guides for the eyes. (b) Influence of the volumetric pyrite content on the low-frequency and high-frequency electrical conductivity. The data are normalized to the conductivity of the background ($\sigma$), for which $\phi_m = 0$. Data are from De Witt (1979, p. 55) [71]. The theory predicts the following trends: $\sigma_\infty = \sigma_b(1 + 3\phi_m)$ and $\sigma_0 = \sigma_b(1 - 3\phi_m/2)$ for spherical grains. Modified from [47].
shows a comparison between the prediction of our model (plain line) and a variety of experimental data.

Figure 13 shows very clearly that the relationship between the chargeability and the volume content of metallic particles is mostly linear and is independent of the grain radius (a) of the metallic particles (at least for concentrations of semiconductors < 22% vol.). For metallic particle concentrations higher than 22%, the chargeability reaches a plateau that is close to 1. Figures 17 and 18 show the effect of the content of metallic particles on the resistivity and chargeability, with fairly good agreement between the prediction of the proposed model and the laboratory data. Finally, the proposed model also predicts that the chargeability is independent of the pore-water conductivity, which is in agreement with the observations that were made by Clavier et al. (1976) [21].

<figure><img src="image_21.png" /></figure>

**Figure 17.** Influence of the volumetric content of pyrite on the amplitude of the resistivity. (a) Influence of the volumetric pyrite content on the resistivity spectra. The resistivity is normalized by the resistivity of the background porous material. The plain lines are just guides for the eyes. (b) Influence of the volumetric pyrite content on the low-frequency and high-frequency electrical conductivity. The data are normalized to the conductivity of the background ($\sigma$), for which $\phi_m = 0$. Data are from De Witt (1979, p. 55) [71]. The theory predicts the following trends: $\sigma_\infty = \sigma_b(1 + 3\phi_m)$ and $\sigma_0 = \sigma_b(1 - 3\phi_m/2)$ for spherical grains. Modified from [47].
<figure><img src="image_22.png" /><figcaption>**Figure 18.** Effect of the presence of semiconductors (pyrite and magnetite) on the chargeability of volcanic materials. For altered rocks without semiconductors, the chargeability of the background reaches a maximum level, which is given by $M_b = R = \lambda/B = 0.10 \pm 0.02$ [95]. Literature data (pyrite and chalcopyrite) are from Pelton et al. (1978) [11], Mahan et al. (1986) [65], Phillips (2010) [82], and Gurin et al. (2015) [96]. The predicted Trends 1 and 2 are based on Equation (38) and use two distinct values for the chargeability for the background material ($M_b$). When the volume fraction of the metal is higher than 1%, it can be easily detected by the induced-polarization method.</figcaption></figure>

## 6. Forward and Inverse Modeling

### 6.1. *Classical Approach*

In time-domain induced polarization, the apparent resistivities are determined by using the electrical potential that is recorded during the injection of the primary current through the current electrodes (A, B). The voltage response ($\psi$), which is due to the known point current source and Sinks A and B (at positions $\mathbf{r}_{(\pm)}$, with current intensity $\pm I$) can be described by the following elliptic partial differential equation that is obtained by combining Ohm's law with a conservation equation for the charge in the quasistatic limit of the Maxwell equations:

$$ \nabla \cdot (\sigma_{\infty} \nabla \psi_{\infty}) = -I[\delta(\mathbf{r} - \mathbf{r}_{(+)}) - \delta(\mathbf{r} - \mathbf{r}_{(-)})] \quad (53) $$

where $\delta(\cdot)$ denotes the Dirac (delta) function (in m⁻³).
<figure><img src="image_22.png" /><figcaption>**Figure 18.** Effect of the presence of semiconductors (pyrite and magnetite) on the chargeability of volcanic materials. For altered rocks without semiconductors, the chargeability of the background reaches a maximum level, which is given by $M_b = R = \lambda/B = 0.10 \pm 0.02$ [95]. Literature data (pyrite and chalcopyrite) are from Pelton et al. (1978) [11], Mahan et al. (1986) [65], Phillips (2010) [82], and Gurin et al. (2015) [96]. The predicted Trends 1 and 2 are based on Equation (38) and use two distinct values for the chargeability for the background material ($M_b$). When the volume fraction of the metal is higher than 1%, it can be easily detected by the induced-polarization method.</figcaption></figure>

# 6. Forward and Inverse Modeling

## 6.1. *Classical Approach*

In time-domain induced polarization, the apparent resistivities are determined by using the electrical potential that is recorded during the injection of the primary current through the current electrodes (A, B). The voltage response ($\psi$), which is due to the known point current source and Sinks A and B (at positions $\mathbf{r}_{(\pm)}$, with current intensity $\pm I$) can be described by the following elliptic partial differential equation that is obtained by combining Ohm's law with a conservation equation for the charge in the quasistatic limit of the Maxwell equations:

$$ \nabla \cdot (\sigma_{\infty} \nabla \psi_{\infty}) = -I[\delta(\mathbf{r} - \mathbf{r}_{(+)}) - \delta(\mathbf{r} - \mathbf{r}_{(-)})] \quad (53) $$

where $\delta(\cdot)$ denotes the Dirac (delta) function (in m$^{-3}$).
The ground surface ($\partial\Omega$) being insulating, the normal component of the electrical field vanishes, and therefore:

$$ \left. \frac{\partial \psi_{\infty}}{\partial \hat{\mathbf{n}}} \right|_{\partial \Omega} = 0 \qquad (54) $$

where $\sigma$ is the local electrical conductivity (in S m$^{-1}$), and $\hat{\mathbf{n}}$ is the unit vector normal to $\partial\Omega$. The apparent resistivity data can be inverted first by using an iterative Gauss–Newton algorithm by minimizing a cost function that is composed of a data misfit term and a regularizer, as discussed below.

Once the apparent resistivity data have been inverted, we can use the effect of the chargeability on the apparent chargeability, as discussed by Seigel (1959) [81], to invert the chargeability. As discussed above, induced polarization can be described by a macrscopic parameter that is called “chargeability” ($M$) [0, 1] (unit in V/V but usually expressed in mV/V). When the porous material is chargeable, the steady-state conductivity of the ground decreases with respect to the instantaneous conductivity, according to $\sigma_0 = \sigma_{\infty}(1 - M)$ (the charges responsible for polarization do not contribute to conduction anymore). Equation (7) can be used to rerun the conductivity to obtain a new potential value ($\psi$) (Figure 1). We could define the apparent chargeability ($M_a$) as $M_a = (\psi_0 - \psi_{\infty})/\psi_0$. During the inversion, the log of the apparent chargeability is treated as the observed data that characterizes the polarization of the subsurface. The log of the chargeability is treated as the model parameter for each cell.

For the resistivity and chargeability problems, we can optimize an objective or cost function written as:

$$ P^{\lambda}(\mathbf{m}) = \|\mathbf{W}_d(\mathbf{d}(\mathbf{m}) - \mathbf{d}_{obs})\|^2 + \lambda \|\mathbf{W}_m\mathbf{m}\|^2 \qquad (55) $$

where the vector ($\mathbf{m}$) is the model vector by using the log of the resistivity (or conductivity) and the chargeability for each cell, which are used to discretize the subsurface in order to enforce positiveness. The matrix $\mathbf{W}_d$ denotes the observation-data-weighting matrix that contains the uncertainties in the data (apparent resistivities or apparent chargeabilities). We can apply a first-order gradient matrix for the regularization matrix ($\mathbf{W}_m$), or we can look for a smooth model by penalizing the roughness of the model vector. The regularization parameter ($\lambda$) is used to balance the influence between the data misfit and the regularizer in the objective function, and there are various ways to optimize this parameter during the iterative inversion process. The observation vector ($\mathbf{d}_{obs}$) corresponds to the observed logarithm of the apparent resistivities or chargeabilities.

Following the classical algorithm that was developed by Oldenburg and Li (1994) [97], we use a deterministic inversion method with Tikhonov regulation to obtain a tomogram of the resistivity or chargeability. The intrinsic resistivity and chargeability are iteratively updated with:

$$ \mathbf{m}_{k+1} = \mathbf{m}_k + \delta\mathbf{m} \qquad (56) $$

where the perturbation is given by:

$$ \delta\mathbf{m} = -[G^T(\mathbf{W}_d^T\mathbf{W}_d)G + \lambda\mathbf{W}_m^T\mathbf{W}_m]^{-1}[G^T(\mathbf{W}_d^T\mathbf{W}_d)(\mathbf{d}(\mathbf{m}_i) - \mathbf{d}_{obs}) + \lambda\mathbf{W}_m^T\mathbf{W}_m\mathbf{m}_k] \qquad (57) $$

The vector $\delta\mathbf{m}$ provides the update of the vector $\mathbf{m}$ at iteration $k$. In this formulation, the matrix ($\mathbf{G}$) denotes the sensitivity matrix, which is updated at each iteration. A benchmark of such an algorithm is provided in Figure 19. The upper two figures denote the true model (for resistivity and chargeability), while the two lower figures denote the inverted tomograms that are obtained by simulating an acquisition survey at the ground surface.
The ground surface (∂Ω) being insulating, the normal component of the electrical field vanishes, and therefore:

$$ \left. \frac{\partial \psi_{\infty}}{\partial \hat{\mathbf{n}}} \right|_{\partial \Omega} = 0 \qquad (54) $$

where σ is the local electrical conductivity (in S m⁻¹), and **n** is the unit vector normal to ∂Ω. The apparent resistivity data can be inverted first by using an iterative Gauss–Newton algorithm by minimizing a cost function that is composed of a data misfit term and a regularizer, as discussed below.

Once the apparent resistivity data have been inverted, we can use the effect of the chargeability on the apparent chargeability, as discussed by Seigel (1959) [81], to invert the chargeability. As discussed above, induced polarization can be described by a macrscopic parameter that is called “chargeability” (M) [0, 1] (unit in V/V but usually expressed in mV/V). When the porous material is chargeable, the steady-state conductivity of the ground decreases with respect to the instantaneous conductivity, according to σ₀ = σ∞(1 − M) (the charges responsible for polarization do not contribute to conduction anymore). Equation (7) can be used to rerun the conductivity to obtain a new potential value (ψ) (Figure 1). We could define the apparent chargeability (M_a) as M_a = (ψ₀ − ψ_∞)/ψ₀. During the inversion, the log of the apparent chargeability is treated as the observed data that characterizes the polarization of the subsurface. The log of the chargeability is treated as the model parameter for each cell.

For the resistivity and chargeability problems, we can optimize an objective or cost function written as:

$$ P^{\lambda}(\mathbf{m}) = \|\mathbf{W}_d(\mathbf{d}(\mathbf{m}) - \mathbf{d}_{\text{obs}})\|^2 + \lambda \|\mathbf{W}_m\mathbf{m}\|^2 \qquad (55) $$

where the vector (m) is the model vector by using the log of the resistivity (or conductivity) and the chargeability for each cell, which are used to discretize the subsurface in order to enforce positiveness. The matrix W_d denotes the observation-data-weighting matrix that contains the uncertainties in the data (apparent resistivities or apparent chargeabilities). We can apply a first-order gradient matrix for the regularization matrix (W_m), or we can look for a smooth model by penalizing the roughness of the model vector. The regularization parameter (λ) is used to balance the influence between the data misfit and the regularizer in the objective function, and there are various ways to optimize this parameter during the iterative inversion process. The observation vector (dobs) corresponds to the observed logarithm of the apparent resistivities or chargeabilities.

Following the classical algorithm that was developed by Oldenburg and Li (1994) [97], we use a deterministic inversion method with Tikhonov regulation to obtain a tomogram of the resistivity or chargeability. The intrinsic resistivity and chargeability are iteratively updated with:

$$ \mathbf{m}_{k+1} = \mathbf{m}_k + \delta\mathbf{m} \qquad (56) $$

where the perturbation is given by:

$$ \delta\mathbf{m} = -[G^T(\mathbf{W}_d^T\mathbf{W}_d)G + \lambda\mathbf{W}_m^T\mathbf{W}_m]^{-1}[G^T(\mathbf{W}_d^T\mathbf{W}_d)(\mathbf{d}(\mathbf{m}_i) - \mathbf{d}_{\text{obs}}) + \lambda\mathbf{W}_m^T\mathbf{W}_m\mathbf{m}_k] \qquad (57) $$

The vector δm provides the update of the vector m at iteration k. In this formulation, the matrix (G) denotes the sensitivity matrix, which is updated at each iteration. A benchmark of such an algorithm is provided in Figure 19. The upper two figures denote the true model (for resistivity and chargeability), while the two lower figures denote the inverted tomograms that are obtained by simulating an acquisition survey at the ground surface.
<figure><img src="image_23.png" /><figcaption>**Figure 19.** Deterministic induced polarization tomography test using a synthetic model. The simulation domain covers a volume of 18 m × 14 m × 6 m. A total of 6 profiles are used at the top of this volume to simulate a survey (15 electrodes separated by 1 m from each other along a profile, and 2 m between profiles). We use a dipole-dipole configuration (372 apparent resistivity data and 372 apparent chargeability data, and 5% white noise added to the data). (a) True conductivity distribution. (b) True chargeability distribution. (c) Inverted conductivity distribution. (d) Inverted chargeability distribution. Modified from [98].</figcaption></figure>

## 6.2. Source Current Approach

In time-domain induced polarization, the secondary voltage is associated with a secondary current, as discussed above in Section 3.2. The occurrence of the source of the current in the absence of an external applied electrical field is similar to a self-potential problem. Indeed, in geophysics, the self-potential method is a passive geophysical method in which the electrical field is passively measured at the surface of the Earth and is associated with the existence of current sources in the ground. After the shutdown of the primary current, a total current density (J) (in A m⁻²) can be described by a generalized Ohm's law:

$$J = \sigma E + J_s \qquad (58)$$

where E denotes the secondary electrical field (in V m⁻¹), and Js is the secondary (source) current density that is responsible for the observed electrical field (in A m⁻²). Since the time variation of the current density is very slow, we can use the quasistatic limit of the Maxwell equations, for which the conservation of the charge takes the form:

$$\nabla \cdot J = 0 \qquad (59)$$
<figure><img src="image_23.png" /><figcaption>**Figure 19.** Deterministic induced polarization tomography test using a synthetic model. The simulation domain covers a volume of 18 m × 14 m × 6 m. A total of 6 profiles are used at the top of this volume to simulate a survey (15 electrodes separated by 1 m from each other along a profile, and 2 m between profiles). We use a dipole-dipole configuration (372 apparent resistivity data and 372 apparent chargeability data, and 5% white noise added to the data). (a) True conductivity distribution. (b) True chargeability distribution. (c) Inverted conductivity distribution. (d) Inverted chargeability distribution. Modified from [98].</figcaption></figure>

## 6.2. Source Current Approach

In time-domain induced polarization, the secondary voltage is associated with a secondary current, as discussed above in Section 3.2. The occurrence of the source of the current in the absence of an external applied electrical field is similar to a self-potential problem. Indeed, in geophysics, the self-potential method is a passive geophysical method in which the electrical field is passively measured at the surface of the Earth and is associated with the existence of current sources in the ground. After the shutdown of the primary current, a total current density (J) (in A m⁻²) can be described by a generalized Ohm's law:

$$J = \sigma E + J_s \qquad (58)$$

where E denotes the secondary electrical field (in V m⁻¹), and Js is the secondary (source) current density that is responsible for the observed electrical field (in A m⁻²). Since the time variation of the current density is very slow, we can use the quasistatic limit of the Maxwell equations, for which the conservation of the charge takes the form:

$$\nabla \cdot J = 0 \qquad (59)$$
By combining Equations (59) and (60), the scalar electrical potential (ψ) is the solution of a Poisson equation:

$$ \nabla \cdot (\sigma \nabla \psi) = \nabla \cdot \mathbf{J}_s \qquad (60) $$

Therefore, Equation (61) is both the classical self-potential equation and the equation that describes the secondary electrical potential fluctuations in time-induced polarization after the shutting down of the primary current or the electrical field. The electrical potential solution of Equation (61) can also be written in an integral form, which, once discretized, can, in turn, be written as:

$$ \mathbf{d} = \mathbf{K} \mathbf{m} \qquad (61) $$

The (N×M)-matrix (K) denotes the kernel, which is the collection of the Green's function for the electrical potential problems, and which is obtained by putting the element-source current densities of the unit strength into the three main directions for each cell that is used to discretize the subsurface; d is an N-vector in the collection of recorded electrical potentials at the N-electrodes; and m is an M-vector that corresponds to the source current-density times the volume of each cell (therefore, each element of m is expressed in A m) for the M-cells that are used to discretize the volume of interest. Alternatively, the integration volume could be accounted for in the elements of the kernel. Since the elements of the kernel should account for the boundary conditions and the complicated electrical-conductivity distribution, we cannot use analytical solutions for the Green's functions. We compute them numerically with the finite element software Comsol Multiphysics 4.0 by placing elementary dipoles in the three principal directions: x, y, and z.

The secondary voltage decays are treated as an N-observed data vector (dobs) at each time following the shutdown of the primary current. The objective function to minimize is defined as:

$$ P^{\lambda}(\mathbf{m}) = \|\mathbf{W}_d(\mathbf{K}\mathbf{m} - \mathbf{d}_{obs})\|^2 + \lambda \|\mathbf{W}_m\mathbf{m}\|^2 \qquad (62) $$

In this equation, Wm denotes a matrix that is used to regularize the inverse problem. Equation (63) can be minimized, and the optimized model vector can be obtained in a single iteration as:

$$ \mathbf{m} = [\mathbf{K}^T(\mathbf{W}_d^T\mathbf{W}_d)\mathbf{K} + \lambda(\mathbf{W}_m^T\mathbf{W}_m)]^{-1} \cdot [\mathbf{K}^T(\mathbf{W}_d^T\mathbf{W}_d)\mathbf{d}_{obs}] \qquad (63) $$

Various approaches are used to focus the source current densities and to properly normalize the model vector to avoid favoring sources that are close to the ground surface (see details in Soueid Ahmed and Revil, 2018 [99]). Figure 20 shows a tomography of the Cole-Cole parameters for an anomalous polarizable body in a tank [99]. This approach is very robust and fast, and especially because we solve a linear problem.

Figure 21 shows, for instance, that the normalized chargeability and the quadrature conductivity can be related to each other. This linear relationship allows for the conversion of the time-domain acquired data to complex conductivity data. If time-domain measurements are performed at different periods, we can even reconstruct complex conductivity spectra in a very simple way.
By combining Equations (59) and (60), the scalar electrical potential (ψ) is the solution of a Poisson equation:

$$ \nabla \cdot (\sigma \nabla \psi) = \nabla \cdot \mathbf{J}_s \qquad (60) $$

Therefore, Equation (61) is both the classical self-potential equation and the equation that describes the secondary electrical potential fluctuations in time-induced polarization after the shutting down of the primary current or the electrical field. The electrical potential solution of Equation (61) can also be written in an integral form, which, once discretized, can, in turn, be written as:

$$ \mathbf{d} = \mathbf{K} \mathbf{m} \qquad (61) $$

The (N×M)-matrix (K) denotes the kernel, which is the collection of the Green's function for the electrical potential problems, and which is obtained by putting the element-source current densities of the unit strength into the three main directions for each cell that is used to discretize the subsurface; d is an N-vector in the collection of recorded electrical potentials at the N-electrodes; and m is an M-vector that corresponds to the source current-density times the volume of each cell (therefore, each element of m is expressed in A m) for the M-cells that are used to discretize the volume of interest. Alternatively, the integration volume could be accounted for in the elements of the kernel. Since the elements of the kernel should account for the boundary conditions and the complicated electrical-conductivity distribution, we cannot use analytical solutions for the Green's functions. We compute them numerically with the finite element software Comsol Multiphysics 4.0 by placing elementary dipoles in the three principal directions: x, y, and z.

The secondary voltage decays are treated as an N-observed data vector (dobs) at each time following the shutdown of the primary current. The objective function to minimize is defined as:

$$ P^{\lambda}(\mathbf{m}) = \|\mathbf{W}_d(\mathbf{K}\mathbf{m} - \mathbf{d}_{obs})\|^2 + \lambda \|\mathbf{W}_m\mathbf{m}\|^2 \qquad (62) $$

In this equation, Wm denotes a matrix that is used to regularize the inverse problem. Equation (63) can be minimized, and the optimized model vector can be obtained in a single iteration as:

$$ \mathbf{m} = [\mathbf{K}^T(\mathbf{W}_d^T\mathbf{W}_d)\mathbf{K} + \lambda(\mathbf{W}_m^T\mathbf{W}_m)]^{-1} \cdot [\mathbf{K}^T(\mathbf{W}_d^T\mathbf{W}_d)\mathbf{d}_{obs}] \qquad (63) $$

Various approaches are used to focus the source current densities and to properly normalize the model vector to avoid favoring sources that are close to the ground surface (see details in Soueid Ahmed and Revil, 2018 [99]). Figure 20 shows a tomography of the Cole-Cole parameters for an anomalous polarizable body in a tank [99]. This approach is very robust and fast, and especially because we solve a linear problem.

Figure 21 shows, for instance, that the normalized chargeability and the quadrature conductivity can be related to each other. This linear relationship allows for the conversion of the time-domain acquired data to complex conductivity data. If time-domain measurements are performed at different periods, we can even reconstruct complex conductivity spectra in a very simple way.
a. True model

<figure><img src="image_24.png" /></figure>

Results of the inversion

<figure><img src="image_25.png" /></figure>

<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /></figure>

<figure><img src="image_28.png" /></figure>

**Figure 20.** True and estimated Cole–Cole parameter distributions using a source-current-density approach. We use 24 electrodes and only 2 current injections to reconstruct the distribution of the Cole–Cole parameters. (a) True distribution of the Cole–Cole parameter models. (b) Estimated electrical conductivity (σ) at Iteration 4 (RMS = 0.8%, 5% Gaussian noise added to the synthetic data). (c) Estimated intrinsic chargeability at time 0 (M). (d) Estimated relaxation time (τ). (e) Estimated frequency of Cole–Cole exponent (c). Modified from [99].
a. True model

<figure><img src="image_24.png" /></figure>

Results of the inversion

<figure><img src="image_25.png" /></figure>

<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /></figure>

<figure><img src="image_28.png" /></figure>

**Figure 20.** True and estimated Cole–Cole parameter distributions using a source-current-density approach. We use 24 electrodes and only 2 current injections to reconstruct the distribution of the Cole–Cole parameters. (a) True distribution of the Cole–Cole parameter models. (b) Estimated electrical conductivity (σ) at Iteration 4 (RMS = 0.8%, 5% Gaussian noise added to the synthetic data). (c) Estimated intrinsic chargeability at time 0 (M). (d) Estimated relaxation time (τ). (e) Estimated frequency of Cole–Cole exponent (c). Modified from [99].
<figure><img src="image_29.png" /></figure>

**Figure 21.** Normalized chargeability (1 Hertz-1 kHz) versus the quadrature conductivity for the volcanic samples from Hawaii (low salinity, NaCl, 25 °C). The “partial” normalized chargeability is determined as the difference in the in-phase conductivity between two frequencies (1 Hz and 1 kHz), while the quadrature conductivity is determined at the geometric mean frequency of 32 Hz. The plain line corresponds to a linear trend. The theory predicts that this trend is independent of the presence of magnetite (or pyrite) in the rocks. There are 93 samples in total for the shield volcanoes in Hawaii, and 32 additional samples from La Soufrière and Papandayan volcanoes. The value of the coefficient (α) can be predicted from the theory (unpublished results).

### 6.3. Tomography of the Relaxation Time

Our goal in this section is to show that we can easily perform a tomography of the main relaxation time by performing two tomographies of two distinct chargeabilities. First, it is worth mentioning that equivalent electrical circuits have been also broadly used to model the induced-polarization voltage-decay patterns in time-domain induced polarization. For the sake of simplicity, we consider the simplest circuit, which is, namely, an RC circuit (R and C stand for resistance and capacitance, respectively, see Figure 22).
<figure><img src="image_29.png" /></figure>

**Figure 21.** Normalized chargeability (1 Hertz-1 kHz) versus the quadrature conductivity for the volcanic samples from Hawaii (low salinity, NaCl, 25 °C). The “partial” normalized chargeability is determined as the difference in the in-phase conductivity between two frequencies (1 Hz and 1 kHz), while the quadrature conductivity is determined at the geometric mean frequency of 32 Hz. The plain line corresponds to a linear trend. The theory predicts that this trend is independent of the presence of magnetite (or pyrite) in the rocks. There are 93 samples in total for the shield volcanoes in Hawaii, and 32 additional samples from La Soufrière and Papandayan volcanoes. The value of the coefficient (α) can be predicted from the theory (unpublished results).

### 6.3. Tomography of the Relaxation Time

Our goal in this section is to show that we can easily perform a tomography of the main relaxation time by performing two tomographies of two distinct chargeabilities. First, it is worth mentioning that equivalent electrical circuits have been also broadly used to model the induced-polarization voltage-decay patterns in time-domain induced polarization. For the sake of simplicity, we consider the simplest circuit, which is, namely, an RC circuit (R and C stand for resistance and capacitance, respectively, see Figure 22).
<figure><img src="image_30.png" /><figcaption>Figure 22. Representation of the equivalent circuit of a mineralized rock in classical textbooks and papers. (a) The rock is composed of two types of pores or channels: one with ions free to move (Channel 1), and a second one with a pyrite grain and a blocking effect of the charge carriers at the interface between the metallic particle and the free pore water (which is an electrolyte with cations and anions). (b) Equivalent electric circuit that facilitates the response calculation in the time and frequency domains. For historical reasons, the polarization that is associated with the blocking of the charge carriers at the surface of the metallic grains is named the “electrode polarization”, which is characterized by a capacitance (C) with an impedance that is written as “$Z_C$”.</figcaption></figure>

Charging this circuit first with the normalized unit voltage, and then switching off the circuit leads to a decay in the secondary voltage, as:

$$ \psi(t) = \frac{R_1}{R_1 + R_2} \exp\left(-\frac{t}{(R_1 + R_2)C}\right) = M e^{-\frac{t}{\tau}} \quad (64) $$

The quantity $\psi(t)$ describes the difference in the potential between the voltage electrodes ($M$ and $N$); $C$ is the capacitance; $R_1$ and $R_2$ are the two resistances shown in Figure 22; and the chargeability is given by $M = R_1/(R_1 + R_2)$ (dimensionless). The relaxation time of the circuit, which is shown in Figure 22, is given by $\tau = (R_1 + R_2)C$. We have the following properties: $0 \le M \le 1$ and $0 \le \tau < +\infty$. The spectral response relative to this model can be obtained by injecting a harmonic current into the circuit, as shown in Figure 22. The associated impedance can be easily determined, and the associated complex resistivity corresponds to the so-called “Debye response”, which is written as:

$$ \rho^* = \frac{1}{\sigma^*} = \rho_0 \left[ 1 - M \left( 1 - \frac{1}{1 + i\omega\tau} \right) \right] \quad (65) $$

where $M = (\rho_0 - \rho_\infty)/\rho_0$ is the chargeability, and $\rho_0$ and $\rho_\infty$ denote the apparent resistivities at the null frequency (continuous or direct current (DC condition)) and at an infinite frequency, respectively, and which is also called the “instantaneous resistivity” ($\sigma_0 = 1/\rho_0$ and $\sigma_\infty = 1/\rho_\infty$). The null frequency means that the current would be applied for an infinite time and that all of the polarization phenomena in the material have reached steady-state conditions in a thermodynamic sense. To enforce that the chargeability ($M$) in Equation (66) is the same as in Equation (65), we first write Equation (66) as:
<figure><img src="image_30.png" /><figcaption>Figure 22. Representation of the equivalent circuit of a mineralized rock in classical textbooks and papers. (a) The rock is composed of two types of pores or channels: one with ions free to move (Channel 1), and a second one with a pyrite grain and a blocking effect of the charge carriers at the interface between the metallic particle and the free pore water (which is an electrolyte with cations and anions). (b) Equivalent electric circuit that facilitates the response calculation in the time and frequency domains. For historical reasons, the polarization that is associated with the blocking of the charge carriers at the surface of the metallic grains is named the “electrode polarization”, which is characterized by a capacitance (C) with an impedance that is written as “$Z_C$”.</figcaption></figure>

Charging this circuit first with the normalized unit voltage, and then switching off the circuit leads to a decay in the secondary voltage, as:

$$ \psi(t) = \frac{R_1}{R_1 + R_2} \exp\left(-\frac{t}{(R_1 + R_2)C}\right) = M e^{-\frac{t}{\tau}} \quad (64) $$

The quantity $\psi(t)$ describes the difference in the potential between the voltage electrodes ($M$ and $N$); $C$ is the capacitance; $R_1$ and $R_2$ are the two resistances shown in Figure 22; and the chargeability is given by $M = R_1/(R_1 + R_2)$ (dimensionless). The relaxation time of the circuit, which is shown in Figure 22, is given by $\tau = (R_1 + R_2)C$. We have the following properties: $0 \le M \le 1$ and $0 \le \tau < +\infty$. The spectral response relative to this model can be obtained by injecting a harmonic current into the circuit, as shown in Figure 22. The associated impedance can be easily determined, and the associated complex resistivity corresponds to the so-called “Debye response”, which is written as:

$$ \rho^* = \frac{1}{\sigma^*} = \rho_0 \left[ 1 - M \left( 1 - \frac{1}{1 + i\omega\tau} \right) \right] \quad (65) $$

where $M = (\rho_0 - \rho_\infty)/\rho_0$ is the chargeability, and $\rho_0$ and $\rho_\infty$ denote the apparent resistivities at the null frequency (continuous or direct current (DC condition)) and at an infinite frequency, respectively, and which is also called the “instantaneous resistivity” ($\sigma_0 = 1/\rho_0$ and $\sigma_\infty = 1/\rho_\infty$). The null frequency means that the current would be applied for an infinite time and that all of the polarization phenomena in the material have reached steady-state conditions in a thermodynamic sense. To enforce that the chargeability ($M$) in Equation (66) is the same as in Equation (65), we first write Equation (66) as:
$$ \frac{\rho^*(\omega)}{\rho_0} = (1 - M) + M \left( \frac{1}{1 + i\omega\tau} \right) = \frac{\tilde{\psi}^*(\omega)}{\psi_0} \qquad (66) $$

where $\tilde{\psi}^*$ denotes the complex potential that is proportional to the complex resistivity, through Ohm's law, in the Fourier domain. Taking the inverse Fourier transform provides the impulse response of this filter, precisely as:

$$ h(t) = \text{FT}^{-1} \left[ \frac{\tilde{\psi}^*(\omega)}{\psi_0} \right] \qquad (67) $$

$$ h(t) = \text{FT}^{-1} \left[ (1 - M) + M \left( \frac{1}{1 + i\omega\tau} \right) \right] \qquad (68) $$

$$ h(t) = (1 - M)\delta(t) + M\tau e^{-\frac{t}{\tau}} \qquad (69) $$

where $\text{FT}^{-1}$ denotes the inverse Fourier transform. The response to a unit step is obtained by the integration of this impulse response, and it would correspond to the time that the current is just sent, which passes from 0 to +1, and from $t = 0$. We obtain, in the time domain, the following voltage-decay response:

$$ \frac{\psi(t)}{\psi_0} = \int_0^t h(t)dt = 1 - Me^{-\frac{t}{\tau}} \qquad (70) $$

When discharging (that is, after switching off the injected current), we just obtain the decaying-potential ratio: $\psi(t)/\psi_0 = M \exp(-t/\tau)$.

Another chargeability was defined in the realm of the mining industry and it is called the "integral chargeability" (expressed in s), which is defined as:

$$ \hat{M} = \frac{1}{\psi_0} \int_0^\infty \psi(t)dt \qquad (71) $$

$$ \hat{M} = \int_{0}^{\infty} M e^{-\frac{t}{\tau}} dt = M\tau \qquad (72) $$

In other words, in the context of the electrical circuit that is shown in Figure 22, the relaxation time ($\tau$) is the ratio of the integral chargeability by the chargeability. This relation is also valid for a Cole-Cole model instead of a Debye model (Florsch, personal communication 2019). In terms of the tomography, we can therefore perform a tomography of the integral chargeability according to one of the approaches underlined above, and a second chargeability tomography by using the classical apparent chargeability. Then, cell by cell, the relaxation time is obtained by dividing the inverted integral chargeability by the chargeability. We will show one example of the relaxation-time tomography that was obtained by using this approach in the next section in the context of a metallic tank (Section 7.4).

# 7. Applications

## 7.1. Treasure Quest: Finding Pyrite in a Sandbox

The first test of the method is to use electrical-resistivity tomography and induced polarization to find an area that is rich in pyrite cubes in a sandbox experiment. The sandbox is saturated with tap water. The results that are shown in Figure 23 indicate that the electrical-conductivity distribution is not sensitive enough to show the position of the disseminated pyrite cubes, while the induced polarization clearly shows where the ore is located (Figure 24). This is consistent with our model, for which the presence of ore has only a second-order effect on the conductivity distribution, while the volume content of the ores directly controls the chargeability of the material. Furthermore, at low frequencies, the pyrite cubes would appear to be insulators that increase the resistivity of the material (Figure 17). The observation of resistive bodies inside a more conductive background is
$$ \frac{\rho^*(\omega)}{\rho_0} = (1 - M) + M \left( \frac{1}{1 + i\omega\tau} \right) = \frac{\tilde{\psi}^*(\omega)}{\psi_0} \qquad (66) $$

where $\tilde{\psi}^*$ denotes the complex potential that is proportional to the complex resistivity, through Ohm's law, in the Fourier domain. Taking the inverse Fourier transform provides the impulse response of this filter, precisely as:

$$ h(t) = \mathrm{FT}^{-1} \left[ \frac{\tilde{\psi}^*(\omega)}{\psi_0} \right] \qquad (67) $$

$$ h(t) = \mathrm{FT}^{-1} \left[ (1 - M) + M \left( \frac{1}{1 + i\omega\tau} \right) \right] \qquad (68) $$

$$ h(t) = (1 - M)\delta(t) + M\tau e^{-\frac{t}{\tau}} \qquad (69) $$

where $\mathrm{FT}^{-1}$ denotes the inverse Fourier transform. The response to a unit step is obtained by the integration of this impulse response, and it would correspond to the time that the current is just sent, which passes from 0 to +1, and from $t = 0$. We obtain, in the time domain, the following voltage-decay response:

$$ \frac{\psi(t)}{\psi_0} = \int_0^t h(t)dt = 1 - Me^{-\frac{t}{\tau}} \qquad (70) $$

When discharging (that is, after switching off the injected current), we just obtain the decaying-potential ratio: $\psi(t)/\psi_0 = M \exp(-t/\tau)$.

Another chargeability was defined in the realm of the mining industry and it is called the "integral chargeability" (expressed in s), which is defined as:

$$ \hat{M} = \frac{1}{\psi_0} \int_0^\infty \psi(t)dt \qquad (71) $$

$$ \hat{M} = \int_{0}^{\infty} M e^{-\frac{t}{\tau}} dt = M\tau \qquad (72) $$

In other words, in the context of the electrical circuit that is shown in Figure 22, the relaxation time ($\tau$) is the ratio of the integral chargeability by the chargeability. This relation is also valid for a Cole-Cole model instead of a Debye model (Florsch, personal communication 2019). In terms of the tomography, we can therefore perform a tomography of the integral chargeability according to one of the approaches underlined above, and a second chargeability tomography by using the classical apparent chargeability. Then, cell by cell, the relaxation time is obtained by dividing the inverted integral chargeability by the chargeability. We will show one example of the relaxation-time tomography that was obtained by using this approach in the next section in the context of a metallic tank (Section 7.4).

# 7. Applications

## 7.1. Treasure Quest: Finding Pyrite in a Sandbox

The first test of the method is to use electrical-resistivity tomography and induced polarization to find an area that is rich in pyrite cubes in a sandbox experiment. The sandbox is saturated with tap water. The results that are shown in Figure 23 indicate that the electrical-conductivity distribution is not sensitive enough to show the position of the disseminated pyrite cubes, while the induced polarization clearly shows where the ore is located (Figure 24). This is consistent with our model, for which the presence of ore has only a second-order effect on the conductivity distribution, while the volume content of the ores directly controls the chargeability of the material. Furthermore, at low frequencies, the pyrite cubes would appear to be insulators that increase the resistivity of the material (Figure 17). The observation of resistive bodies inside a more conductive background is
inherently difficult. Figure 25 shows the results of another set of sandbox experiments that used a target that was composed of sand mixed with pyrite and bentonite. The targets were immersed in a sand background. The sandboxes were saturated with tap water. These experiments confirm the previous findings and show that, for bentonite, it is better to use the normalized chargeability to show up the target because of its high cation exchange capacity (CEC).

<figure><img src="image_31.png" /></figure>

**Figure 23.** Treasure hunt in a sandbox. The goal of this test is to see which of resistivity or induced polarization tomography performs better to localize the ore in the sandbox without prior knowledge of its position (a) Sandbox with the position of the electrodes and the cylindrical area in which the pyrite is located. For each acquisition, we inject and retrieve the current using two electrodes A and B and all the other electrodes are voltage electrodes $M_i$ with respect to a reference electrode N = 1 (Ref, see Figure 23). (b) Electrical resistivity tomography. The area encapsulated by the white line corresponds to the true position of the target (pyrite grains mixed with sand). The resistivity tomography is unable to show up the position of the ores because the effect of the ores has only a minor effect on the DC conductivity of the mixture. The induced polarization tomography is shown in Figure 24.
inherently difficult. Figure 25 shows the results of another set of sandbox experiments that used a target that was composed of sand mixed with pyrite and bentonite. The targets were immersed in a sand background. The sandboxes were saturated with tap water. These experiments confirm the previous findings and show that, for bentonite, it is better to use the normalized chargeability to show up the target because of its high cation exchange capacity (CEC).

<figure><img src="image_31.png" /><figcaption>**Figure 23.** Treasure hunt in a sandbox. The goal of this test is to see which of resistivity or induced polarization tomography performs better to localize the ore in the sandbox without prior knowledge of its position (a) Sandbox with the position of the electrodes and the cylindrical area in which the pyrite is located. For each acquisition, we inject and retrieve the current using two electrodes A and B and all the other electrodes are voltage electrodes $M_i$ with respect to a reference electrode N = 1 (Ref, see Figure 23). (b) Electrical resistivity tomography. The area encapsulated by the white line corresponds to the true position of the target (pyrite grains mixed with sand). The resistivity tomography is unable to show up the position of the ores because the effect of the ores has only a minor effect on the DC conductivity of the mixture. The induced polarization tomography is shown in Figure 24.</figcaption></figure>
<figure><img src="image_32.png" /><figcaption>**Figure 24.** Use of induced polarization to localize an ore body in a sandbox. (a) Chargeability tomography from the classical Gauss–Newton inversion of the induced polarization data and using the setup shown in Figure 23. Inversion of the chargeability with 12 bipoles [A, B] (A and B are current electrodes). All the other electrodes are voltage electrodes $M_i$ with respect to a reference electrode N. These electrodes are used to measure the secondary voltages after shutting down the primary electrical current. (b) Inversion of the source-current density ($J_s$) with four bipoles [A, B]. The area encapsulated by the white line corresponds to the position of the target (pyrite grains mixed with sand). The black dash line shows the position of the array of electrodes as shown in Figure 24 (only four bipoles used for the source-current-density inversion). The white line underlines the area where the ore is truly located. Adapted from [100].</figcaption></figure>
<figure><img src="image_32.png" /><figcaption>**Figure 24.** Use of induced polarization to localize an ore body in a sandbox. (a) Chargeability tomography from the classical Gauss–Newton inversion of the induced polarization data and using the setup shown in Figure 23. Inversion of the chargeability with 12 bipoles [A, B] (A and B are current electrodes). All the other electrodes are voltage electrodes $M_i$ with respect to a reference electrode N. These electrodes are used to measure the secondary voltages after shutting down the primary electrical current. (b) Inversion of the source-current density ($J_s$) with four bipoles [A, B]. The area encapsulated by the white line corresponds to the position of the target (pyrite grains mixed with sand). The black dash line shows the position of the array of electrodes as shown in Figure 24 (only four bipoles used for the source-current-density inversion). The white line underlines the area where the ore is truly located. Adapted from [100].</figcaption></figure>
<figure><img src="image_33.png" /></figure>

**Figure 25.** Comparison of the resistivity, chargeability, and normalized chargeability tomograms between localized pyrite and localized bentonite in a sandy background. (a) Resistivity tomography in the case of disseminated pyrite in a sandy background. (b) Chargeability tomography in the case of disseminated pyrite in a sandy background. (c) Normalized chargeability tomography in the case of disseminated pyrite. (d) Resistivity tomography in the case of bentonite. (e) Chargeability tomography in the case of disseminated pyrite. (f) Normalized chargeability tomography in the case of bentonite. Adapted from [100].

## 7.2. Treasure Quest: Finding Slag at an Archeological Site

Dumped slag heaps that have been used over the history of humankind can be seen as valuable deposits to feed the growing need for metal resources. As shown in the previous section, the spatial delineation of slag heaps cannot easily be obtained from the resistivity distribution. The induced-polarization method can be used to fulfill this goal. A least-square 3D induced-polarization inversion based on the least-square technique is used to obtain the chargeability distribution that delimits a slag heap at the archeological site of Saint-Vincent sur L'Isle, Dordogne, France (Figure 26). From the linear relationship that is determined between the chargeability and the volumetric metal content, the metallic volume of the slag heaps can be directly determined. A volume of high chargeability (above the background value of 0.08) indicates the position of the slag heap. The linear relationship between the chargeability and the metallic volumetric content (see Figure 14) is used to provide an image of the distribution of the metals that are present at this archeological site (Figure 27).
<figure><img src="image_33.png" /></figure>

**Figure 25.** Comparison of the resistivity, chargeability, and normalized chargeability tomograms between localized pyrite and localized bentonite in a sandy background. (a) Resistivity tomography in the case of disseminated pyrite in a sandy background. (b) Chargeability tomography in the case of disseminated pyrite in a sandy background. (c) Normalized chargeability tomography in the case of disseminated pyrite. (d) Resistivity tomography in the case of bentonite. (e) Chargeability tomography in the case of disseminated pyrite. (f) Normalized chargeability tomography in the case of bentonite. Adapted from [100].

## 7.2. Treasure Quest: Finding Slag at an Archeological Site

Dumped slag heaps that have been used over the history of humankind can be seen as valuable deposits to feed the growing need for metal resources. As shown in the previous section, the spatial delineation of slag heaps cannot easily be obtained from the resistivity distribution. The induced-polarization method can be used to fulfill this goal. A least-square 3D induced-polarization inversion based on the least-square technique is used to obtain the chargeability distribution that delimits a slag heap at the archeological site of Saint-Vincent sur L'Isle, Dordogne, France (Figure 26). From the linear relationship that is determined between the chargeability and the volumetric metal content, the metallic volume of the slag heaps can be directly determined. A volume of high chargeability (above the background value of 0.08) indicates the position of the slag heap. The linear relationship between the chargeability and the metallic volumetric content (see Figure 14) is used to provide an image of the distribution of the metals that are present at this archeological site (Figure 27).
<figure><img src="image_34.png" /><figcaption>**Figure 26.** Induced-polarization survey performed at an archeological test site in France. A total of 21 profiles (24 electrodes each) are used to perform 3D tomograms of chargeability and electrical resistivity. The color scale represents the topography of the area (altitude in meters (m) above sea level). The present survey comprises 3150 dipole–dipole measurements. Modified from [98].</figcaption></figure>

<figure><img src="image_35.png" /><figcaption>**Figure 27.** Tomography of the volume fraction of metal (in %) in a slag heap that corresponds to an archeological site. The tomogram is directly determined by using the chargeability tomogram. We use a threshold chargeability for the background of the chargeability: $M_b = 0.08 \approx R$ (from experimental data performed in the laboratory by using materials from the archeological test site). Modified from [98].</figcaption></figure>
<figure><img src="image_34.png" /><figcaption>**Figure 26.** Induced-polarization survey performed at an archeological test site in France. A total of 21 profiles (24 electrodes each) are used to perform 3D tomograms of chargeability and electrical resistivity. The color scale represents the topography of the area (altitude in meters (m) above sea level). The present survey comprises 3150 dipole–dipole measurements. Modified from [98].</figcaption></figure>

<figure><img src="image_35.png" /><figcaption>**Figure 27.** Tomography of the volume fraction of metal (in %) in a slag heap that corresponds to an archeological site. The tomogram is directly determined by using the chargeability tomogram. We use a threshold chargeability for the background of the chargeability: $M_b = 0.08 \approx R$ (from experimental data performed in the laboratory by using materials from the archeological test site). Modified from [98].</figcaption></figure>
### 7.3. Graphitic Schists and Permafrost

This section demonstrates the effects of the temperature below the freezing point on rocks that contain metals, with an application to graphitic schists. Spectral induced polarization spectra performed with graphite graphitic show a very strong effect of temperature below the freezing point (see Figure 28). Since the temperature has a strong effect on the instantaneous conductivity, it is not surprising, in the framework of the model discussed above, that the relaxation time is also dependent on the instantaneous conductivity and on the temperature (Figure 29). As long as the temperature decreases, the electrical conductivity decreases, and the relaxation time increases. Figures 30 and 31 demonstrate the effects of both graphite and temperature on the electrical conductivity and chargeability tomography.

<figure><img src="image_36.png" /><figcaption>**Figure 28.** Results of spectral induced polarization measurements with graphite immersed in water-saturated sand. Note the effect of temperature below the freezing point. Below the freezing point, the relaxation times (from the peak of the quadrature conductivity) seem to be very dependent on the temperature. Modified from [95].</figcaption></figure>
### 7.3. Graphitic Schists and Permafrost

This section demonstrates the effects of the temperature below the freezing point on rocks that contain metals, with an application to graphitic schists. Spectral induced polarization spectra performed with graphite graphitic show a very strong effect of temperature below the freezing point (see Figure 28). Since the temperature has a strong effect on the instantaneous conductivity, it is not surprising, in the framework of the model discussed above, that the relaxation time is also dependent on the instantaneous conductivity and on the temperature (Figure 29). As long as the temperature decreases, the electrical conductivity decreases, and the relaxation time increases. Figures 30 and 31 demonstrate the effects of both graphite and temperature on the electrical conductivity and chargeability tomography.

<figure><img src="image_36.png" /><figcaption>**Figure 28.** Results of spectral induced polarization measurements with graphite immersed in water-saturated sand. Note the effect of temperature below the freezing point. Below the freezing point, the relaxation times (from the peak of the quadrature conductivity) seem to be very dependent on the temperature. Modified from [95].</figcaption></figure>
<figure><img src="image_37.png" /><figcaption>**Figure 29.** Relationship between the instantaneous conductivity and the low-frequency relaxation time for samples with metallic grains (galena, chalcopyrite, graphite, magnetite, pyrite) immersed in sand as background. The low-frequency relaxation time is inversely proportional to the instantaneous conductivity of the mixture. In the equation, *A* is a fitting constant for each material. Modified from [95].</figcaption></figure>

<figure><img src="image_38.png" /><figcaption>**Figure 30.** Large-scale electrical conductivity and normalized chargeability tomograms crossing an area that is associated with a rock glacier with permafrost below it. (a) Electrical conductivity tomogram showing position and depth of the permafrost (about 70 m). (b) Normalized chargeability tomogram (high normalized chargeability zones correlated with graphitic schist formation). The frozen area is in blue. Modified from [43].</figcaption></figure>
<figure><img src="image_37.png" /><figcaption>**Figure 29.** Relationship between the instantaneous conductivity and the low-frequency relaxation time for samples with metallic grains (galena, chalcopyrite, graphite, magnetite, pyrite) immersed in sand as background. The low-frequency relaxation time is inversely proportional to the instantaneous conductivity of the mixture. In the equation, *A* is a fitting constant for each material. Modified from [95].</figcaption></figure>

<figure><img src="image_38.png" /><figcaption>**Figure 30.** Large-scale electrical conductivity and normalized chargeability tomograms crossing an area that is associated with a rock glacier with permafrost below it. (a) Electrical conductivity tomogram showing position and depth of the permafrost (about 70 m). (b) Normalized chargeability tomogram (high normalized chargeability zones correlated with graphitic schist formation). The frozen area is in blue. Modified from [43].</figcaption></figure>
<figure><img src="image_39.png" /></figure>

**Figure 31.** Normalized chargeability versus electrical conductivity at different temperatures for inverted field data (red- and black-filled circles) and laboratory experiments (open squares and circles). The data above the lower plain line indicate that the sample contains semiconductor or metallic particles. The normalized chargeability is calculated from the equation: $M_{II} = M \sigma_{\infty}$, where $\sigma_{\infty}$ and M (= M1) come from the fitting with a double Cole-Cole model, which corresponds to Equation (42). The grey area corresponds to the graphite-rich zone. The high and low values of the normalized chargeability and conductivity are also controlled by the lithology. Modified from [43].

### 7.4. Localization of a Metallic Tank

We illustrate now the use of a relaxation-time tomography to perform a tomography of a metallic tank underground (see Section 6.3 above). This case study is illustrated in Figure 32. We see that both the chargeability and the relaxation-time tomography can be used to detect a metallic tank. The tomography of the relaxation time provides a very clear image of the tank.
<figure><img src="image_39.png" /></figure>

**Figure 31.** Normalized chargeability versus electrical conductivity at different temperatures for inverted field data (red- and black-filled circles) and laboratory experiments (open squares and circles). The data above the lower plain line indicate that the sample contains semiconductor or metallic particles. The normalized chargeability is calculated from the equation: $M_{II} = M \sigma_{\infty}$, where $\sigma_{\infty}$ and M (= M1) come from the fitting with a double Cole-Cole model, which corresponds to Equation (42). The grey area corresponds to the graphite-rich zone. The high and low values of the normalized chargeability and conductivity are also controlled by the lithology. Modified from [43].

### 7.4. Localization of a Metallic Tank

We illustrate now the use of a relaxation-time tomography to perform a tomography of a metallic tank underground (see Section 6.3 above). This case study is illustrated in Figure 32. We see that both the chargeability and the relaxation-time tomography can be used to detect a metallic tank. The tomography of the relaxation time provides a very clear image of the tank.
<figure><img src="image_40.png" /><figcaption>**Figure 32.** Tomography of a dam with the position of a metallic tank underground. (a) Tomography of the electrical conductivity. (b) Tomography of the chargeability. (c) Tomography of the relaxation time (in s). Note how the tomography of the relaxation time can be used to assess the position of the metallic tank with a precision that is better than using chargeability tomography. Unpublished results.</figcaption></figure>

### 7.5. Application to a Deep Cu-Pb-Zn Deposit

Although the traditional induced polarization method in the time and frequency domains [28,52–55] has achieved great successes, it cannot necessarily meet the needs of large-scale induced-polarization exploration in the presence of strong electromagnetic interferences. The idea that underlies spread-spectrum communication [101] was introduced into induced-polarization instrumentation by Xi et al. (2013, 2014) [102,103] and their corresponding data processing. Liu et al. (2016, 2017) [104,105] proposed a spread-spectrum induced-polarization (SSIP) method that is based on ZigBee wireless sensor networks and GPS timing to realize synchronous real-time full-waveform acquisition and the real-time data-quality monitoring of massive SSIP receivers [103,106]. At the same time, the GPS timing and synchronization of the SSIP transmitter and the spread-spectrum IP receiver, as well as the ZigBee wireless networking of the spread-spectrum IP transmitter, the receiver array, and the personal computer (control center) are realized [103]. On the basis of these technologies, exploration depths greater than 1000 m, even in the presence of strong electromagnetic interferences, could be achieved [105,107,108].
<figure><img src="image_40.png" /><figcaption>**Figure 32.** Tomography of a dam with the position of a metallic tank underground. (a) Tomography of the electrical conductivity. (b) Tomography of the chargeability. (c) Tomography of the relaxation time (in s). Note how the tomography of the relaxation time can be used to assess the position of the metallic tank with a precision that is better than using chargeability tomography. Unpublished results.</figcaption></figure>

### 7.5. Application to a Deep Cu-Pb-Zn Deposit

Although the traditional induced polarization method in the time and frequency domains [28,52–55] has achieved great successes, it cannot necessarily meet the needs of large-scale induced-polarization exploration in the presence of strong electromagnetic interferences. The idea that underlies spread-spectrum communication [101] was introduced into induced-polarization instrumentation by Xi et al. (2013, 2014) [102,103] and their corresponding data processing. Liu et al. (2016, 2017) [104,105] proposed a spread-spectrum induced-polarization (SSIP) method that is based on ZigBee wireless sensor networks and GPS timing to realize synchronous real-time full-waveform acquisition and the real-time data-quality monitoring of massive SSIP receivers [103,106]. At the same time, the GPS timing and synchronization of the SSIP transmitter and the spread-spectrum IP receiver, as well as the ZigBee wireless networking of the spread-spectrum IP transmitter, the receiver array, and the personal computer (control center) are realized [103]. On the basis of these technologies, exploration depths greater than 1000 m, even in the presence of strong electromagnetic interferences, could be achieved [105,107,108].
In order to showcase a proof of concept of such an exploration tool, it was applied to a Cu-Pb-Zn polymetallic mine in Baiyin City, Gansu Province, China. The silver-copper-lead-zinc polymetallic deposit belongs to the VMS (volcanogenic massive sulfide) deposits, which was once the largest copper mine in China and which has a mining history of hundreds of years [109]. The exploration section is located west of the ore body that is being mined underground. Because of the existence of strong electromagnetic interferences, such as underground mineral mining, mineral-processing plants, and smelters around the profile, it is extremely challenging to achieve an IP exploration depth of 1000 m. In order to realize the high-resolution exploration of SSIP at depths higher than 1000 m, 25 sets of four-channel SSIP receivers were used for synchronous observation, which achieved a minimum measurement dipole spacing of 20 m, and which covered a profile length of 2000 m. In order to achieve high-resolution exploration, the observation starts from the minimum distance of 20 m from the supply electrode, and the center of the supply dipole is the center of the survey line. Then, increase the power-supply dipole distance by 80 m each time inside the survey line. To the outside of the survey line, increase the power-supply dipole distance by 400 m each time until the maximum power-supply dipole distance is 5200 m (Figure 33). Each time that power is supplied, the 25 receivers start observing simultaneously.

<figure><img src="image_41.png" /><figcaption>**Figure 33.** Nonconventional gradient array protocol for the survey line. Current was supplied by using 34 current electrode groups with a maximum current electrode spacing of 5200 m, and a minimum spacing of 20 m. The spacing between the potential electrodes is 20 m (modified from [105]).</figcaption></figure>

There are various electromagnetic interferences in the measured waveforms. Considering the short-term and random nature of electromagnetic interferences, robust superposition and correlation detection are used to enhance the signal-to-noise ratio [104,105]. Considering that the maximum power-supply dipole distance is 5200 m, there is also electromagnetic coupling induction in the measurement results. Therefore, the relative phase spectrum is used to suppress the influence of electromagnetic coupling interference [110].

Figure 34 shows the resistivity-and-chargeability (relative phase is converted into chargeability in the inversion) profile after the inversion using ZondRes2D, as well as the Cu-Pb-Zn-ore-body morphology that was verified by drilling. From the resistivity profile in the figure, the existence of a steeply standing fault can be roughly inferred, but the existence of an ore body cannot be inferred. Fortunately, an oblong polarizability anomaly can be found from the polarizability profile of the figure. The location of the chargeability anomaly is consistent with the fault location, which is implied in the resistivity profile, so that it can be reasonably inferred that the IP anomaly is caused by the VMS ore body. Subsequent drill holes have verified the above inferences, and the IP anomaly that was discovered this time is highly consistent with the verified VMS ore body.
In order to showcase a proof of concept of such an exploration tool, it was applied to a Cu-Pb-Zn polymetallic mine in Baiyin City, Gansu Province, China. The silver-copper-lead-zinc polymetallic deposit belongs to the VMS (volcanogenic massive sulfide) deposits, which was once the largest copper mine in China and which has a mining history of hundreds of years [109]. The exploration section is located west of the ore body that is being mined underground. Because of the existence of strong electromagnetic interferences, such as underground mineral mining, mineral-processing plants, and smelters around the profile, it is extremely challenging to achieve an IP exploration depth of 1000 m. In order to realize the high-resolution exploration of SSIP at depths higher than 1000 m, 25 sets of four-channel SSIP receivers were used for synchronous observation, which achieved a minimum measurement dipole spacing of 20 m, and which covered a profile length of 2000 m. In order to achieve high-resolution exploration, the observation starts from the minimum distance of 20 m from the supply electrode, and the center of the supply dipole is the center of the survey line. Then, increase the power-supply dipole distance by 80 m each time inside the survey line. To the outside of the survey line, increase the power-supply dipole distance by 400 m each time until the maximum power-supply dipole distance is 5200 m (Figure 33). Each time that power is supplied, the 25 receivers start observing simultaneously.

<figure><img src="image_41.png" /><figcaption>**Figure 33.** Nonconventional gradient array protocol for the survey line. Current was supplied by using 34 current electrode groups with a maximum current electrode spacing of 5200 m, and a minimum spacing of 20 m. The spacing between the potential electrodes is 20 m (modified from [105]).</figcaption></figure>

There are various electromagnetic interferences in the measured waveforms. Considering the short-term and random nature of electromagnetic interferences, robust superposition and correlation detection are used to enhance the signal-to-noise ratio [104,105]. Considering that the maximum power-supply dipole distance is 5200 m, there is also electromagnetic coupling induction in the measurement results. Therefore, the relative phase spectrum is used to suppress the influence of electromagnetic coupling interference [110].

Figure 34 shows the resistivity-and-chargeability (relative phase is converted into chargeability in the inversion) profile after the inversion using ZondRes2D, as well as the Cu-Pb-Zn-ore-body morphology that was verified by drilling. From the resistivity profile in the figure, the existence of a steeply standing fault can be roughly inferred, but the existence of an ore body cannot be inferred. Fortunately, an oblong polarizability anomaly can be found from the polarizability profile of the figure. The location of the chargeability anomaly is consistent with the fault location, which is implied in the resistivity profile, so that it can be reasonably inferred that the IP anomaly is caused by the VMS ore body. Subsequent drill holes have verified the above inferences, and the IP anomaly that was discovered this time is highly consistent with the verified VMS ore body.
<figure><img src="image_42.png" /><figcaption>**Figure 34.** VMS Cu-Pb-Zn ore bodies imaged by SSIP. (a) The 2D resistivity tomogram. (b) The 2D IP chargeability tomogram. The chargeability ($M = \eta$) is here expressed in % rather than in a fraction from 0 to 1. Both tomograms are derived from the apparent complex resistivity data at a frequency of 0.41 Hz. The thick black lines show the shape and location of the Cu-Pb-Zn ore bodies defined by the eight wells, which are projected onto the plane of the tomogram (modified from [105]).</figcaption></figure>

**8. Future Trends in Marine Exploration**

In this section, we want to discuss the possible integration of the electrical conductivity/inducedpolarization methods with the self-potential technique. The self-potential method is a passive geophysical technique that measures the electrical field that is associated with the natural source of the currents in the subsurface [111]. Induced-polarization surveys can be performed at the seafloor [112]. Negative self-potential (SP) anomalies near seafloor massive sulfide deposits are associated with redox-active environments, in which metallic bodies are located [113].

Recently, self-potential surveys were used to investigate SMS (seafloor massive sulfides) deposits at the seafloor. The inversion of the self-potential data can be used to obtain a 3D tomogram of the causative-source-current density, which, in turn, can underline the ore body [114]. Meanwhile, a self-potential survey uses very simple equipment (nonpolarizing electrodes and a sensitive voltmeter), which can be easily mounted together with active (galvanometric or induction-based) sensors to image the resistivity and the induced polarization [115].

In order to show the new capabilities of the self-potential method, we discuss here two case studies. Zhu et al. (2020) [115] conducted a combined self-potential and transient
<figure><img src="image_42.png" /><figcaption>**Figure 34.** VMS Cu-Pb-Zn ore bodies imaged by SSIP. (a) The 2D resistivity tomogram. (b) The 2D IP chargeability tomogram. The chargeability ($M = \eta$) is here expressed in % rather than in a fraction from 0 to 1. Both tomograms are derived from the apparent complex resistivity data at a frequency of 0.41 Hz. The thick black lines show the shape and location of the Cu-Pb-Zn ore bodies defined by the eight wells, which are projected onto the plane of the tomogram (modified from [105]).</figcaption></figure>

## 8. Future Trends in Marine Exploration

In this section, we want to discuss the possible integration of the electrical conductivity/inducedpolarization methods with the self-potential technique. The self-potential method is a passive geophysical technique that measures the electrical field that is associated with the natural source of the currents in the subsurface [111]. Induced-polarization surveys can be performed at the seafloor [112]. Negative self-potential (SP) anomalies near seafloor massive sulfide deposits are associated with redox-active environments, in which metallic bodies are located [113].

Recently, self-potential surveys were used to investigate SMS (seafloor massive sulfides) deposits at the seafloor. The inversion of the self-potential data can be used to obtain a 3D tomogram of the causative-source-current density, which, in turn, can underline the ore body [114]. Meanwhile, a self-potential survey uses very simple equipment (nonpolarizing electrodes and a sensitive voltmeter), which can be easily mounted together with active (galvanometric or induction-based) sensors to image the resistivity and the induced polarization [115].

In order to show the new capabilities of the self-potential method, we discuss here two case studies. Zhu et al. (2020) [115] conducted a combined self-potential and transient
electromagnetic survey on the southwest Indian ridge. Six horizontal AgCl electrodes were attached to the transient electromagnetic system to passively record the horizontal electrical-field component (Figure 35). The natural electrical field was extracted from the raw data and was integrated to obtain the self-potential signals. Negative self-potential anomalies (~-27 mV) were observed at the known Yuhuang hydrothermal field. The inversion of these self-potential data reveals localized SMS deposits with a thickness of ~65 m (Figure 36).

The second case study is reported by Kasaya et al. (2020) [116], who conducted a self-potential and marine direct-current-resistivity survey by using two autonomous underwater vehicles (Figures 37 and 38). The self-potential and the apparent resistivity were obtained along a survey line that crossed the known hydrothermal mounds of sulfide ore in the mid-Okinawa Trough. Similar to the previous case study, the raw self-potential data were integrated to obtain a potential field profile. The distribution of the negative self-potential anomalies obtained and the low apparent resistivity (0.2 ohm-m or less) show good agreement with the locations of the known hydrothermal deposits (Figure 39).

Recent studies have demonstrated that self-potential tomography can play an important role in the exploration of SMS deposits. However, the electrical-conductivity distribution of SMS deposits is an important ingredient in the inversion of the recorded self-potential data. Indeed, the self-potential governing equation is given as [111]:

$$ \nabla \cdot (\sigma \nabla \varphi) = \nabla \cdot \mathbf{J}_s \qquad (73) $$

where $\sigma$ denotes the electrical conductivity, $\varphi$ is the self-potential field, and $\mathbf{J}_s$ denotes the source-current density that is associated with the ore body. We can develop the left-hand side to obtain [117]:

$$ \nabla^2 \varphi = \frac{\nabla \cdot \mathbf{J}_s}{\sigma} + \nabla \ln \sigma \cdot \mathbf{E} \qquad (74) $$

where $\mathbf{E} = -\nabla \varphi$ is the electrical field. Equation (74) shows that, in order to obtain the source-current density, knowledge of the conductivity distribution ($\sigma$) below the seafloor is required [111,118]. Future analysis must combine induced-polarization tomography with the inversion of the self-potential data in order to better quantify SMS deposits.

<figure><img src="image_43.png" /><figcaption>Figure 35. Sketch of the deep-towed geophysical instrument that can also be used to measure the self-potential and possibly the resistivity and the induced polarization. The green line denotes the cable-link-control terminal and the time-domain electromagnetic (TEM) coil. The numbers on the cable show the location of 6 Ag/AgCl nonpolarizing electrodes. The data logger is mounted inside of the control terminal (modified from [115]).</figcaption></figure>