Since the 1966 M7.2 Xingtai earthquake, permanent observatories have been set up in some of the major earthquake zones in China to investigate earthquake precursor of geoelectric resist-vity and current field^[1]. Precursors to the 1976 M7.8 Tangshan earthquake have been found, e.g. decrease of geoelectric resistivity^[2~5], geo-resist-vity impulse related to tidal force M₂ wave^[6], and geo-current field changes related to tidal force MS_f^{[7, 8]}. A model of short-term earthquake precursor of tidal force has been established in accordance^{[9, 10]}. The PS100 anti-interference high-precision geoelectric measurement system has been developed to overcome the instability of electrode polarization, the interference of the natural electromagnetic field, and especially the interference of random current due to industrialization, which will prevail on the globe. This system utilizes the Code Division Multiple Access (CDMA) technique and the concept of black box to the earth-resistivity measurement instrument. It applies convolution and deconvolution to the encoded input time series and the observed output time series, and identifies the value of georesistivity with eliminated electromagnetic interference^[11]. Specification tests have been conducted by China National Institute of Metrology regarding the anti-interference ability of the PS100 system. APS100 system had measured the resistance of a 0.01% measuring resistance, on which interference voltage from a signal generator (as source of the interference) was applied. The results¹⁾ indicate that when the signal-to-noise ratio was between 20dB and 0dB, both the difference between the detected and true values and the standard deviation of multiple measurements were better than 0.1% (between 0.01 % at 20 dB and 0.1% at 0 dB).

1) No. of certificate: DCcm 2004-0012. and the calibration/measurement described in this certificate, as on item of the Calibration Measurement Capabilities (CMC) of National Institute of Metrology (NIM) in DC and low frequency (IMHz and below) of electromagnetic metrology field, has been listed in the Appendix C of MRA and can be found in the BIPM's website: http://kcdb.bipm.fr/BIPM-KCDB/AppendixC/search.asp?MET=EM

Four sets of PS100 systems has been recording at 4 stations in Sichuan and Yunnan (MN station at 28.532°N, 102.169°E, HG station at 26.544°N, 101.928°E, YM station at 25.690°N, 101.861°E, and LJ station at 26.837°N, 100.245°E) since April 2004(Fig. 1)

图 1

Fig. 1

Fig. 1 PS-100 HRT station deployment and the distribution of epicentres of earthquakes, before which the HRT wave were recorded (The case number in brackets is the same in Table 2)

In our two internal reports to the China Earthquake Administration in 2003, we pointed out that along the middle segment of the North-South Earthquake Zone (within 200~300 km from the MN station), it is likely that a huge earthquake occur in the next two to three years^[12]. Because of this judgment, we deployed four PS-100HRT stations under an emergency funding from the China Earthquake Administration to confirm the HRT earthquake prediction method. Since the April of 2004, these four stations have been running continuously to capture the HRT waves^[13]. The most northern station MN is 288 km away from the epicentre of the Wenchuan earthquake (Fig. 1). The station and earthquake epicentre are within the boundary zone (about 105°E) between the East and the West Interaplate Blocks in China. The four PS100HRT system recorded geo-resistivity with a standard deviation of less than 0.1%. For some points with signal to noise ratio < 0 dB, the measurement precision is better than 0.17%. This demonstrate that our network can capture the tiny HRT precursor waves^{[14, 15]}. The four geo-electric stations have been running continuously since the beginning of 2004. Baseline behaviour of the data is well understood. At each station, we took a special measure to eliminate the influence of the surface layer. To remove shallow ground-water seasonal variations^{[16, 17]}, the detection depth is ensured to reach 300 m with current electrodes separation of more than 1000 m. To reduce the polarization instability, the electrodes are made of stable lead plank. At each station, there are three completely independent resistivity monitoring systems (instrumentation, electrodes, wires, etc.) and in each direction, there are two synchronous channels to measure the electric resistivity (ρ_x) and telluric potential (ΔV). The electric resistivity and telluric potential are measured every live minutes and variance of the data is given every hour. Data are transferred to a control center in Beijing through a public network in real time. The time base in each station is synchronized. The standard deviation of the hourly geo-resistivity measurements at the four stations are less than 0.1%.

In the focal zone, the earth material is porous with fluids. The earthquake forming process will create disturbances to the fluids in the pore. The characteristics of the fluid low are an important factor to control the mechanical behaviours of the Earth medium. Small disturbances of the acoustic pressure superimposed on the balanced value of the pore pressure plays an important role in the earthquake forming processand precursors^[18~21]. The geo-electric properties of rocks (geo-resistivity and geo-current) are mostly related to pore fluids. Because of this relationship, scientists in many countries (Japan, Russia, USA, Greece, etc.) have paid great attentions to the geo-electric (resistivity and geo-current) precursor studies. In recent years, much of the attention has been paid to geo-current studies^[22~31].

For earthquake prediction, it is necessary to establish a model. In addition, all parameters involved changes with time and its dynamics need to be characterized. The difficulty of earthquake prediction comes from that the earthquake forming process and precursor mechanism so far are not understood by people. Motivated by the observations before the 1976 M_s7.8 Tangshan earthquake and in consideration that the tidal force is well understood, it appears everywhere and its period is far larger than the natural period of the fault system, we choose the tidal force as the "input signal" and geo-resistivity and telluric current as "output signal" due to its sensitivity to pore fluid disturbance. By using of process recognition^[32] technology, we are able to establish the HRT wave model for earthquake prediction and quantitatively estimate the mathematical model and its parameter variations.

In acoustics^[33], we shall make a start by looking at the simplest case, that of an idealized fluid, uniform and continuous in its properties, at rest in thermodynamic equilibrium, except for the motion caused by the sound waves themselves, and with this acoustic mode motion small enough in magnitude so that many nonlinear effects are negligible and the equation of motion takes on its simplest form. In the presence of a sound wave the pressure becomes P + p. Because a acoustic pressure (p) gradient produces an acceleration of the fluid, we have the wave equation^[33]:

(1)

The equation can be brought into this form by defining an effective fluid density ρ, compressibilit κ, and wave velocity C, for the fluid in the pores, which is the wave equation of porous media for simple-harmonic waves. We start by considering the propagation of a plane wave in a uniform porous media of infinite extent, at rest in the absence of sound. We have already known^[33] that a general solution of the wave equation (1), traveling in the positive x direction's given by the function p(x-C_t).This represents a wave of constant form p(z), traveling with constant speed C. And then we consider that the effect of the change (p) of confining pressure (P), which is coming from the source on the pore pressure change (p_p) in medium under a station. We have already known^[34] that discussion of the general situation is quite complicated and only two kinds of simply special cases are discussed here: one is undraind case, in which p_p is proportional to p; and the other is draind case, in which p_p keeps no changing, unrelated to p.

In the earthquake forming process, the coupling of the solid rock matrix and its pore fluid will modify the wave velocity of both parts of the porous crust. Here, as an example, let us have a look at the plane-wave solution and show that there are two possible waves with two veloci-ties in the porous rocks^[33]. This particular example shows that the wave with velocity near C_f has most of its energy in the fluid motion (the HRT wave); that with velocity near C_s carries most of its energy in motion of the solid network (seismic wave).

Considering viscosity of the pore fluids, pore fluid motion should be expressed with the full expression of each component of the viscosity-stress tensor. From the Stokes-Navier equation of motion of a viscous, compressible fluid, it can be shown^[33] that the lluid velocity u can always be uniquely separated into a longitudinal part u₁, for which curl u₁=0, related to p and a transverse part u_t, for which divu_t=0, unrelated to p and therefore we know that the RT wave always has longitudinal and transverse waves simultaneously.

On the whole the geo-electric properties of rocks (geo-resistivity and geo-current) can directly or indirectly reflect this displacement (strain) or velocity (strain rate) and can be very sensitive to small strain or its rate. Therefore the electric parameters observed at stations can be very sensitive to small stress-strain waves.

To circumvent the difficulty of the earthquake prediction, solving the prediction practical problem and finding the operational way of earthquake prediction, Biru Zhao suggested treating the earthquake process as a "black box". Using the system recognition technology^[32] in the modern control theory, we can calibrate this "black box" mathematically with different input and output parameters. Similar to the "probing method" in industrial engineering such as in order to determine whether there exists a crack, people listens to the sound from that a hammer beats a wheel & axle of train; or from that a bowl is touched by another bowl. In order to determine whether there exists an earthquake source, we analyze the HRT wave recorded before an earthquake and consider earthquake prediction as a systematizing engineering and can establish a mathematical model for this "black box". This will enable us to better understand the earthquake forming and precursor mechanism.

From the rock mechanic experiments^[35], right before an earthquake, the stiffness of the epicentre medium will approach zero or even goes to negative (F is force, u is displacement, σ is stress and ε is strain). At present, people are unable to measure the complete constitutional relationship of the earth medium (rock matrix and its pore fluids) in-situ. In our model of an impending-earthquake precursor of geoelectricity driven by tidal force, we consider that a system of two geological masses suffers under quasi-steady initial loading F₀ and on which a sine wave of a small tidal forces is superimposed^[9]. We obtain the increment of displacement and its geo-electric response is , the output signal of the black box is inversely proportional to media stiffness λ in the epicentre area medium (as shown in Fig. 2). Appliance of some theories of acoustics^[33] to the HRT precursor waves, we obtain the velocity of the medium in the epicentre area is inversely proportional to the mechanical resistance R_m.

图 2

Fig. 2

Fig. 2 Schematic diagram showing the model (mechanism) of HRT wave (a, b, c) Schematic diagram showing the model (mechanism) of HRT wave; (d) The measured HRT wave cases prior to the Ludian earthquake on 10th of August, 2004 (cf. No.16 in Table 2), which can be compared with the model.

We can see that the quasi-steady initial loading F₀(in our modal) is analogous to the equilibrium pressure P (in acoustics) and the amplitude of tidal force is analogous to acoustic pressure (p). When λ continues to decrease and as the period T₀ (the natural period of the fault system) equals to the period of the input signal.

From Fig. 2c we can see that physical indicators to short-term and/or impending stages of an earthquake forming process and we can see that in the early stage (medium term stage) of an earthquake forming process, the epicentre medium is at an elastic stage and rock stiffness λ is a large, positive number. In this case, the geo-electric response due to tidal forces is very small (between 0.05 % and 0.1%) as the curve E-t shows. When the stiffness approaches zero (short-term stage), the displacement of the medium (mostly the pore fluids and electric response) in the epicentre region follows the variation of that of tidal force (the period equal to that of the tidal force) while the amplitude is inversely proportional to the stiffness λ (Although Δu is proportional to A, but normally A changes little. Therefore the electric response is inversely proportional to λ.) and increases dramatically, i.e., when , the harmonic oscillation wave (HT wave) happens, which follows an earthquake in the period from two months to several days^[9] and what which is the physical indicator of short-term stages, as shown by curve D'-t₁; When λ continues to decrease and as the period T₀(the natural period of the fault system) equals to the period of the input signal (impending stage), i.e., when , we are in the condition of resonance and the resonance oscillation wave (RT) appears in the epicentre region and what which is the physical indicator of impending stages, which follows an earthquake in the period, normally, from several days to one or two hours. The RT waves are emitted in pairs with opposite polarity. The velocity of the medium in the epicentre region (and its electric response) is inversely proportional to impedance. When T=T₀ (natural period of the fault), the imaginary part of the impedance is zero. The velocity of the medium in the epicentre area is inversely proportional to the mechanical resistance R_m. The mechanical resistance R_m of the pore fluid is normally very small, and thereby the electric response increases suddenly and is at its maximum^{[14, 15, 33]}. Two RT waves which travel with different velocities (V_fast, V_slow), propagate outward from the epicentre at the same time. In the same station, the difference of arrival time of two RT waves is proportional to the distance from the station to the epicentre. The length of time required for the amplitude of the motion to reduce to (1/e) which occurs in one cycle, of its original value is the natural period T₀ which is related to the scale of the fault system (or earthquake magnitude), as shown by curve D'''-t₃.

From the HRT wave model, we can see that the negative stiffness value (i.e. the descending slope of a constitutional relationship curve) is a necessary condition for the instability (earthquake) of the fault, i.e., the medium, mostly due to the pore fluid, in the epicentre area has a stiffness of λ≤0. In this stage, the medium in the epicentre area will produce harmonic oscillation (the period equal to that of the tidal force) HT wave. Right before the earthquake when the period of the driving force equals to the natural period T₀, the resonance oscillation RT wave will be emitted resonantly. The weak, but clear HRT wave is propagating from the epicentre and being recorded by our anti-interference geo-electric monitoring instrument (PS-100HRT) network. When a network of monitoring stations received such waves, we can predict the upcoming earthquake from the spacial and temporal variations of the HRT waves.

A network of PS100 stations in Sichuan and Yunnan, China has recorded the HRT precursor wave to the M_w9.0 Sumatra tsunami earthquake on December 26^th2004. This is the first record of precursor wave. Thereafter, we have captured HRT waves from more than twenty strong earthquakes, which are well-matched and show repeatability, consistency and regularity.

Fig. 3a shows that three month before the Sumatra earthquake, there are half-month period HT waves recorded at MN station with increasing amplitude (maximum of 2.4%). This is 20 times of the standard deviation of the measurement. Half a month after the earthquake, the HT wave disappeared and there is no appearance of the HT wave from January 2005 to May 2005. Therefore the half-month period HT waves prior to three month before the earthquake are short-term precursor wave to the Sumatra earthquake (cf. No.1 in Table 2).

图 3

Fig. 3

Fig. 3 HT waves to Sumatra Mw9.0, Wenchuan Ms8.0 and Sumatra Ms8.7 earthquakes (a) HT wave prior to the Sumatra Mw9.0 earthquake [3.15°N, 91.79°E; 08:58:51.1(Beijing time)]on 26th of December, 2004; (b) HT wave prior to the Sumatra Mw9.0 earthquake recorded from 2004-12-01 to 2005-01-06; (c) HT wave prior to the Sumatra Ms8.7 earthquake [2.07°N, 97.01°E; 00:09:36] on 29th of March, 2005; (d) HT wave 7 and half day variations prior to the Wenchuan Ms8.0 earthquake [31.0°N, 103.4°N; 2008-05-12 14:28:04]; (e) HT half-month variations prior to the Wenchuan Ms8.0 earthquake; (f) Before the Wenchuan Ms8.0 earthquake, HT wave in NW channel recorded every 5 minutes, corresponding to stiffness λ=0 in the epicentre region.

Comparing Fig. 3b with Fig. 3c, we can see that there are similar patterns of half-month period HT waves, but different amplitudes and that of the wave prior to the Sumatra M_s8.7 earthquake (cf. No.2 in Table 2) is smaller.

From Figs. 3d and 3e, we can see that before the Wenchuan M_s8.0 earthquake, the NW channel of HG station shows the HT half-month and three 7-day period variations with increasing amplitude (cf. No. 10 in Table 2). About 10 hours before the earthquake, the geo-resistivity decreased for about 1.0%. This is 10 times of the standard deviation of the measurement. From Fig. 3f we can see two half-day period variations indicating that the medium stiffness λ≤0. That is, from 4 May 2008 to 6 May 2008, there are two half-day period variations indicating that the stiffness of medium in the epicentre area changes from "+" to "-" and the resonance oscillation RT wave and the earthquake is coming up.

Fig. 4a is the hourly measurements of geocurrent at MN station. Fig. 4b is the 5-minute measurements of geo-resistivity at MN station (Δ=2903 km, A_z=-12.8°). Fig. 4c is the 5-minute measurements of geo-resistivity at HG station (epicentral distance to Sumatra earthquake Δ=2686 km, azimuth A_z=-13.5°) during December 24 to 29, 2004. The two stations were in the first P-wave arrival "+" zone of the M_w9.0 Sumatra earthquake. According to the results of in situ experiments, geo-resistivity decreases when com-pressing^[36~38]. There were two obvious minima of geo-resistivity after the earthquake atHGand MN stations, a, b and a', b', respectively. Corresponding to a, a', there was a remarkable extreme value of geo-current a'' at MN station. Comparing the time of the minima at the two stations, a propagation phenomenon of significance can be found. (Time service of all the stations was offered by the control center, error from CST within 1 s.) The minima a, b at HG station were earlier than a', b' at MN station, respectively. Accordingly, phase velocity of the RT wave was about 306 km/h. Multiplying the travel times of the minima a, a from the M_w9.0 earthquake with the phase velocity obtained above, we have, source of the electric signal was 2.64 × 10³ km away from HG station, and 2.89 × 10³ km from MN station, close to the real epicentral distances of the two stations, respectively. This indicates that the electric signal was generated at the same location and time as the M_w9.0 Sumatra earthquake, i. e. the first minima of geo-resistivity (a, a') and the first extreme value of geo-current (a'') after the earthquake are the coseismal effects of the Sumatra earthquake. The signal did not arrive at MN station until 9 hours later than the earthquake only because of its low propagation velocity.

图 4

Fig. 4

Fig. 4 Comparison of the RT waves (a) Hourly measurements of geo-current at MN station on December 24th~28th, 2004;(b)5-minite measurements of geo-resistivity at MN station on December 24th~28th, 2004; (c) 5-minite measurements of geo-resistivity at HG station on December 24th~28th, 2004; (d) RT waves 2 days before the Mw9.0 Sumatra earthquake; (e) Comparison of the RT waves prior to 5 earthquakes.

Fig. 4d gives the resonance waves (RT waves) 2 days before the M_w9.0 Sumatra earthquake at MN station. At resonance T=T₀, where the wattles component of the impedance is equal to zero, the amplitude of fluid velocity is inversely proportional to the acoustic impedance R_m (i. e. watt component). That RT wave controlled by acoustic impedance appears in a fault system indicates an earthquake is impendent, e.g. within 1~3 days, according to the experience of a few earthquakes. If the acoustic impedance is small during resonance, which is usual, the peak will be very sudden and sharp, with magnitude as large as 2.9%, 17 times larger than the standard deviation of the observations, resulting in a remarkable anomaly.

In Fig. 4d there are 2 RT waves, an earlier decreasing one and a later increasing one (different from RT waves of the Jiujiang earthquake, increasing then decreasing), with arrival time difference about 14 h at MN station, which is of significance (cf. No.1 in Table 2).

Fig. 4e compares RT waves of several earthquakes. The curves from top to bottom are the records of RT waves before the M_s7.6 Banda Sea earthquake at MN station, the M7.8 Pakistan earthquake of 2005-10-08 at YM station and at MN station²⁾, the M_s5.7 Jiujiang earthquake at HG station and at MN station, and the M5.6 Ludian, Yunnan earthquake of 2004-08-10 at MN station, respectively. The arrival time of the primary RT waves before the earthquakes were aligned. It shows that the arrival time of the second RT wave was a function of epicentral distance of the station. In general, the arrival time of the second RT wave was delayed more greatly when the epicentral distance was larger, i.e. the arrival time difference between the primary and the second RT waves increased with the epicentral distance.

2) Three distant specification tests were conducted by the control center in Beijing on 2005-08-12, 2005-10-04, and 2005-10-10, for the working state of the PS100 system at MN station. The result indicated that during the specified time period (August 12 to October 10, 2005), the PS100 system was stable. Therefore, the observed geo-resistivity anomalies before the Pakistan earthquake are real.

It must be noted that the PS100 stations were in the first P-wave arrival "-" zone of the M_s5.7 Jiujiang earthquake (According to the results of in situ experiments, increasing of geo-resistivity is due to dilatation^[36]), so that the RT waves at the three stations increased at first and then decreased one day before the earthquake; whereas for the other earthquakes the stations were in the "+" zone so that the RT waves decreased before increased 1~3 days before the earthquakes.

Fig. 5a shows the relation of arrival time and epicentral distance from several earthquakes to the four stations. The points representing Sumatra earthquake (YN) to MN station, to HG station, the M_s5.7 Jiujiang, Jiangxi earthquake (JJ) of 2005-11-26 to MN station and points of some other cases nicely fall onto a straight line of velocity of 307 km/h.

图 5

Fig. 5

Fig. 5 Cases of RT wave and its mechanism (a) Velocity curve of geo-electric RT wave (longitudinal) and the propagation velocity related to pore pressure changes calculated from the data of large and small earthquakes; (b) The geo-resistivity measurement at MN and telluric field measured at HG two days before and one day after the Sumatra earthquake; (c) RT longitudinal and transverse waves of geo-resistivity before the Mw9.0 Sumatra earthquake and its ellipse of principal stress; (d) Vector graph of the RT longitudinal and transverse waves of geo-current field during the Mw9.0 Sumatra earthquake.

From the opinion that shallow strong earthquakes can result in significant transient changes of pore pressure, while pore pressure changes can cause some of the moderate and small earthquakes^[21], we plotted the points of the distance difference and the time difference between large and small earthquakes representing the M_w9.0 Sumatra earthquake (2004-12-26, 08:58: 51. 1 CST, 3.15°N, 95.79°E), and the M_s5.0 Chuxiong (CX) earthquake (2004-12-26, 15 : 30:13.9 CST, 24.76°N, 101.75°E), and the M_s7.6 Banda Sea earthquake (2006-01-28, 00:58:50.5 CST, 5.4°S, 128.1°E), and the M_s4.8 Yongsheng (YS) earthquake (2006-01-28, 15:02:01.4 CST, 26.3°N, 100.8°E) on Fig. 5a (cross-triangles), and found that the points of large and small earthquakes were all close to the straight line of velocity 307 km/h obtained from geo-electricity data. This indicates that for both large and small earthquakes, propagation velocity of waves related to pore pressure was consistent with the phase velocity of RT longitudinal wave calculated from geoelectricity data. That is, a type of acoustic waves with velocity much lower than that of seismic waves does exist in the crust, energy propagating mainly in pore fluids in the porous rocks of the crust. It has been proved^[33] that there are two possible types of waves in elastic porous media: one is seismic waves, with velocity close to the wave velocity in the solid framework, and most of the energy propagating during the motion of the framework; the other type of waves is acoustic waves with most of its energy propagating in the motion of pore fluids, and velocity close to the wave velocity in the pore fluids. Coupling of the two types of waves makes the seismic waves in the solid framework becoming faster, and the acoustic waves in the pore fluids becoming slower. The RT wave belongs to the latter type.

Fig. 5b's upper panel shows the geo-resistivity measurement every five minutes at MN and lower panel is telluric field measured every five minutes at HG. From Fig. 5b, we can see that two days before and one day after the Sumatra earthquake, i.e., on 25^th, 26^th and 27^th of December, 2004, both geo-resistivity and telluric field start to show step-like variations (basically one step a day), and MN station shows retard compared to HG station. This is because MN station has a larger epicentre distance Δ=2903 km than the HG station Δ=2686 km, thereby it takes longer time for the HRT wave to propagate to the MN station. This measurement confirms the existence of the HRT wave.

From Fig. 5b, the most interesting phenomenon is that the width of the step becomes longer daily. It is 4.5 hour on 25^th; 9 hour on 26^th and 13.5 hour on 27^th. From this figure, we can clearly see the Doppler effect, i.e., the natural period T₀ becomes longer before the earthquake (25^th and 26^th) due to the fluid movement in the fault (from a pipe model) before earthquake rapture^[33]. We think this is due to the pore fluid movement along the fault during pre-rapture fault.

According to the results of in situ experiments of geo-resistivity changes^[36~38] and the method to calculate the anomaly parameters of geo-resistivity heterogeneity^[39]. Fig. 5c gives the effective resistivity ellipses before the M_W9.0 Sumatra earthquake calculated from the three-component geo-resistivity values at the stations, to describe the directions of the principal stress. We can see, the primary wave was longitudinal wave, and the second wave was transverse wave. Left panel in Fig. 5c shows first motion of P-wave positive (+). a-a' resistivity decrease (compression) and the propagation direction of the RT wave is in parallel to the maximum principal stress direction. This is the compression mode wave propagation. Middle panel shows b-b' resistivity increase (expansion) and the propagation direction of the RT wave is perpendicular to the maximum principal stress direction. This is the shear mode wave propagation. Right panel shows the three component georesistivity values at the MN station. The distance between the MN station and Sumatra M_w9.0 earthquake epicentre is 2903 km. The resistivity decrease (a-a') happens 48 hours before the earthquake, which leads b-b' resistivity increase by 14 hours.

Fig. 5d is the vector graph of the geo-current field, according to the three components of the geo-current field at MN station during the M_w9.0 Sumatra earthquake. We can see, as above, the primary wave was longitudinal wave, and the second wave was transverse wave. Relative to the occurring time of the M_w9.0 Sumatra earthquake, using the arrival times of the primary and the second RT waves of geo-current field, propagation velocities of the longitudinal and the transverse RT waves were determined: V₁=346 km/h, V_t=126 km/h, respectively. Taking these values into the equation of virtual wave velocity V_v=(ViVt)/ (V₁V_t)/(V₁-V_t), we have, virtual wave velocity during the earth-quake was V_v=198 km/h, which is almost consistent with the RT virtual wave velocity before the earthquake 207 km/h (Fig. 6e). The arrival time difference of longitudinal and transverse geo-current fields during the earthquake was consistent with the arrival time difference of longitudinal and transverse waves of geo-resistivity before the earthquake. Therefore, it can be concluded that the RT waves (longitudinal, transverse, and virtual waves) appearing before the earthquake were emitted by the early nonseismic failure of the faults before the earthquake in the principal hypocenter area, because the arrival time difference of RT longitudinal and transverse waves and virtual wave velocity during the earthquake were consistent with the RT waves before the earthquake, and the RT waves during the earthquake were emitted by the fault break of the earthquake. That is, at least for the shallow large earthquakes investigated, modulated by tidal forces, non-seismic early fault failures had happened in the hypocenter areas 1~3 days before the earthquakes. Left panel in Fig. 5d shows the a-a' electric field is in the NE direction, parallel to the RT wave propagation (compression mode propagation). Middle panel shows the c-c' electric field is in the ES direction, perpendicular to the RT wave propagation (shear mode propagation). Right panel shows the three components of the geo-current field at MN station during the M_w9.0 Sumatra earthquake. The a-a' electric field at MN rises 8.4 hours after the earthquake. V₁=346 km/h. The c-c' electric field rises 23 hours after the earthquake. V_t=126 km/h. V_v-V₁V_t/(V₁-V_t)=198 km/h.

图 6

Fig. 6

Fig. 6 Cases of RT wave and its propagation (a) Case of RT wave to the two-Taiwan Ms7.0 earthquakes on December 26th, 2006;(b) Case of RT wave to the Kuril Islands Ms8.1 earthquake on November 15th, 2006;(c) Case of RT wave to the Banda Sea Mw 7.6 earthquake on January 28th, 2006;(d) The curve of earthquake magnitude (M) versus the natural period of the fault system (T0); (e) Velocity of RT virtual waves before the earthquakes.

Fig. 6a shows example from the two-Taiwan M_s7.0 earthquakes on December 26^th, 2006 and from which, we see that there are no changes 6~8 days before the earthquake. But since five days before the earthquake, we see RT wave with an increasing amplitude. The time difference for this RT wave is 10 hours, which gives an epicentre distance of about 2070 km. This agrees with the actual epicentre distance of 2019 km (cf. No.6 in Table 2). The amplitude increases to 1.5%, which is significant compared to the standard deviation of the measurements.

Fig. 6b shows example from the Kuril Islands M_s8.1 earthquake on November 15^th, 2006. The epicentre distance estimated agrees with the earthquake tact (cf. No.5 in Table 2).

Example of HRT wave from the Banda sea M_w7.6 earthquake on January 28^th, 2006 is given in Fig. 6c and which shows that MN recorded HRT wave 6 days before the earthquake with increasing amplitude daily but no change from 10 to 7 days before the earthquake. One day before the earthquake, the RT wave with amplitude of 2.3% appeared. The time difference between first sudden drop of RT wave and sudden rise of RT wave is 23 hours. From this, we estimated the epicentre distance is 4680 km. This agrees with the actual epicentre distance of 4683 km (cf. No.4 in Table 2).

The decay rate of the RT wave gives the natural period of the fault (see Fig. 4d). The earthquake magnitude can be predicted by M=4.14logT₀ -1.486 (Fig. 6d), where T₀ is measured in minutes.

Plotting arrival time difference of the two RT waves before a earthquake at each station Δt=t_t -t₁ with respect to the epicentral distance from the station to the coming earthquake Δ into Fig. 6e, we can see all the points are close to a straight line with virtual wave velocity V_v=(V₁V_t)/(V₁-V_t) ≈ 207 km/h, indicating that the virtual wave velocity were the same for those earthquakes. Hence, as long as arrival time difference of the two RT waves at one station is known, it is possible to calculate the epicentral distance of the coming earthquake quantitatively. With a network of appropriately distributed PS100 stations, epicenter of an earthquake about to occur can be determined 1~3 days in advance.

Therefore, it can be concluded that the RT waves (longitudinal, transverse, and virtual waves) appearing before the earthquake were emitted by the early non-seismic failure of the faults before the earthquake in the principal hypocenter area, because the arrival time difference of RT longitudinal and transverse waves and virtual wave velocity during the earthquake were consistent with the RT waves before the earthquake, and the RT waves during the earthquake were emitted by the fault break of the earthquake. That is, at least for the shallow large earthquakes investigated, modulated by tidal forces, non-seismic early fault failures had happened in the hypocenter areas 1~3 days before the earthquakes.

Fig. 7a shows example from Indonesia M_s8.5 earthquake on September 12^th, 2007 and which shows that five days before the earthquake, there are RT waves (shaded area) recorded at MN station. From this, we estimated the magnitude and epicentre distance are in agreement with the actual earthquake parameters (cf. No.8 in Table 2).

图 7

Fig. 7

Fig. 7 RT wave cases and its regularity (a) Example from Indonesia Ms8.5 earthquake on September 12th, 2007;(b) NW channel recordings of HRT waves recorded at HG before Wenchuan Ms8.0 earthquake; (c) NW channel recordings removed from long term changes; (d) Example from the M5.7 Jiujiang earthquake of 2005-11-26 at HG, MN and LJ stations; (e) Example from the M7.8 Pakistan earthquake of 2005-10-08 at YM, HG and MN stations and three distant specification tests.

From Figs. 7b and 7c, we can see he half-day HRT variation with amplitudes several times the standard deviation of the measurements. Several days before the Wenchuan earthquake, the amplitude of the half-day variation starts to increase. Superimposed on top of the half-day variation, there are RT wave with the period of 3~4 hours. We can use these wavelets to obtain the magnitude and epicentre distance information. One day before the earthquake, RT wave reverse its trend indicating that the epicentre medium stiffness λ≤0.About 1 0 hours before the earth-quake, the geo-resistivity suddenly decreases for more than 1%, and there are 3 hour period variations present. This demonstrates the rapture length (earthquake magnitude) and epicentre distance (see Table 1 and No.10 in Table 2), which forms the strong earthquake warning signal. From Fig. 7c, we can see the NW channel with increasing amplitude of the RT variation (May 9^th, 0.3%; May 10^th, 0.5%; May 12^th, 1.0%). Prediction results for epicentre distance and magnitude are shown in Table 1.

表 1(Table 1)

Table 1 Earthquake prediction results from the RT wave in EW channel for epicentre distance and magnitude

Table 1 Earthquake prediction results from the RT wave in EW channel for epicentre distance and magnitude

表 2(Table 2)

Table 2 Cases of HRT wave

Table 2 Cases of HRT wave

The epicentre distance estimated at different times varies between 500 and 800 km. Because the epicentre is the north of HG station, the center of the pore fluid wave might move north-south. The center on May 12^th was 555 km and 810 km, respectively. If we consider this difference to be the length of the rapture L=810 -555=255 km, migrating from south to north and we assume that this represents the two ends of the fault, we can estimate the rapture length. The magnitude can then be estimated by^[40]: M=3.3 + 2.1logL=8.5. This is in good agreement both with the estimation shown in Table 1 and with the result of Zhang et al.^[41].

The four stations in Sichuan and Yunnan are generally in a straight line of NS direction. The M_s5.7 Jiujiang earthquake was to the east of those stations; the Pakistan earthquake was to the west of the stations. The stations were generally on the same P-wave front of the two earthquakes, so the in-phase anomalies of the two earthquakes at three stations were shown in Figs. 7d and 7e respectively.

The M_s5.7 Jiujiang earthquake was to the east of those stations. The epicentre distance estimated agrees with the earthquake fact (cf. No.14 in Table 2).

Fig. 7e shows example from the M7.8 Pakistan earthquake of 2005-10-08 at YM, HG and MN stations and three distant specification tests were conducted by the control center in Beijing on 2005-08-12, 2005-10-04, and 2005-10-10, for the working state of the PS100 system at MN station. The resut indicated that during the specified time period (August 12 to October 10, 2005), the PS100 system was stable. Therefore, the observed the HRT waves before the Pakistan earthquake are real. The magnitude and epicentre distance of prediction are consistent with the facts of the earthquake (see No.3 in Table 2).

Cases of HRT wave and the consistency of prediction with the earthquake fact are shown in Table 2.

The HRT precursor wave has been recorded at four PS-100 stations in South-West China before more than twenty strong earthquakes. The HRT wave short-term and impending earthquake prediction method was developed on the bases of the former geoelectricity earthquake prediction method, which included earth-resistivity, MT, geocurrent field, artificial super low frequency (SLF) waves and seismic electromagnetic signals etc., and has been used in China for more than 40 years^[42~59]. According to our two internal reports to the China Earthquake Administration in 2003^[12], we got an emergency funding from the China Earthquake Administration and deployed four geo-electric monitoring stations in Sichuan and Yunnan provinces from 2004, using the new generation of equipment (PS-100) and technologies to capture the HRT wave earthquake precursor. The first record of HRT precursor wave has been recorded from the 26th/Dec/2004 Sumatra M_w9.0 earthquake with the largest epicentre distance Δ=2900 km. Thereafter, we have captured HRT waves from more than twenty strong earthquakes, which are well-matched and show repeatability, consistency and regularity. Before the Wenchuan M_s8.0 earthquake, we recorded the HRT wave precursor at the only operating station in Hongge (HG, Δ=465 km). We found that significant impending signal has a causal relationship with the Wenchuan earthquake recorded at the HG station in the early morning (0~5AM) of 12th/May/2008. These facts demonstrated that an earthquake, particularly a strong earthquake, has precursors, and it is predictable. Short-term or impending earthquake prediction is feasible in the near future. If we establish the PS-100HRT network of stations to a large scale, we might be able to predict earthquake in advance. The HT wave is generated by pore fluid movement driven by tidal force when the epicentre medium stiffness approa-ches zero before the earthquake. The RT wave is emitted by resonance when the natural period of the fault equals that of the driving tidal force. If we have a global PS100 network of appropriately distributed stations, accurate earthquake prediction will be possible in the near future.

The HRT wave short-term and impending earthquake prediction method was developed on the bases of the former geo-electricity earthquake prediction method, which has been used for more than 40 years. This method was gradually developed with the supervision and effort of many colleagues in China Earthquake Administration. In particular, Gu Gongxu, Qin Xinling, Zeng Rongsheng, Xie Lili, Zhu Chuan-zhen, Luo Zhuoli, Guo Zengjian, Chen Yuntai, Ma Zongjin, ChenZhangii, Wu Zhongliang, Yang Jiansi, Chen Rong, Ma Jin, Xu Shaoxie, Mei Shirong, Che Shi, and many others have shared their research results and ideas with us. This has made a signiticant contribution to the development of the HRT short-term and impending earthquake prediction method. Yan Guliang, Wang Chengmin, Geng Qingguo and Lin Mingzhou offered support and encouragement for the composition of this article. The geoelectricity data are from a program funded by the China Earthquake Administration and its sub-division in Sichuan and Yunnan provinces.Professor Tang Yuxiong is acknowledged for diligently running the observatories. The seismic P-wave data is from the work of Li Baokun, Zhang Liwen and Qiu Haijiang. Thegravity tidal data was from the computation of Liu Duanfa.