Rethinking earthquake-related DC-ULF electromagnetic phenomena: towards a physics-based approach

Numerous electromagnetic changes possibly related with earthquakes have been independently reported and have even been attempted to apply to short-term prediction of earthquakes. However, there are active debates on the above issue because the seismogenic process is rather complicated and the studies have been mainly empirical (i.e. a kind of experience-based approach). Thus, a physics-based study would be helpful for understanding earthquake-related electromagnetic phenomena and strengthening their applications. As a potential physics-based approach, I present an integrated research scheme, taking into account the interaction among observation, methodology, and physical model. For simplicity, this work focuses only on the earthquake-related DC-ULF electromagnetic phenomena. The main approach includes the following key problems: (1) how to perform a reliable and appropriate observation with some clear physical quantities; (2) how to develop a robust methodology to reveal weak earthquake-related electromagnetic signals from noisy background; and (3) how to develop plausible physical models based on theoretical analyses and/or laboratory experiments for the explanation of the earthquake-related electromagnetic signals observed in the field conditions.


Introduction
As one important branch of geo-electromagnetism, the study of earthquake-related electromagnetic (EM) phenomena plays a key role in the study of seismic hazard.The initial interest of seismo-electromagnetism is mainly the relationship between geoelectric (potential or resistivity) changes and earthquakes.One famous case study is the resistivity Correspondence to: Q. Huang (huangq@pku.edu.cn)changes before the 1976 Tangshan earthquake (Zhao andQian, 1978, 1994).After the early report of geoelectric potential changes before the 1969 M s = 7.4 Bohai earthquake in China (Central Earthquake Working Committee, 1969) and the Kamchatka earthquakes in Russia (Sobolev, 1975), three physicists in Greece (Varotsos, Alexopoulos, and Nomicos, abbreviated as VAN hereafter) started to make a series observation test of geoelectric potential changes and developed the so-called VAN method, which is a short-term earthquake prediction method using the Seismic Electric Signals (SES) extracted from the continuous observing data of geoelectric potential changes (Varotsos and Alexopoulos, 1984a, b).Besides geoelectric changes, there are independent reports of geomagnetic changes before earthquakes.More attention has been given to the frequency at Ultra-Low-Frequency (ULF, up to 30 Hz, which is different from the definition of the International Telecommunication Union; ITU) (Fraser-Smith et al., 1990;Kopytenko et al., 1993;Hayakawa, 1999;Hattori, 2004).The magnitude of the possible reported anomalies of ULF geomagnetic changes tends to decrease with time, partially due to the increase of the resolution of magnetometers (Johnston, 1997).The ULF anomalies seem to follow an empirical relation between earthquake magnitude and epicentral distance (Hattori, 2004).Simultaneous changes of EM signals at various frequencies before the 1995 Kobe earthquake were reported (Nagao et al., 2002).It should be mentioned that there are also some interesting reports about the secondary EM phenomena before earthquakes (Hayakawa, 1999).
Because numerous observed EM data at a broad range of frequencies support the existence of earthquake-related EM phenomena (Zhao andQian, 1978, 1994;Fraser-Smith et al., 1990;Kopytenko et al., 1993;Park et al., 1993;Qian and Lin, 1994;Johnson, 1997;Hayakawa, 1999;Nagao et al., 2000Nagao et al., , 2002;;Du et al., 2002;Uyeda et al., 2002;Huang and Liu, 2006;Zhang et al., 2006;Gao et al., 2010;Tang et al., 2010; Q. Huang: Rethinking earthquake-related DC-ULF electromagnetic phenomena Huang, 2011), the related-work has been paid much attention in the study of seismic hazard.Remote EM networks have been established in some countries such as Japan, Greece, and so on (Uyeda et al., 2004;Varotsos, 2005).Recent work has indicated that appropriate EM array observation would provide much more useful information of underground crust activity, compared to that of conventional EM network observation.For example, EM changes possibly related to the volcano-seismic activity in the summer of 2000 were detected by an EM array at Izu peninsula (Uyeda et al., 2002).
However, seismo-EM study up to now is mainly empirical, which means a kind of experience-based approach.Here, the term of "empirical" is a well-accepted term in the community of earthquake forecast and prediction, which means that no matter what the inside physics of the observed phenomenon is, one tends to make a prediction based on the past experience learned from the observed phenomena.Note that such experience is not sufficient enough to ensure that one prediction is right, because the above earthquake-related EM phenomena have not been fully understood due to the complexity of seismogenic process and inhomogeneity of underground resistivity structure.For example, there are some active debates on either the VAN method or other seismo-EM methods; even some successful case studies are reported (Lighthill, 1996;Uyeda and Kamogawa, 2008).The main concerns include the physical generation mechanism of the reported EM signals (Park et al., 1993;Johnson, 1997;Huang, 2002;Varotsos, 2005) and the selectivity phenomenon of SES (Orihara et al., 2002;Uyeda et al., 2004;Varotsos, 2005).Thus, it is important for contempory empirical studies on earthquake-related EM signals to move to a physics-based study, which would be helpful for understanding earthquake-related EM phenomena and strengthening their potential applications.
For simplicity, this study focuses only on the earthquakerelated DC-ULF EM phenomena, although the frequency of the co-seismic EM signals discussed in Sect.2.1.2is at Extremely-Low-Frequency (ELF, 30-300 Hz, which is different from the definition of ITU) band.Thus, unless stated explicitly, earthquake-related EM phenomena in this paper means the frequency of DC-ULF.First, I discuss some controversial or challenging problems in seismo-EM phenomena, e.g. the physical generation mechanism, the propagation and selectivity of seismo-EM signals, the methodology of extracting weak signals from noisy background, etc.Then, following some recent attempts to understand the above problems, as a potential physics-based approach, I present an integrated research scheme, taking into account the interaction between observation, methodology, and physical modeling.
2 Some controversial/challenging problems in seismo-electromagnetics

Generation mechanisms of earthquake-related EM signals
There are some possible physical models of generating earthquake-related EM signals, e.g.electrokinetic effect (Mizutani et al., 1976;Johnson, 1997), piezoelectric effect (Nitsan, 1977;Ikeya et al., 1997;Huang, 2002), stressstimulated charges (Varotsos, 2005), stress-activated positive holes (Freund, 2010), etc.Among the above generation mechanisms, electrokinetic effect would be one of the most favorite mechanisms in terms of literature (Johnson, 1997).However, how to explain the observed phenomena in the field condition is a challenging problem deserving further study.In this section, as an example, I'll discuss the following two controversial issues related to the above problem: one is the effective piezoelectricity in the field condition, and the other is the co-seismic EM signals.

Is the effective piezoelectricity in the field condition sufficient enough to produce observable signals?
Regarding the generation mechanism of earthquake-related EM signals, piezoelectric effect is one of the candidate mechanisms because one can easily connect the EM signals and the stress release of an earthquake.Such an idea has been applied to the explanation of earthquake lightning (Finkelstein and Powell, 1970).However, piezoelectric effect in the earth is controversial (Tuck et al., 1977;Ghomshei and Templeton, 1989).On the one hand, piezoelectric effect may be dismissed because it should be screened by the compensating charges in the conductive earth (Tuck et al., 1977); but on the other hand, a model study considering the role of the compensating charges has indicated that piezoelectric effect could be a plausible explanation of earthquake-related EM signals (Ikeya et al., 1997).Such an idea seems to be supported by some laboratory rock experiments (Yoshida et al., 1997;Sasaoka et al., 1998).The problem is that if one simply extends the laboratory results to the field case, an order of 10 −2∼−3 effective piezoelectricity is required to produce observable signals in laboratory experiments (Sasaoka et al., 1998), which means a very strong preferred orientation of piezoelectric crystals that is difficult to satisfy in the field condition.Therefore, whether or not the effective piezoelectricity in the field condition is sufficient enough to produce observable signals becomes an interesting issue deserving further study.Huang (2002) presented a model based on the role of the compensation charges and the elastic dislocation theory of a fault to explain the electric signals during rupture and evaluate the effective piezoelectricity in the scale of the field condition.Figure 1 gives an example of evaluating effective  piezoelectricity of a rock, which consists of piezoelectric crystals.For simplicity, the effective piezoelectric polarization is estimated without considering the interaction among elementary polarizations of crystals.If the crystals in the rock are completely random, the number of orientated grains would be n 1/2 , considering the statistical fluctuation, where n is the total number of piezoelectric crystals.Assuming the volume fraction of crystal grain is η, and N is the total number of grains, n = N η.Then the effective piezoelectric constant, α eff , would satisfy, where α is the piezoelectric constant of a single crystal.In the scale of the field condition, the number of crystals (n) becomes very huge, and the effective piezoelectric constant, α eff , would tend to zero.Thus, the effective piezoelectricity could be negligible.However, if there are some preferred orientations of piezoelectric crystals in rocks (Fig. 1), the effective piezoelectric constant, α eff , should not satisfy Eq. ( 1).Instead, it should be evaluated using, where ς is the fraction of preferred orientations of crystals.
Combining the above effective piezoelectric constant, α eff , given by Eq. ( 2) and the EM fault model presented by Huang (2002), we found that even a weak preferred orientation of piezoelectric crystal grains could lead to observable signals.This is due to its weak effective piezoelectricity in the field condition, e.g. with an order of 10 −6 , which is much weaker than that required in laboratory experiments (10 −2∼−3 by Sasaoka et al., 1998).

Which co-seismic EM signals can be expected?
Comparing with the numerous reports on pre-seismic EM signals, co-seismic DC-ULF/ELF EM signals are not reported frequently, although there are some related reports (Nagao et al., 2000;Skordas et al., 2000;Tang et al., 2010).However, it was found that all of the observed signals are most likely synchronous with the arrivals of seismic waves.Huang (2002) presented an EM fault model that considers the compensation of piezoelectric effect and the dislocation theory of a fault as a possible explanation of the observed coseismic electric signals.Besides such signals accompanying those of seismic waves, is there any other kind of co-seismic EM signals?
Taking into account electrokinetic effect, Pride (1994) derived a set of macroscopic equations governing the coupled seismic and EM waves in fluid-saturated porous media.Haartsen and Pride (1997) presented a numerical technique to solve the above governing equations and investigated the properties of the coupled seismic and EM waves in a fluid-saturated stratified porous medium.As an improvement of the above numerical technique by using a new straightforward evaluation method of the generalized reflection and transmission coefficients, Ren et al. (2010a, b) extended the method of calculating synthetic seismogram from Chen (1999Chen ( , 2007) ) to the numerical simulation of coupled seismic and EM waves in multi-layered porous media.As an example, Fig. 2 shows the results of applying the above technique to a two-layer model with an explosion source near the free surface.For this model (Fig. 2a), a sand-sandstone interface is located at 100 m depth.Receivers are located at the free surface.The source time function adopted here is Ricker wavelet with a dominant frequency of 100 Hz. Figure 2b indicates that, besides the clear EM signals accompanying with seismic waves, there are also weak P-to-EM conversion waves that are radiated from the interface and travel at EM wave velocity.Note that Fig. 2 gives just a conceptual example to demonstrate the coupled seismic and EM signals in a layered porous medium.Thus, if such P-to-EM conversion waves of natural earthquakes can be detected in the field condition, these signals may have some potential application in earthquake early warning system because they travel much faster than the velocity of P waves.
It should be mentioned that besides the above EM signals accompanying those of seismic waves and/or being radiated from the interface, there might also be EM signals directly associated with the rupture process.However, so far no co-rupture EM signals have been detected for natural earthquakes.Unfortunately, there is no known theoretical work  considering both the coupled mechanisms (e.g.electrokinetic effects) and the dynamic rupture process of earthquakes.If the model of Ren et al. (2010a) can be further improved after considering the dynamic rupture model of earthquakes, it may provide some theoretical feature of the possible corupture EM signals.In fact, we are starting to consider some dynamic rupture modeling of an earthquake in the theoretical coupling model.The related result will be reported independently elsewhere.

Propagation and selectivity of earthquake-related EM signals
The selectivity phenomenon of SES (Orihara et al., 2002;Uyeda et al., 2004;Varotsos, 2005) is an interesting topic attracting some hot discussions (Lighthill, 1996;Uyeda et al., 2004;Uyeda and Kamogawa, 2008).Varotsos et al. (1998) proposed a conductive channel model to explain the SES selectivity.Their numerical simulation indicated that there is some amplification effect of EM signals above the upper end of a conductive channel (Sarlis et al., 1999), indicating that the underground conductive channel could be an explanation for the SES selectivity.This conclusion has been further supported by the recent numerical simulation based on a 3-D finite element method (FEM) on the SES selectivity of the 2000 Izu earthquake swarm (Huang and Lin, 2010a).It should be mentioned that the 3-D FEM simulation indicated that the above amplification effect above the upper end of a conductive channel exists only for ULF band.
As another possible explanation of the SES selectivity, some laboratory analog experiments on a geographical scale model and a waveguide model were developed to simulate the propagation of seismic EM signals (Huang andIkeya, 1998, 1999;Huang, 2005).These experimental results showed that the geographical effect, such as the distribution of ocean and land, may lead to some aspect of the selectivity phenomenon.Such experimental approaches (e.g.Huang and Ikeya, 1998) have also demonstrated the feasibility of long range EM signals in the band from very-low-frequency (VLF) to ELF, as well as enhanced the reliability of field observations in the VLF-ELF band as commented by Tzanis and Vallianatos (2001).
After investigating the possible influences in numerical simulation of the SES selectivity, Huang and Lin (2010b) showed that the selectivity phenomena of SES may be caused either by surface inhomogeneity of resistivity or by underground conductive channel, depending on the detailed electric structure and contrast of resistivity in the model.Figure 3 shows a sketch model, which consists of ocean-land distribution, land resistivity inhomogeneity, and underground conductive channel.Figure 4a gives the electric field distribution on the surface obtained from the 3-D FEM calculation based on the model given in Fig. 3. Figure 4b shows the results after adding to Fig. 3 a horizontal conductive channel at a depth of 10 km, connecting the dipole current source and the vertical underground conductive channel.One can see from Fig. 4b that the electric field can be enhanced significantly around both the surface resistivity boundary and the upper end of the vertical conductive channel.Thus, one may conclude from the above model results that the double contributions from surface inhomogeneity of resistivity and underground conductive channel may play an important role in the selectivity phenomena of SES.It should be mentioned that besides the above model calculations, the evidence from the field observations of the real underground structure is also very important.Fortunately, some field results from the detailed EM sounding study around the Ioannina SES sensitive station in Greece (Eftaxias et al., 2002;Balasis et al., 2002Balasis et al., , 2005) ) seem to support the above conclusion from the model results.

Extracting weak earthquake-related EM signals from noisy background
At the current empirical stage, it would be important to make the reliability and statistic analyses of the above phenomena, similarly to those that have already been done for seismicity precursors (Huang, 2006).Another important issue is how to reveal earthquake-related EM signals from background noise, which requires further study on methodology.
As an example of methodology of earthquake-related EM phenomena, we provide a new seismomagnetic data processing method by combining principal component analysis (PCA) (Hattori, 2004(Hattori, , 2006;;Han et al., 2009;Huang, 2011) and geomagnetic diurnal variation analysis (Han et al., 2009(Han et al., , 2011;;Huang, 2011).After applying the PCA method to the geomagnetic diurnal data, one can investigate the temporal variations of the contribution from each principal component (Huang, 2011).As a case study of the 2000 Izu earthquake swarm, Fig. 5 shows the temporal variation of contribution of each principle component to the geomagnetic diurnal variation (see Huang, 2011 for the principle of the methodology) at Seikoshi (SKS), Mochikoshi (MCK), and Kamo (KAM) stations.A clear enhancement of the contribution from the second principle component was revealed about one month before the Izu earthquake swarm.The result indicated that the methodology of applying PCA to geomagnetic diurnal variation is effective in revealing possible earthquake-related signals from background disturbances (Han et al., 2010).The detailed PCA analyses of the geomagnetic diurnal variation associated with the Izu earthquake swarm will be given in a separate paper elsewhere.
It should be mentioned that sometimes the wavelet analysis can be also useful for extracting weak signals from noisy background.For example, the wavelet analysis of geomagnetic data at SKS, MCK, and KAM stations indicated that an unusual enhancement of geomagnetic signals would be a kind of phenomenon associated with the 2000 Izu Islands earthquake swarm (Han et al., 2011).

A physics-based approach
As discussed in the previous section, there are some controversial and/or challenging problems in seismo-EM phenomena due to the complicated seismogenic process.The recent attempts have also provided some useful information for understanding the above problems.However, because the seismo-EM study up to now is mainly empirical, a physicsbased study would definitely be helpful for understanding earthquake-related EM phenomena and strengthening their applications.As a potential physics-based approach, I present an integrated research scheme, taking into account the interaction among observation, methodology, and physical model.The basic idea is given in Fig. 6.The main approach includes the following key problems: 1. How to perform a reliable and appropriate observation with some clear physical quantities.Similar to other branches in geophysics, seismo-electromagnetism is also an observation-based discipline.Thus, in order to obtain convincing EM precursors, the fundamental observation data must be credible enough.Also, the observation should have clear physical meaning, which is useful for either repeatable testing or theoretical investigation.An appropriate observation not only requires the past experience data, but also should take into account some results and/or requirements of the physical model and the methodology studies.Sometimes, it would be important to design and perform some controlled observations to test the physical model or the proposed methodology.
2. How to develop a robust methodology to reveal weak earthquake-related EM signals from noisy background.As discussed previously, contemporary study on earthquake-related EM phenomena is mainly at an empirical stage.Past experience has shown that the earthquake-related EM signals are generally weak, especially compared with the background noises.Thus, developing robust methodology becomes a useful and challenging problem.If there is an appropriate observation, adopting some existed techniques in signal processing may be effective in extracting weak signals from noisy background, as has already been discussed in the previous section.Taking into account the empirical characteristics of earthquake-related EM signals, as well as some insights from the physical model study, developing new robust methodology should also be expected.Of course, in order to enhance the reliability of the revealed precursory signals, another important issue is performing the reliability analysis and rigorous statistic test, similar to those that have already been done for seismicity precursors (Huang, 2006).
3. How to develop plausible physical models based on theoretical analyses and/or laboratory experiments for the explanation of the earthquake-related EM signals observed in the field conditions.The empirical features of earthquake-related EM signals obtained from observational and methodological approaches can provide fundamental constraints and field validations for the physical models.On the other hand, the theoretical features predicted from the physical models may be useful for designing observation and developing methodology.In order to strengthen the understanding of the multi-physics process of earthquakes, it is important to promote theoretical work considering both the coupled mechanisms (e.g.electrokinetic effects) and the dynamic rupture process of earthquakes.Furthermore, designing and performing some controlled field testing with multi-geophysical observations may provide an opportunity to test the integrated research scheme given in Fig. 6.
Being a kind of potential precursor of earthquakes, the earthquake-related DC-ULF EM phenomena have been paid considerable attentions and play an important role in the study of seismic hazards.However, at the current empirical stage, there are some active debated problems deserving further study.After reviewing some recent works in an attempt to understand the above controversial problems, I present a potential physics-based research approach that takes into account the interaction among observation, methodology, and physical model.It is expected that such an integrated research scheme should be effective in understanding the physics of earthquake-related EM phenomena, as compared to the previous empirical approach.
(a) A sketch figure of the crystals in the rock with a weak preferred orientation.The bold arrow indicates the orientation of each crystal.(b) An equivalent orientation of tals in the rock, which is the same as that given in (a).The thick bold arrow gives the ivalent orientation of crystals in (a).

Fig. 1 .
Fig. 1.(a) A sketch figure of the crystals in the rock with a weak preferred orientation.The thin bold arrow indicates the orientation of each crystal.(b) An equivalent orientation of crystals in the rock, which is the same as that given in (a).The thick bold arrow gives the equivalent orientation of crystals in (a).
Figure 2 (a) A two-layer model with an explosion source near the free surface.(b) 6 Seismograms and electro-seismograms of the two-layer model.Direct waves and free surface 7 reflections are omitted in order to show coupled waves clearly.The P-to-EM conversions 8 were amplified by a factor of 500.9 10

Fig. 2 .
Fig. 2. (a) A two-layer model with an explosion source near the free surface.(b) Seismograms and electro-seismograms of the two-layer model.Direct waves and free surface reflections are omitted in order to show coupled waves clearly.The P-to-EM conversions were amplified by a factor of 500.

8 Fig. 3 .
Figure 3 A sketch model consists of ocean-land distribution, land resistivity inhomogeneity, 3 and underground conductive channel (from Huang and Lin, 2010b).(a) Model distribution on 4 the surface.The star and circle indicates the projection on the surface of the dipole current 5 source and the underground conductive channel, respectively.(b) A cross-section including 6 the dipole current source and the underground conductive channel.7 8

Figure 4 7 Fig. 4 .
Figure 4 Electric field distribution on the surface (modified after Huang and Lin, 2010b).(a) 3 Results of the model given in Fig. 3. (b) Results of the model after adding in Fig. 3 a 4 horizontal conductive channel at a depth of 10 km, which connects the dipole current source 5 and the vertical underground conductive channel.6 7

Fig. 5 .Fig. 6 .
Fig. 5. Temporal variation of contribution (percentage) of each principle component to the geomagnetic diurnal variation, which is obtained by applying Principle Component Analysis (PCA) to the geomagnetic diurnal fitting data of X component of geomagnetic field at Seikoshi (SKS), Mochikoshi (MCK) and Kamo (KAM) stations.The number of earthquakes around the stations and the magnetic K p index are also given for reference.The vertical dashed lines indicate two M > 6.0 earthquakes in the investigated spatio-temporal range.