Abstract

Seismic waves are generally observed through the measurement of undulating elastic ground motion. We report the remote detection of the Earth's electric field variations almost simultaneously with the start of fault rupturing at about 100 km from the fault region using a special electric measurement. The rare but repeated detection indicates that the phenomenon is real. The characteristic time of diffusion is almost instantaneous, that is, less than 1 second to travel 100 km, more than ten times faster than ordinary seismic P wave propagation. We suggest that the measured electric field changes are produced by the electrokinetic effect through increased pore water pressure of the seismic pulse. It is also suggested that the long range propagation is due to the surface wave mode confined near the interface of the different conductivity. The length scale of the finite strength of the electric field is 16 km, 160 km for electric conductivity of 0.01, 0.001, Sm−1, respectively. This phenomenon suggests a new seismic sensing method and a new earthquake early warning system providing more seconds of lead time.

1. Introduction

Seismic waves are emitted from the slipping part of rupturing faults during an earthquake, and seismometers sense the induced ground elastic motion. The seismic wave speed is about 6 km/s for the primary (P) wave and 3.5 km/s for the secondary wave [1, 2] in the upper crust. There have been many reports on the electromagnetic phenomena of earthquakes and volcanic eruptions (e.g., [37]), but only a few papers have dealt with the field observation of electromagnetic variation that started just after the rupture of the fault [823]. Belov et al. [8] detected associated electromagnetic fields about 10 seconds before the seismic wave arrival and 4 seconds after the origin time. The ULF band field changes were detected some 0.3–0.8 seconds before P wave arrival and about three seconds after the origin time during the Izmit earthquake [13, 24] and between the P wave arrival and the origin time during the 1995 Hyogoken-Nanbu earthquake [11]. The VLF band pulse-like signals were detected about 1,000 km from the fault at about the origin time of the great earthquake east of Hokkaido on October 4, 1994 [10].

Recent application of seismoelectronics to exploration surveying has resulted in finding and confirming the coseismic electromagnetic field [1417]. However, we need more data to understand the phenomenon to build a consistent generation model and to develop application technology based on this phenomenon. The most interesting problem is to investigate whether the signal starts at the beginning of the fault rupture and to understand the generation and the propagation mechanism. The propagation speed of the signal is another important point to be revealed, especially from the viewpoint of application to the earthquake early warning (EEW) system [2832]. The greater speed of this electric field variation is expected to provide several more critical seconds to issue the warning before the arrival of the hazardous S wave, which could contribute to decreasing the damage to life and property during an earthquake.

2. Observation Method

We have been observing electric field changes at 13 ground electric field monitoring stations in the central parts of Japan in relation to the research and development of the prediction method for earthquakes and volcanic eruptions (Figure 1) since 1989. The detection sensor is a special antenna made of a vertical casing pipe in boreholes of depths ranging  m [20, 21]. The observation frequency bandwidth is DC (0–0.7 Hz), ULF (0.01–0.7 Hz), and ELF/VLF (1–9 kHz). For frequencies lower than the ULF band the antenna works as an immersed ultralong electrode. We used the three-channel analog recorder to monitor the raw waveform of the DC and ULF bands and the envelop of waveform of the VLF band. The running speed of the recorder is 1.5 cm/h. The system is very sensitive to crustal activities such as earthquakes and volcanic eruptions and is highly robust to both meteorological and urban noises. The evaluation of the system was based on observations at different geological situations at 13 sites in central Japan and for more than 13 years (e.g., [11, 25, 33]).

To investigate waveform characteristics of precursory anomalous signals of earthquakes and volcanic eruptions, we made temporal high-speed analog recordings. The recorder’s original role was as a continuous monitoring seismograph with a chart speed of 3 mm/s. The temporal electric field observation was conducted for about 27 months during June 1994 to December 2002 at two sites of the borehole measurement network, Hasaki and Chikura (Figure 1).

3. Coseismic Changes

A moderate earthquake of magnitude 4.9 occurred at a depth of 57 km in one of the most active seismic zones in the south of Ibaraki Prefecture, north of Tokyo, Japan, on 14 June 2002 (Figure 1). Figure 2(a) shows the electric monitoring record of low speed at Hasaki, on the day of the earthquake. We can see that at three time periods in the DC and ULF bands (black and blue) there were three groups of anomalies, A1, A2, before the occurrence of the earthquake, and B at the time of occurrence with very different indications from the normal state record (Figure 2(b)).

We fixed the recording parameters in our analog recording after a testing period of approximately two months. We analyzed the level and its variation by taking account of several background noises of natural and human-made origin during the normal state in terms of crustal activity. One of the most important efforts in electromagnetic observation is to decrease cumbersome noises to assure high value of the S/N or to find a method to discriminate signals embedded in large background noises (e.g., [34]).

As to our observation system, the unique noise is the surge of electric current at the time of atmospherics. The occurrence of very large atmospherics is rare and easy to discriminate from normal noises, so it is not difficult to select amplification parameters and make a zero adjustment to mask the usual noises to start the observation in search of anomalies because of the high robustness in comparison to environmental noises. It was easy for us to pick up the present anomalous signals A and B shown here because of their peculiar pattern of evolution.

One of the anomalies (A1), lasting about 20 minutes, occurred about 1 hour and 40 minutes before the earthquake, and the second anomaly (A2), which lasted about 25 minutes, occurred about 30 minutes before the earthquake. The time scales of these variations that appeared on the two channels were 20 seconds ~2 minutes. These kind of anomalies have been argued that they are not the result of meteorological or geomagnetic causes but rather are precursory phenomena related to earthquakes [35, 36].

Precursory anomalous electric variations were also detected at the site in Chikura (Figure 1). The records in the DC and ULF bands are, however, quite different from those at Hasaki. None of the individual anomalies was simultaneous with other anomalies, and none was coherent; moreover the anomalous appearance persisted for several weeks at Chikura, with maximum ULF band anomalies appearing on July 6, 2002. Later in the present work, we will touch upon these phenomena when we discuss the propagation channel of the remotely detected electric signals.

At the first stage of the analysis of the phenomenon, the B anomalies (Figure 2(a)) at Hasaki were thought to be well-known coseismic electric variations observed and reported in numerous papers [9, 10, 12, 1423]. They correlated in general with the ground surface fluctuation induced by the seismic wave arrival. Moreover, seismograms and electrograms at the observation site correlated well with each other [10, 1517]. Theoretical and observational studies indicated that there were linear relations between the electric field and seismic velocity [10] or acceleration [17], with the coefficient (transfer function) depending on the kinetic property on the site.

During this observation period, we simultaneously made a continuous high-speed analog recording. The time marker pulse is overwritten on the signal trace every second. The clock of the recorder was calibrated every minute using the radio clock with a resulting accuracy to within 1 ms. The high-speed recording is shown in the upper part of Figure 3(a).

Signal B was quite different from those coseismic electric variations correlated with the arrived seismic fluctuation of the ground. The initial phase of B can be clearly seen from 42 m 50 s ( in Figure 2(a)) to 43 m 5 s (P) in the record. is the origin time of the focal parameter for the earthquake determined by the National Research Institute for Earth Science and Disaster Prevention (NIED). The peak of B at around 43 m 15 s in Figure 2(a) is not recorded due to signal saturation. The signal started at the origin time of the earthquake, not at the time of seismic wave arrival. Moreover, the frequency band is quite different: in the case of signal B, the band is in the 0.01–0.7 Hz range, while the seismogram shows several Hz, as suggested in the strong motion record described below.

We also compared the electric signal waveform at Hasaki with the seismogram there. The lower part of Figure 3(b) is the assumed seismic record at Hasaki using the seismograph at the nearest site of the strong motion network (K-NET) of NIED [37], which was Choshi. We took into account the difference of about 11 km in the hypocentral distance between Choshi and Hasaki when building the assumed seismogram at Hasaki. The horizontal time axis was shifted by the P wave traveling time difference,  s, by adopting 6.1 km/s for P wave speed. Anomaly B is different from the ordinary coseismic variation moving with the seismic wave as recorded and studied so far [9, 10, 12, 1423]. This is the first report of the detection of this kind of electric variation, that is, correlation with the seismic wave starting from the beginning of fault rupture, in our long history of electric field measurement by using a borehole antenna.

Comparison of the electric field variation and the observed seismic waveform leads us to the following points. (1)The starting time of the anomaly (11 h 42 m 49.0 s (±1 s)) is close to the origin time of the earthquake (=11 h 42 m 50.0 s (±1 s)).(2)The strength of the electric field increased gradually until P wave arrival, when it peaked.(3)The electric field strength decreased after P wave arrival.(4)The electric field strength increased very sharply from S wave arrival and decreased very slowly as the seismic wave diminished.(5)The phenomena did not occur coincidently with the earthquake, but has an intrinsic relation to the seismic wave. (6)It is indicated by point (1) that the effective characteristic time of diffusion from source to receiver of some 100 km is less than 1 second, corresponding to an effective spreading rate to be greater than 100 km/s.

4. Statistical Characteristics

Another similar variation was observed for the earthquake (M4.9, 13 h 04 m, 35.815 N, 140.914 E, depth 32.2 km) that occurred on October 16, 2002, near Hasaki (E2 in Figure 1), as shown in Figure 4. In this case the electric field started at October 16 13 h 04 m 34.2 s (±1 s) compared with the origin time of October 16 13 h 03 m 34.5 s (±1 s), corresponding with the beginning of fault rupturing, but the signal waveform did not increase as the P wave approached the site. At around the time of the S wave arrival the electric field strength increased sharply. There are many similarities between the two cases, and thus we investigated the reappearance of the phenomena using all available high-speed recording data.

We had access to high-speed recordings covering the 27 months from June 1994 to December 2002, giving data for 20 months at Hasaki and 7 months at Chikura. The whole set of data was checked for the reappearance of coseismic signals for moderately large earthquakes at the two observation sites. The seismic events to be compared were tentatively limited to those of magnitude M ≥ 4 and an epicentral distance less than 100 km from Hasaki or Chikura. There were 55 such earthquakes in this time period on the basis of the seismic catalog of NIED.

There were 3 cases of corupture anomalies in 42 earthquake events at Hasaki and 2 cases in 24 such events at Chikura. In the analysis we checked whether there were similar anomalies just around the origin time of each earthquake. The epicentral distance to detect the corupture anomalies ranged from 37 to 98 km. Five cases had the largest signals at S wave arrival. However, the correlation is not so high compared with those on October 16, 2002, as shown in Figures 2 and 4: there are no peaks at P wave arrival or there are no discernible changes in between the origin time and P arrival. The changes in each of the signals are found to be related to neither the seismic source mechanism nor the seismic magnitude under the limited volume of the earthquake ensemble. Further research is needed to clarify this point.

5. Discussion

Here we discuss generation mechanism and propagation passage for the observed electromagnetic signal associated with seismic events.

5.1. Generating Source

An extensive discussion on the generation mechanism on seismic electromagnetic signals (SES) (e.g., [34]) indicates that there are three dominant mechanisms. The first mechanism involves active sources such as crustal deformations, crack generation in the crust, and pore water movement (group 1); the second involves passive sources such as changes in the physical parameters of the crust, that is, electrical conductivity, magnetic susceptibility, and electric permittivity (group 2). The third involves apparent sources, that is, the observation instruments and the environment [34]. As a matter of fact the third source is generally dominant in the SES observation, except special measurement systems such as ours [33].

The active sources are of piezomagnetic, piezoelectric, nonclassical piezoelectric, electrokinetic, and induction effects in nature. The piezoelectric effect is induced by mechanical stress changes of rocks containing minerals, such as quartz. The stress changes cause electric polarization of the monocrystalline quartz grain without a preferred orientation. The piezoelectricity is very weak outside the radio frequency [34]. The piezomagnetic effect is induced by the deformation of minerals with residual magnetism through ferromagnetic inclusions. The orientation of inclusions is changed in the course of the rock. The random and anisotropic features of the orientation result in the small magnetic effect of the seismic wave scale. There is also another kind of electric polarization change through the dislocation (e.g., [34]). In dielectric materials, dislocations are charged, and the dislocations move when the materials undergo stress changes. These effects seem to play a minor role (see the appendix).

The present phenomena do not seem to be induced by passive causes, because the properties of the crust do not change quickly, for example, over the maximum several tens of seconds of a large seismic wave’s amplitude. External causes are also discarded on large robustness to environmental noises of the present observation [33]. Moreover, the order estimation of the SES strength due to several causes [34, 38] suggests that the dominant causes are the piezomagnetic effect [11, 13, 14, 39] and/or the electrokinetic effect [7, 1417, 4045] for small stress changes such at the seismic wave frequency (see the appendix).

5.2. The Electrokinetic Effect and a Simple Model

The strong correlation of electric field variation with P wave propagation leads us to infer that the generation zone is limited to around the seismic wave front not to the rupture zone. In the following two subsections, we present simple model to discuss the signal propagation path semiquantitatively, that is, using a homogeneous diffusion model [33, 45] and the induced model of seismic electrical signal through the seismic pulse using the two-layer model [34].

5.2.1. Homogeneous Diffusion Model

The diffusion type Equation (A.1) of the electric field in unbounded conductor deduced from the Maxwell equation was used to infer the electric field evolution [33, 45]. Assuming the conductivity and the streaming potential coefficient, electric field strength is expressed by the use of the spherical symmetric Green function. Using possible values of the parameters for the electrokinetic effect, we inferred the source time function for the characteristic signal waveforms at the time of seismic swarm and volcanic eruption using observed signal waveform [33, 45].

A calculation at any point indicated that the sharp initial rise of the signal is conserved even at a distance of 30 km from the source (see [33, Figure  13]) by assuming electrical conductivity  Sm−1 (the conductivity is around  Sm−1 in central Japan [46]). The amplitude decreases by two orders but still can be detected by measurement system having high enough resolving power and ratio, which is the case for the present special borehole measurement system.

On the other hand, at the long distance of 100 km from the source, the estimated time of diffusion [7] after the start of the fault rupture is 10 s and 100 s for the electric conductivities of and 0.01 Sm−1, respectively, which is much larger than the detection time of less than 1 second, so that we cannot think the propagation in the homogeneous media as a most possible one.

5.2.2. Two-Layered Model

A Russian group treated an electromagnetic field induced by a plane seismic wave by assuming a half space of layered media of different conductivity [43] in order to discuss the precursory SES observed by [8].

In the model of a finite angle of seismic wave incidence for the magnetic induction interaction, they indicate that the induced variation can be observed also around the ground surface of the conductivity contrast with the characteristic detection distance of 1.6 km–160 km [43] (see the appendix). The surface mode dissipates by power law and dominate at far distance. The treatment can be extended to other mechanoelectric processes such as piezomagnetic, piezoelectric, and electrokinetic effects.

In the appendix we formalize the case of electrokinetic effects to show the equation for the seismoelectric field with the source induced by the induction effect can be used for the case of the electrokinetic effect. We can see that the electric field far from the rupture source is attributed to the surface wave confined just at the ground surface induced at the propagating wave front. We can expect fast and long distance detection as 160 km for the value of conductivity of 0.001 S/m.

In the case of the shallow earthquake the ground surface can be assumed as the wave guide to explain very short time delay of detection compared with the origin time. And for deep event we can assume underground contrast as interface between the upper crust and deeper crust.

As to the waveform of the forerunner the wave signal is shown to have porality reversal at the moment that the initial pulse balances to the following surface mode. The point is contrary to the observed profile waiting further studies.

5.3. Non-Newtonian Effects

Recent developments in the field of seismo electromagnetics started from the Frenkel-Biot model [47, 48] treating the charged material of linear poroelastic characteristics (e.g., [14, 15] and references cited there in). Assuming a two-phase system of poroelastic media, Biot deduced that there is a third elastic seismic wave (a slow compressive wave) in addition to the ordinary compressive P and shear S waves. Pride completed the theory of two-phase media comprised of porous elastic materials and Newtonian fluid as a solvent. He used an averaging procedure to get macroscopic equations to show that conversion of seismic wave energy to electromagnetic variation at the interface of different rocks is most likely due to the mechanism of electrokinetic coupling [1517, 47, 48]. There are electromagnetic fields accompanying a seismic body and surface waves limited to the neighborhood of the seismic wave. The electromagnetic disturbances generated at the material interface can diffusive very fast, comparable to electromagnetic wave speed [17, 47, 48]. Revil and Jardani extended the theory by assuming elasticoviscous rheology of the generalized Maxwell type [15, 16]. In particular, they found that the seismoelectric variation converted from a seismic wave at the interface depends on the frequency and that it can have very large amplitude for the resonant seismic frequency with the relaxation frequency of the solvent of the generalized Maxwell law. It is also noteworthy that the DC components of the streaming coefficient may be large (see the appendix).

5.4. Interpretation of the Phenomena

On the basis of these research results, we are going to interpret the present phenomena. The very similar features of evolution of the electromagnetic field to those of the elastic compressive wave suggest that the present electric field variation is caused by seismic coupling induced by some kind of sources. The very early arrival of the variation led us to infer that the direct or reflected wave makes a large contribution to that variation. The most probable mechanism is that the electromagnetic field is induced by the conversion from a seismic wave originating at the substance interface caused by electrokinetic coupling at the interface. The coupling can be amplified by the non-Newtonian effect of the solvent. The interface includes the air/ground surface.

We suppose there are additional propagation mechanisms awaiting further research. Regarding the unstable detection of the signal, we assume the heterogeneous crustal property determines the streaming potential parameter along the path of the seismic propagation to the receiver, but the rarity is difficult to explain.

In this frame of reasoning, we may interpret that the precursory anomalies A1 and A2 are possibly due to a preslip at the fault [49, 50] several tens of minutes before the earthquake occurrence: the preslip of the porous media causes an electromagnetic field through the electrokinetic effect, and the electromagnetic field is propagated by diffusion to the observation point via some kind of efficient propagation path.

6. Possible Application

From a practical viewpoint, the phenomena suggest the possibility of a new measurement method in seismology using electromagnetic sensors in addition to the usual seismometers. For example, we may be able to approximate that the first detection time is the earthquake origin time, thus eliminating one of the unknowns in the focal parameters determination.

In the earthquake early warning (EEW) system available now in Japan [29, 30, 32, 51], the focal parameters are determined by the use of P wave arrival data at a small number of points. The practical use of this information started in October 2007 [29], with the implementation of particular countermeasure systems such as the emergency shutdown of gas fittings and home electrical devices [51]. One of the shortcomings of the EEW is that it cannot provide very much time for alerting people of an impending earthquake in the case of an earthquake of shallow depth just under a city. If we could find a new electromagnetic instrument to sense the seismic wave front using the present phenomena, we could gain several seconds before the S wave arrival, which could decrease the level of earthquake disaster considerably.

7. Concluding Remarks

Here we have presented observational data showing that there are seismoelectric field variations possibly converted at the physical interface of the different geology in the path of seismic wave propagation. The effective detection speed on the order of 100 km/s or possibly electromagnetic speed is much larger for the diffusion spreading for the realistic conductivity of 0.1–0.001 Sm−1. This phenomenon may be used in the future development of the EEW technology, producing more lead time before the arrival of the big shock, compared with the usual monitoring that relies on the elastic wave detection.

Further studies are needed before we can develop instruments able to detect the electric field variations almost every time when an earthquake occurs. It has been suggested that the ultralong electrode measurement (ULEM) [33, 35, 52] method may provide a useful sensor to solve this problem. Future efforts will be aimed at increasing the reproducibility of the observations as well as at understanding the generation and propagation mechanisms of the phenomenon.

Appendix

Assuming quasistationary cases of the frequency band of interest for phenomena of the earth crust, we can neglect the displacement current. In this case, the electric field obeys the following diffusion equation: (see, e.g., [7, 15, 33, 34, 38, 43, 45]). Here, , is the magnetic permeability and the electrical conductivity. is the source electric current and can be decomposed into where , , , and are the electric current due to the piezomagnetic effect, piezoelectric effect, and electrokinetic effect and the induced current by the motion of the medium in the geomagnetic field, respectively [34].

The ratio is 0.07 (conductivity  S/m, ), 0.7 ( S/m, ), 0.7 ( S/m, ), and 7 ( S/m, ). Two effects are comparable. The ratio is 30 ( S/m), 3 ( S/m) denoting that the piezomagnetic effect is lager for larger conductivity and comparable for the smaller conductivity.

In the case of the electrokinetic effect as a driving source, the electrokinetic current is written as where and are the coefficient of the streaming potential and the pore water confining pressure, respectively [41, 42].

 The phenomenological parameters of the streaming potential had been known by field observation or laboratory experiment to be 10−6 to 10−8 VPa−1 (e.g., [15, 34]). Pride [14] started to make a kinematic approach to express the frequency-dependent parameter in terms of the electrochemistry of mineral/water interface such as zero potential and surface conductivity and textural parameters of poroelastic rock and solvent [7, 15, 16, 47, 48]. The dynamic streaming coefficient is expressed as where, , , , and are electrical permittivity, zeta potential, fluid viscosity, and fluid conductivity, density of pore water. And are textural parameters (see [15]) of the porous media. Critical frequency is of order of 10 kHz. For the low frequency limit, In case of a porous material saturated with the generalized Maxwell fluid, the coupling coefficient is where is the peak relaxation time of the Cole-Cole distribution of the relaxation in a generalized Maxwell fluid and is a constant in the form of the distribution (see [15]). The value for crude oil is shown to be several orders larger than the Newtonian fluid of order 10−5 V/Pa [15]. Moreover the coefficient is also several orders larger for the excitation by resonance frequency of the viscoelastic fluid.

We use the analytical method for the present purpose of gaining an intuitive understanding of electromagnetic field evolution via remote observation. Gershenzon et al. [43] used a very simple two-layered model to explain the magnetic precursor preceding seismic wave arrival observed by Belov et al. [8]. They assumed the induction mechanism for the electromagnetic field generation. Here we show the mechanism of the electrokinetic effect. The pore water pressure for the porous media for the slow process is expressed as follows: where , , the porosity, and the volume strain. And where , , and , are the bulk moduli of the rock matrix, dry porous rock, and pore water, respectively, is the porosity, and is the volume strain [34, 47, 48].

We treat a two-dimensional problem with the wave propagating vector situated in the plane normal to the -plane with the incidence angle (see Figure 5). We take the excess pressure pulse, , associated with the seismic wave, as follows: where is the pulse height, is the width, and is the seismic P wave velocity. The governing equation of the fluctuation part of the electric field, , in the th layer is where , , is the velocity of the electromagnetic wave in the Earth, and is the conductivity of the th layer.

The expression is the same as those of the magnetic field [43], except for the form of the driving force term, . The boundary condition in both directions from infinity is The condition at the boundary is that the tangential component of and the normal component of electric displacement of are continuous. For the sake of simplicity, the permeability is assumed to be the same in the two layers. In this case the boundary condition is The equation for the electric field is reduced to almost the same as that for the magnetic field to treat the induction driving force [43]. The final solutions are where , is the Heaviside function, , and are integrals of functions including , . The , , and are expressed in similar ways in [43]. The right-hand term is composed of two terms, and . The first term is composed of and exp . These two terms are limited to the wave front. On the other hand, the term cos is confined to the layer boundary propagation in the forward direction. The second term has a finite value at the far field.

The spatial length, , which shows the finite strength of the electric field, is to correspond to as . The length is 1.6, 16, and 160 km for the conductivity of 0.1, 0.01, and 0.001 S/m. Observation far from the wave front corresponds to , which extends in the forward direction of wave propagation and in the opposite direction. At the far point the electric field strength is shown to decay in the power law as , .

Acknowledgments

The authors thank R. Ikeda, Y. Hagiwara, H. Takahashi, S. Uehara, T. Katayama, S. Yamane, N. Fujii, and M. Ukawa for support and encouragement and E. Yamamoto and Y. Yamauchi for technical support. They also thank NEC Corporation, Nippon Steel Welding Products & Engineering Co., Ltd., and Tics. Inc., Sankousha Corporation, Franklin Japan Co., Ltd., Izu-Oshima Onsen Hotel, Teikoku Oil Co., Ltd., Nakazato Village of Gifu Prefecture, Kofu City, and Shizuoka Prefecture for support. They would like to express their sincere thanks for the careful reading and valuable suggestions regarding the original and revised paper by Professors Andre Revil, Naum I. Gershenzon and the Editor Professor Gary Egbert.