Short-term Earthquake Prediction Based on Seismic Precursory Electric  Signals Recorded on Ground Surface.

2. VARIOUS TOPICS IN SEISMOLOGY – TECTONICS PERTAINING TO EQ PREDICTION

Introduction

   There is a general notion, in the seismological community, that earthquakes can occur in a random way in every place in a large area that exhibits seismicity. This is corroborated by the existing theories with regard to the generation mechanisms of earthquakes and the strain charge conditions held in such places. The Self Organized Criticality (SOC) theory, proposed by Bak (1966); Bak; Tag; and Wiesenfeld (1988), indicates that physical systems, like the seismogenic areas, are at critical state or instability and therefore, any small earthquake has some probability to cascade into a large earthquake. Other researchers (Otsuka 1972a, b; Vere-Jones 1976; Bak and Tang 1989; Ito and Matzutaki 1990) concluded that the generation of a strong earthquake depends on very small variations of stored elastic energy, fault strength or variations of the elastic properties of the seismogenic region, therefore, statistically, a strong earthquake can take place anywhere.
   Since all the latter theories are combined to the fact that the Earth, as it is deduced from the geological, tectonic maps, compiled by the geologists, is highly fractured, then it is justified to conclude that strong earthquakes can, in a random way, occur almost anywhere. However, as the following analysis will demonstrate, the localization of strong earthquakes follows a scheme which is prescribed by the deep tectonics and the intense fracturing of the lithosphere, as far as it concerns the shallow earthquakes which occur in the crust.


2.1 Spatial distribution of strong EQs.

   Before any attempt is to be made for any earthquake prediction, it is required to be known at least the general tectonic setting for the area under study. This is justified by the fact that EQs are the results of intense, tectonic processes which take place in the lithosphere. In this work a different approach, in relation to EQs, was adopted. Applied geophysical (gravitational) methods have been used in an attempt to explain the spatial distribution of strong EQs and to prove that these occur only on predefined narrow faulting - fracture zones of the lithosphere, detectable in a deterministic scheme. The theory of lithospheric plates and their motion is well known and has been adopted by the seismological community. A generalized picture of the Earth's lithospheric plates is presented in figure (2.1.1).

Fig. 2.1.1. The Earth's lithospheric plates (after USGS).
 

   The main lithospheric Earth's plates are subdivided in smaller ones. McKenzie (1972, 1978) presented a small-scale lithospheric model for the particular case of the Greek territory. This is presented in the following figure (2.1.2).

Fig. 2.1.2. McKenzie's (1972, 1978) lithospheric, plate model, proposed for the Greek territory.
 

   As these plates are in continuous motion, due to applied forces from the surrounding larger plates, extension or compression is applied on them as a result. Therefore, extensional or compressional fracture zones are formed. The kinematics of the Greek territory tectonic plates was studied by Papazachos et al. (1996) and is presented in the following figure (2.1.3). It is understood that strong seismic events are expected to occur where large tectonic activity takes place.

 
     

Fig. 2.1.3. Kinematics of the Greek territory tectonic plates proposed, by Papazachos  et al. (1996).
 

   The main features of this map are the areas which are characterized by compressional, extensional and strike-slip motion. As it will be explained later on, this plays an important role in the kinematics of the Aegean area and moreover facilitates the calculation of its angular, rotational velocity.
 

2.1.1 Mapping of major, seismic fracture zones – faults.

   Mapping of fracture zones - faults is of great importance in seismological studies. These are the places where earlier EQs took place and most probably future EQs will occur. Therefore, the accurate mapping of these fracture zones - faults is of highest priority. The problem that arises is: what kind of faults - fracture zones must be mapped and is it always possible?
   Earlier studies have revealed rather broad, seismic zones of increased, seismic risk, based, on the Greek seismicity history. In next figure (2.1.1.1) De Bremaecker et al. (1982) divided the regional, Greek area into smaller blocks which generally move towards southwest with different velocities in respect to the plate of Eurasia.

Fig. 2.1.1.1. Finite element network with boundary condition velocities, regarding European plate as proposed by De Bremaecker et al. (1982) for the Greek regional area.

   In another study, Papazachos (1988) defined the most active seismic areas in Greece, based on the spatial distribution of the shallow earthquakes. The latter is shown in next figure (2.1.1.2).

Fig. 2.1.1.2. Seismic zoning of shallow earthquakes in the Greek and surrounding areas (Papazachos, 1988).


   Furthermore, the Greek territory was divided in specific zones according to the seismicity as it is expressed by the yearly released, seismic moment per 10000 Km2 (Papazachos, 1989). This is presented in the figure (2.1.1.3).

Fig. 2.1.1.3. Seismicity of the seismic zones of the Greek territory as it is expressed by the yearly, released, seismic moment per 10000 Km2 (Papazachos, 1989).


   In a more recent study, Hanus and Vanek (1993) presented another mapping of the active, seismic zones in the Greek territory which is presented in the following figure (2.1.1.4).

Fig. 2.1.1.4 Mapping of the active, seismic zones in the Greek territory as presented by Hanus and Vanek, (1993).


   Comparing the latter three figures (2.1.1.2, 3, 4) it is shown that the mapping of the active, seismic zones in the Greek territory (concerning shallow seismicity) is very subjective to each researcher. The characterization of an area, as an active, seismic one, is based on subjective criteria of each scientist, who presents his work on this topic. The only common, broad feature between them is the absence of active, seismic zones in the central Aegean area (Cyclades). An objective, active, seismic zoning of the Greek territory can be obtained by taking into account the deep lithospheric fracture zones / faults.
   As far as it concerns the detailed mapping of seismic fracture zones and faults, classical geological surface observations cannot detect all of them. This is especially true when no surface trace has been produced by any deep, subsurface, tectonic activity. Furthermore, it is difficult to assign a seismological significance, to a detected on surface fault, without any other deep tectonic knowledge. The previous problem has been overcome, partly, by studying the delineation of the seismic events over a certain area. An example of this procedure is presented in the following figure (2.1.1.5).

Fig. 2.1.1.5. Delineation of a deep fracture zone – fault from its seismic activity, mapped  over  a certain period of time (June, 2001). The seismic activity for the assigned, spatial window (dashed square) was mapped in June 2001. It is obvious the delineation of the seismically activated fault. As a result of this activity a strong seismic event occurred at the lower left part of this fault (Psara EQ, Ms = 5.6R, June 10th, 2001).

2.1.2. Mapping of major, seismic fracture zones - faults from the study

         of the Earth's gravity field.

   The creation of fracture zones and faulting of the Earth is due to the present stress field conditions which are applied in a particular area. Fracturing takes place in various modes (Mattauer, 1973), when extensional or compressional forces are applied on it. In the following figure (2.1.2.1) the black arrows indicate the applied forces, while the lines in the solid block indicate the various generated modes of fracturing. In case (A) there is block movement which creates a typical thrust fault combined with internal micro-fracturing typically normal and parallel to the direction of the applied stress field. In case (B) the generated, internal, micro-fracturing is interconnected and delineated diagonally in respect to the applied stress field direction, while in case (C), a characteristic strike-slip fault has been formed, due to exceptional type (B) stress field conditions.

Fig. 2.1.2.1. Fracturing – faulting which is generated by applied stress (Mattauer, 1973). The arrows indicate the applied stress direction (compressional, extensional, shear-stress).

   The form of a typical fault is presented in the following figure (2.1.2.2), regardless of its cause of origin. Typically, the ground is split into two parts through the fault plane, which exhibit a differential mode of relative movement to each other. The slip-vector of the generated fault, characterizes the latter.

Fig. 2.1.2.2. Generalized sketch of a down-throw fault. The thick vector represents the slip vector; the thin vectors indicate the relative movement of each ground block to the other one. The down-throw block has slipped along the fault plane.

   The slip-vector, that characterizes a fault, can be generally considered as the combination of its orthogonal components, when it is projected to an orthogonal coordinating system X, Y, Z. This is presented in the following figure (2.1.2.3).
   In each case, the interrelation of the magnitudes between these components can generate the various types of faulting which are observed in nature. Moreover, by accepting that the fault plane coincides with the Z-X plane, then the SV-Y component equals to zero. In such a case, (when SV-Z = 0), no thrust, neither normal faulting can exist. Since the slip-vector lays in the fault plane, then the only possible, ground block movement is along the strike of the fault. The latter is the case of a strike-slip fault and is characterized as sinister or dextral one depending on the direction of the block movement. In case SV-X = 0, then normal faulting or thrust faulting is generated, depending on the sign of the SV-Z component. In nature, faulting that exhibits both types of faulting modes, is, typically, generated. The type of generated faulting is characterized by the dominant and much larger component of the fault slip-vector.

Fig. 2.1.2.3. Analysis of a vector into three orthogonal components. The fault slip-vector is analyzed into its orthogonal components SV-X, SV-Y and SV-Z.

   In a most general case, the slip-vector is analyzed into orthogonal components, along the X-Y horizontal plane and the vertical to the ground Z-axis. In such a case the vertical component VS-Z plays a very important role in the overall tectonic and the stratigraphic conditions which characterize an area. A relative, vertical movement of the two ground blocks of the fault has, as a consequence, the stratigraphic modification of the overlaying geological formations. The latter has as an effect, the generation of ground lateral density discontinuities, and therefore, are generated changes in the intensity of the Earth’s gravity field. This is the physical basis upon which the gravity methodologies for the study of the geology tectonics and stratigraphy of an area were founded.
   Large tectonic fault / fracturing zones of the Earth are mainly reflected on the morphology of its gravity field. The wavelength of these features depends on the dimensions of the corresponding fault / fracture zones. Therefore, by studying the large wavelengths of the gravity field, we can obtain information for the deep large faulting / fractures zones of the Earth.
   Moreover, by studying the horizontal gradient of the corresponding gravity field, we can locate the precise location of a fault / fracture zone quite easily. In the following figure (2.1.2.4), (a) represents a geological fault of the Earth, (b) is the corresponding gravity field, while in (c) is presented the corresponding horizontal gradient of it. Details on this topic can be found in any geophysical gravity textbook.

Fig. 2.1.2.4. (a) = Geological fault model, (b) = Corresponding gravity field, (c) = Corresponding gravity horizontal gradient.

   Generally, the gravity method of Applied Geophysics is capable of tracing lateral density discontinuities of the underground, geological formations.
 

2.1.3 Application of the methodology in Greece.

   On 13th May, 1995 a devastating earthquake (M=6.1R, Lat=400.18, Lon=210.71) occurred in "Grevena" area, northern Greece and caused many damages at the nearby towns and villages (some deaths were reported, too).
   What is interesting, in this earthquake, is the fact of total absence of any statistical seismo-logical indication that a strong EQ could take place, in an area like this, which was characterized as of very low seismic risk (zone 1, out of 4) and of total absence of any known, large, geological tectonic features to justify any possible strong seismic event.
   The earthquake which occurred, the absence of statistical and geological - tectonic data which could reveal the seismically dangerous character of the "Grevena" area and the required, deterministic cause of this earthquake motivated the analysis to follow.
   The very same scenario was replicated on 7th September 1999 in Athens, the capital city of Greece. Although Athens is located in an area which is characterized as "zone 2", as far as it concerns its seismic risk, the complete absence of knowledge of deep tectonics had postulated the notion that the capital of Greece is built over "safe" ground. Unfortunately, for the people (over a hundred) who died, apart from the very large damages in buildings, during this earthquake, this was proved to be totally wrong, after this earthquake.

   So the question that arises is: Where do the strong earthquakes really occur?

   The following analysis may be the answer to the above question.

   It is well known that strong earthquakes occur when a fracture zone / fault releases its stress load and rupture occurs along a large part of it. The longest the fault / fracture zone has been activated, the largest the magnitude of the earthquake is.
   The problem therefore, is formulated as follows: Is it possible to map the deep fracture zones, where strong EQs occur? The study of strong earthquakes in Greece, statistically treated, may vaguely reveal the large, deep, tectonic features that are prerequisites for its occurrence. The same problem, of the detection of deep tectonics of an area, faced from the view of Applied Geophysics, is rather simple, as it will be explained later on.
   The key point of this study is that strong earthquakes must coincide in location with large tectonic faulting / fracturing systems of the lithosphere. Consequently, the specific target of this analysis is the location of these tectonic features, which locations can't be traced by other geological methods.
   The following map (fig. 2.1.3.1), where is presented the location of the strong, (M>=6R), seismic events of the Greek territory, which occurred during the period 1900 – 1997, has been compiled, in order to have a clear view on this topic. The circle, around each location, is the area defined by the (conservative) location estimation error (+/- 25Km). It must be pointed out that the location of each earthquake can be anywhere in the error circle, without decreasing the validity of the map.

Fig. 2.1.3.1. Earthquakes with Ms>=6R which occurred in Greece during the period 1901 - 1997. The red circle indicates the error, accepted, for the location determination that equals to +/- 25Km. The map correlates in size with the Bouguer Anomaly map of Greece (Makris, 1981).

   The gravity field transformation approach into its gradient can be applied in two-dimensional form on the gravity field of Greece. A simplified presentation of the Greek gravity field, Bouguer anomaly (Makris, 1981) is presented in the following figure (2.1.3.2).

Fig. 2.1.3.2. Simplified (Thanassoulas, 1998) gravity (Bouguer anomaly) map of Greece (Makris, 1981). The scale is in mgals.

   The following scheme has been applied in order to utilize the horizontal gradient transforma-tion.

                                                     Grad_x,y(G) = G(x,y)*TL(x,y)                                                    (2.1.3.1)


   Where: G(x,y) is the gravity field at point x, y

                 TL(x,y) is the operator, applied, for the transformation to be utilized

                 Grad_x,y(G) is the resulted, horizontal gradient

   TL has been calculated analytically by the use of a polynomial surface of second order, fitted, to a progressively sliding window over the gravity field data.
   In order to avoid near surface tectonic features, and targeting to the deeper structure of the Earth, a 2-D "window" of 20 x 20Km was used in all these transformations. The results of this operation are presented as photo relief map, illuminated, from NE direction. In the next figure (2.1.3.3) is shown the result of the transformation of the Bouguer anomaly map of Greece into a gradient one.

Fig. 2.1.3.3. Horizontal gradient gravity map of Greece (Thanassoulas, 1998).

   The latter map is presented in 3-D photo relief (fig. 2.1.3.4), so that the high gradient value areas are better visualized.
   Deep fracture zones express (horizontal gravity gradient ridges) themselves as the bounda-ries between dark black and white zones.

Fig. 2.1.3.4. 3-D photo relief compiled gradient map of Greece. Each scale unit is 10Km. Lighting from NE (Thanassoulas, 1988).

   It is very interesting to note the deep fracture zone, crossing in a NW - SE direction the Attiki area (X=32, Y=34), where Athens is located and the strong EQ of Athens took place. This feature, unfortunately, reveals the dangerous and seismically risky position, where the capital city of Greece is located.
   The detailed fracture zones-faults (thick black lines) locations, determined by the latter method, are presented in the following figure (2.1.3.5).

Fig. 2.1.3.5. Deep, lithospheric fracture zones (thick black lines) determined from the horizontal gradient map of Greece (Thanassoulas, 1998). The solid dots indicate the locations of strong (Ms>6R) EQs, which occurred during the period 1901 – 1997.

2.1.4 Strong EQs location map of Greece (1901- 1997).

   It was stated earlier that strong earthquakes have their origin in deep and large faults / fracture zones of the Earth. In the following figure (2.1.4.1), is made a comparison between the location of the strong earthquakes, determined by seismological methods, and the location of the deep faults / fracture zones determined by the transformation of the gravity field into its horizontal gradient. The EQs catalog, which is used, is the one provided by the National Observatory of Athens (NOA).
   Therefore, it is expected that at least a fault zone crosses the error location circle of each strong earthquake, which is mapped. This fault zone is, presumably, the actual fault zone that was activated and generated the corresponding earthquake. Consequently, the "most probable true" location, of the corresponding earthquake, is not the one suggested by seismological methods, but the one which resulted by the projection of the "seismological" location of the epicenter, to the nearest fault / fracture zone.

Fig. 2.1.4.1. Location of strong, seismic events (Ms>=6R, solid red dots) in Greece, for the period 1901 – 1997, in relation to the location of the fracture zones - faults (thick black lines), determined by the transformation of the corresponding gravity field.

   A detailed examination of this map indicates that the vast majority of the earthquakes are located in a distance shorter than 25Km from the nearest fracture zone. Therefore, if the location of each earthquake is moved within its location error radius towards the collocated fracture zone, then the validity of the EQ location map does not change, while at the same time there is a very good coincidence of the EQ locations with the location of the fracture zones (fig. 2.1.4.2). A very small number of EQs, mainly in the western part of Greece, do not follow this rule. This could be attributed either to missing existent fracture zones / faults which the transformation of the gravity field didn't succeed to identify or to the wrongly calculated, by the seismological methods, location of the deviating earthquakes. In either case the degree of correlation (statistically) is very high.

Fig. 2.1.4.2. Corrected location, of the presented EQs, according to the calculated fracture zones - faults from gravity map transformation. The EQs data file contains data from 1901 to 1997, when the original study was presented (Thanassoulas, 1997). The EQs, outside the gravity map boundaries, have been excluded.

   A detailed example, showing the validity of the methodology, is presented in the following figure (2.1.4.3).

Fig. 2.1.4.3. Delineation of a deep fracture zone - fault from its mapped, seismic activity over a certain period, compared with the one, calculated, by the gravity field transformation. A (left) = fracture zones, calculated, through gravity field transformation, B (right) = fracture zone, mapped, from its seismicity, which took place during June 2001. Comparison is made between the inset frames.

   It is very interesting to see how the seismically, mapped, fracture zone in fig. (2.1.4.3 - B) coincides with the one which is mapped by the gravity method and is presented in figure (2.1.4.3 - A). Further more the specific fracture zone extends towards NE in the area of Western Turkey. This fracture zone has already produced two strong EQs with Ms >6R in the past.
 

2.1.5. Verification of fault zones by strong EQs which occurred during the period 1998-2006.

   Accepting the validity of the presented methodology, it is applied as a test on the earthquakes which took place during the period 1998-2006 in Greece. These earthquakes must coincide with the corresponding, deep, lithospheric fracture zones / faults.
   In the following figure (2.1.5.1), are presented for validation purposes the strong (Ms>=6R) earthquakes (solid red circles) which took place in Greece, after the completion of this study (Thanassoulas, 1998), that is the period 1998-2006, along with the determined fault / fracture zones.

Fig. 2.1.5.1. Correlation of the location of recent strong (1998 - 2006, Ms>=6R) seismic events (red circles) to the location of the already calculated fracture zones - faults by the transformation of the Greek gravity field.

   It has been made clear that all the latest, strong, seismic events coincide with the location of the fracture zones - faults which have been already calculated.
   Two more examples, related to Lefkada EQ (2003/08/14, Ms = 6.4R) and Kythira EQ (2006/01/08, Ms = 6.9R), follow. Lefkada EQ caused large damages on Lefkada Island and is the most destructive recent one. The above mentioned EQ is presented in the following figure (2.1.5.2).

Fig. 2.1.5.2. Lefkada EQ, red circle (2003/ 08/14, Ms = 6.4R) in relation to the deep, lithospheric fracturing (brown lines), calculated by the gravity field transformation (Thanassoulas, 1998).

   Kythira EQ is the strongest one which occurred in Greece, recently. It was felt almost allover Greece and the neighbor countries, too. Kythira EQ is presented in figure (2.1.5.3).

Fig. 2.1.5.3. Kythira EQ, red circle (2006/01/08, Ms = 6.9R) in relation to the deep, lithospheric fracturing (brown lines), calculated by the gravity field transformation (Thanassoulas, 1998).

   The presented analysis of strong earthquakes’ location which earthquakes occurred in Greece during the last 100 years (more or less) and the detailed examples, too, show clearly that, contrary to seismological notions that strong EQs can occur, statistically, in any place, these happened in specific, deep, lithospheric fracture zones / faults. These can be mapped by the use of applied geophysical methods, and in particular, through the analysis of the Earth’s gravity field that is modified by these strong, tectonic, lithospheric anomalies. The methodology which is presented could be used to prepare more detailed maps of seismic risk for wider areas, even though seismological data of such areas are missing or are of limited extent. This is an advantage, compared to seismological methods which require a statistically large number of earthquakes to take place in order to evaluate the seismic risk. The case of the EQ in Grevena (seismic risk zone I) and the EQ in Athens (seismic risk zone II) are two characteristic examples, which justify the validity of this methodology.
   In conclusion, taking into consideration the latest, increased, seismic activity of the Greek territory, it is justified to expect, strong earthquakes to occur along this faulting / fracture-zoning network in the future. Therefore, the State Authorities must take into consideration the fact that the seismic risk is enhanced along these specific fracturing / faulting zones. The fact that, in some of them we have not experienced a strong earthquake yet, does not justify, at all, the notion that it is an aseismic zone. On the contrary, it must be considered as a seismogenic area, potentially to be activated in the future. What is open, as a question, that is only the timing of the future strong EQ.
 

2.2. Depth distribution of strong earthquakes.

A parameter of an earthquake of great importance is its depth of occurrence. The vast majority of the seismicity of a seismic region takes place in the crust (fig. 2.2.1) of the Earth. The crust, as it is shown in the sketch presentation of the next figure (2.2.1), is very thin compared with the other discrete parts of the Earth’ interior (mantle, outer-inner core).

Fig. 2.2.1. The structure of the Earth (www.physicalgeography.net). The crust, mantle, outer – inner core of the Earth’s interior (not to scale).

   A more detailed presentation of the lithospheric structure is presented in the following figure (2.2.2). The continental and oceanic crust is shown, in relation to a larger geological unit, referred, as “the lithosphere”, which overlays the “plastic asthenosphere” and the “upper mantle”.

Fig. 2.2.2. The structure of the lithosphere (www.physicalgeography.net). The continental – oceanic crust, the lithosphere, the plastic asthenosphere as long as their representative depth extents are shown (not to scale).

   The characteristic depth values that are assigned to each of them have been identified by specific studies upon the change of their physical properties. The analysis of the P velocity of the seismic waves as a function of depth (Mueller and Landisman, 1966) has revealed the existence of a velocity decrease zone at a depth around 10Km, while an abrupt increase is observed at depth of 30Km (Moho discontinuity). The latter is presented in figure (2.2.3).

Fig. 2.2.3. Model of P velocity distribution within the crust, that shows the postulated low velocity layer in the upper crust (Mueller and Landisman, 1966) and its abrupt increase at Moho depth.

   Toksoz et al. (1967) studied the shear wave velocity as a function of depth. The results of this study are presented in the following figure (2.2.4). The main feature of this study, concerning the depth of occurrence of earthquakes, is the shear wave velocity decrease, observed, at a depth of 100Km in tectonic regions.

Fig. 2.2.4. Upper mantle shear velocity models, for oceanic, continental shield and tectonic regions, based on dispersion of surface waves. The profiles are uncertain below 500Km depths, because of insufficient data (Toksoz et al. 1967).

   The majority of the EQs occur in the seismogenic zone of the lithosphere that is where earthquakes mostly take place. The base of the seismogenic zone is the top of the more ductile asthenosphere. A simplified model of an earthquake and its associated rupture surface is shown in figure (2.2.5).

Fig. 2.2.5. A simplified model of an earthquake is shown and its associated rupture surface (S). The dimensions of the rupture are denoted as W and L for its depth and length extent accordingly. S = rupture surface and D = direction of rupture.

   The detailed study of a seismogenic area and particularly the study of the depth distribution of the earthquakes which occur in that area, for a rather large period of time, reveal the depth extent of its corresponding, seismogenic zone.
   This is illustrated in the following figure (2.2.6), as an example that refers to the depth distribution of the foci of EQs, registered, at the Nicoya peninsula, Costa Rica (Avants et al. 2001, Newman et al. 2002).

Fig. 2.2.6. Depth distribution of the foci of EQs registered at Nicoya peninsula, Costa Rica that illustrates along-strike variability in the seismogenic zone seismicity (Avants et al. 2001, Newman et al. 2002).

   The seismogenic zone of the Nicoya Peninsula revealed from the up-dip events extends from 10 – 20Km, while the analysis of the total events (around 2500) indicates that the seismogenic zone starts from almost 10Km and extends up to 50Km.
   The physical parameters of the lithosphere, generally, change along depth and therefore, the mechanisms that generate earthquake precursory phenomena depend upon the depth of occurrence of each EQ. Consequently, earthquakes of different magnitude which occur in different depths, probably will generate different precursory phenomena and of variable intensity.
   In the Greek territory, the depth distribution of earthquakes is analyzed in an introductory form as follows:

     Number of EQs, as a function of occurrence depth, for a fixed magnitude, in steps of .5R

            (5.5 – 7.5 R).

     Number of EQs as a function of magnitude for a fixed depth, in steps of 5 Km (0 – 10Km).

   In both cases the EQ catalog of NOA, which was used, spans from 1901 to 2006. However, it must be pointed out that the first part of this catalog (1901 – 1950) refers to the early period of earthquake parameter estimation and therefore, because of technological lack during this period, a depth of 0Km was assigned to each registered EQ. Therefore, EQs which refer to a depth of 0Km will be ignored.
 

2.2.1. Number of EQs as a function of occurrence depth for a fixed magnitude of

          5.5 R

   Depth distribution of the EQs in Greek territory with magnitude of 5.5 R is shown in the following figure (2.2.1.1).

Fig. 2.2.1.1. Depth distribution of EQs in Greek territory, with magnitude of 5.5R.

   Four distinct peaks are characteristic in this diagram. The first is located at 0Km depth and will be ignored as was mentioned earlier, the second one is located at a depth of 10Km, the third one spans from 30 to 40Km, while a forth one is located at 100Km. Single EQs extend to a depth of 150Km.
 

2.2.2. Number of EQs as a function of occurrence depth for a fixed magnitude of

          6.0 R

   In the next drawing (2.2.2.1) has been used instead a magnitude of 6.0 R.

Fig. 2.2.2.1. Depth distribution of the EQs in Greek territory, with magnitude of 6.0R

   For the magnitude of 6.0 R the depth distribution of the EQs is similar to the 5.5R one, but with lower peak values. The depth zones of 10Km, 30-40Km, and 100Km are characteristically the same but an extra peak is shown, at a depth of 120Km.
 

2.2.3. Number of EQs as a function of occurrence depth for a fixed magnitude of

          6.5 R

   For the next magnitude level of 6.5 R (fig. 2.2.3.1) the No of EQs is less than the previous case, but the very same depth zones are revealed and one more at a depth of 65Km, as well.

Fig. 2.2.3.1. Depth distribution of the EQs in Greek territory, with magnitude of 6.5R

   It must be taken into account that EQs of such a magnitude are rare in the Greek territory and therefore, even a single seismic event of such a magnitude is very important for the iden-tification of these seismic zones.
   The latter is validated by the depth distribution of EQs with magnitude of 7.0 R and 7.5 R, shown in the following figures (2.2.4.1 – 2.2.5.1).
 

2.2.4. Number of EQs as a function of occurrence depth for a fixed magnitude of

         7.0 R

Fig. 2.2.4.1. Depth distribution of the EQs in Greek territory, with magnitude of 7.0R

2.2.5. Number of EQs as a function of occurrence depth for a fixed magnitude of

          7.5 R.

Fig. 2.2.5.1. Depth distribution of the EQs in Greek territory, with magnitude of 7.5R

   The study of the depth distribution of the EQs for the range of magnitudes from 5.5 R to 7.5 R indicates the existence of preferential depths in the lithosphere, where the EQs take place.
 

2.2.6. Number of EQs as a function of magnitude for a fixed depth of 0, 5, 10Km.

Fig. 2.2.6.1. Number of EQs as a function of magnitude for fixed depths of 0, 5 and 10Km.

   In this graph is made a comparison between the No of EQs which occur at an intermediate depth of 5Km and the depth of 10Km. It is evident that, in their majority, the EQs present a tendency to occur at the depth of 10Km of the seismic zone.
 

2.2.7. Number of EQs as a function of magnitude for a fixed depth of 100Km.

   Finally, a similar analysis has been made for the depth of 100Km.

Fig. 2.2.7.1. Number of EQs as a function of magnitude, for a fixed depth of 100Km.

   Although the number of EQs is rather small, at the depth of 100Km, it is evident that mainly (more than 50%) strong events occur at this depth.
   The seismic zones, identified by the previous diagrams, coincide with the depth of P seismic (Fortsch and Moho) velocity discontinuities, observed, by Mueller and Landisman (1966) and with the S shear, seismic wave velocity discontinuity, observed, at a depth of around 100 Km, (Toksoz et al. 1967).
   In conclusion, it is expected that strong EQs which take place in a seismogenic area behave, probably, differently, as far as it concerns the generation of precursory phenomena. Therefore, it is possible that some times precursory phenomena may not be generated or these are undetectable with the present state of physical knowledge and technology, available to the scientists.
 

2.3. Seismic energy density at distance (x) from the focal area.

   During the occurrence of a strong EQ a large amount of strain energy, accumulated and stored in the regional focal area, is released. That energy propagates from the focal area outwards, in a more or less spherical surface mode, until it reaches a seismic wave velocity discontinuity. In such a case it can be reflected, refracted or diffracted, depending on the specific physical - tectonic conditions, met in the area of incidence.
   The energy, which is released through an earthquake, obeys the simple laws of physics. The application of these laws in seismology explains, in a very simple way, some well known notions as: a) a rocky ground is safer than a sedimentary one, b) different ground acceleration is observed in very close distances. These phenomena and observations can be explained by energy manipulation in terms of physical laws, as it will be explained in the text to follow.
   The actual, local, seismic intensity felt, due to a distant earthquake, depends, in its simplest approach, on the energy / surface unit that “in-flows” at this area. Assuming homogeneous Earth, the largest the distance between the focal area and the affected area is, the less the seismic intensity is felt in this area. This is demonstrated in the following figure (2.3.1).

Fig. 2.3.1. Seismic energy that outflows at a solid angle from the focal area and passes through S1 and S2 remote surfaces of different locations and referring to the same solid angle.

Assuming that the outflow energy in a solid angle is (E) and


                                                       S1 = πR1², S2 = πR2²                                                               (2.3.1)

   Then:

                                                       I₁ = E/S1 and I₂ = E/S2                                                              (2.3.2)


   Where I_{1, 2} are the seismic energy per surface unit, in other words the intensity felt at each area 1, 2.

   Substituting S1, 2 in equation (2.3.2) by the expressions of equation (2.3.1) we get

                                                           I₂/I₁ = (R₁/R₂)²                                                                        (2.3.3)

   And by taking into account that R = (solid angle/2) * x

   Where: x stands for the distance of each area from the hypocenter of the EQ, then the equation (2.3.3) results into:

                                                           I₂ = I₁*(x₁/x₂)²                                                                         (2.3.4)

   Equation (2.3.4) indicates that the seismic intensity, which is felt at two different sites from the same EQ, is inversely proportional to the square of the ratio of their distances from the hypocenter of the EQ.
   The same is true for the acceleration which is observed in an area, registered, during the occurrence of an earthquake. However, sometimes the accelerogram indicates larger acceleration than what is normally expected. This is due to local geological conditions and is especially caused by lateral, geological discontinuities, when the seismic energy propagates through material of different densities.
   A simple physical - mathematical explanation is as follows:
   It is assumed that a discontinuity exists between a geological formation of d1 density and another one of d2. Adjacent the discontinuity, we consider the same unit volume M in both geological formations (fig. 2.3.2).

Fig. 2.3.2. Lateral discontinuity sketch model. Seismic energy flows from medium of density d2 towards medium of density d1. At the discontinuity an equal volume M is activated in both sides. It is assumed that there is no energy loss during this process.

   As long as seismic energy propagates, along the geological formations, the same quantity of energy flows through the adjacent unit volumes at the discontinuity. The latter is expressed by the equation (2.3.5).

                                              Seismic energy = ½ d₂M(V₂)²                                                       
                                                                                                                                                             (2.3.5)
                                              Seismic energy = ½ d₁M(V₁)²

   Where V1, 2 are the motion velocities of the two (M) volume masses of the different density media.

   By equating the second part of equations (2.3.5) it results into:

                                              V₁ = V₂*(d₂/d₁)^1/2                                                                            (2.3.6)

   Indicating that, for the same amount of energy transfer through a lateral discontinuity, the relative motion velocity, in the two different media, depends on the square root of the ratio of their densities. This is the physical explanation of the “amplification” of the seismic waves, observed, over loose materials. It is nothing more than the utilization of the conservation of energy which flows from one geological formation to the other, under the assumption that there are no other significant energy losses of any kind.
   Furthermore, the seismic energy is transferred through an oscillating seismic wave of angular velocity:

                                                          V_ω = V₀*sin(ωt)                                                                      (2.3.7)

   Therefore, the ratio V_ω1 / V_ω2 of the angular velocities of the two media is expressed by the equation:

                                             V_ω1 / V_ω2 = [V₀₁* sin(ωt)] / [V₀₂* sin(ωt)]                            (2.3.7)

                                                     V_ω1 / V_ω2 = V₀₁ / V₀₂                                                               (2.3.8)

   Where: V₀₁ and V₀₂ are the maximum values of the velocities, observed, in the two geological formations.


   Since:

                                                               V₀₁ ≠ V₀₂                                                                            (2.3.9)

   Then:

                                                     dV₀₁/dt ≠ dV₀₂/dt                                                                       (2.3.10)

   The latter equation (2.3.10) indicates that the acceleration which will be observed in the two geological media will be different for the same activating seismic wave. Moreover, the difference in the observed acceleration will depend on the ratio of their densities, as it is indicated by equation (2.3.6).
   The already presented topic is related to the seismological notion which more or less suggests: “a large earthquake which occurs in the sea is not so dangerous for the nearby towns, since most of the energy is absorbed by the seawater”. The latter is presented in the following figure (2.3.3).

Fig. 2.3.3. Sketch diagram that represents the seismic energy propagation from hypocenter, through its traveling wave fronts, towards a nearby town.

   The seismic energy that arrives to the town region is larger than what was absorbed by the seawater. This is demonstrated in figure (2.3.3). The energy that travels along the path EQ hypocenter – town and is carried by the seismic wave fronts is larger than what can be absorbed by the thin (compared to the EQ – town distance) seawater layer. Moreover, the S waves that carry part of the total energy, cannot propagate in the seawater. The vast majority of the seismic energy propagates along the referred EQ – town line path through solid ground. Therefore, all the latter physical laws are valid and applied.
In conclusion, a strong earthquake can be destructive and always dangerous for the nearby towns, regardless of its epicentral location in seawater or on solid ground.
 

2.4. Electrical resistivity lithospheric model.

   In section (2.2) the physical parameters of P and S seismic wave velocity were presented as a function of depth in the crust and upper mantle. Although these parameters prescribe more or less the physical behavior of the medium, another property of it, the resistivity, prescribes its electrical behavior.
   It is well known, in the geological – geophysical sciences, that each geological formation presents a range of resistivity.
   Therefore, it is possible, by obtaining knowledge of the depth distribution of the resistivity of a geological area, to assign to it the most probable geological formations. In areas which are more tectonically complicated it is possible to model complex tectonics, just by using the appropriate electrical methodology. Details on this topic can be found in any textbook for Applied Geophysics in the corresponding chapter of the application of electrical methods.
   The depths, involved for the study of the crust and the upper mantle, are of the order of some decades of Kilometers. Therefore, the magnetotelluric methods were used, initially, in order to achieve large depth penetration of the current, used, for this kind of studies. The key feature of these methods is the depth penetration (skin depth) of an electromagnetic wave in the ground which depends on its frequency and the resistivity of the ground. The latter is expressed by the following equation (2.4.1).

                                                                   Z = 500*(ρ/f)^1/2                                                               (2.4.1)

   Where (Z) is the penetration depth in meters, (ρ) is the ground resistivity in Ohm*m and (f) is the frequency of the electromagnetic wave in cycles/second. A main drawback of these methodologies is the limited band pass (T<1hour) of the filters, used for registering the electric and magnetic field, while at the same time, for larger periods signals, a very long time is required for an effective registration. An alternative to this methodology is the traditional Schlumberger electrical sounding of the AMNB electrode configuration. In such a kind of electrical arrays the length AB of the electrical dipole which induces the DC current into the ground, controls the current depth penetration and therefore, the depth of the investigation.
Results obtained from such a kind of operations were reported by Blohm and Flathe (1970) after having used a 150Km AB electrical dipole, Blohm (1972) and Homilius and Blohm (1973) in the Rhine Graben area. Other researchers, using a larger AB dipole length (up to AB = 600Km) in a different tectonic province (Southern Africa) investigated the resistivity distribution in the Earth to larger depths. In this frame of work Van Zijl (1969), Van Zijl et al. (1970) studied the crustal conductivity structure of South Africa, Van Jijl and Joubert (1975) for the same purpose used AB length of 450Km and 400Km. In 1976 Van Jijl reported results, obtained, from 30 deep electrical soundings (AB spacing = 40Km) carried out at the Umtali – Pietersburg area (South Africa) while Blohm, Worzyk and Scriba (1977) presented results, obtained, from a Shlumberger sounding with AB = 1250Km which is the longest used, so far, and known in the literature that deals with this geophysical methodology.
   The unique graph of the apparent resistivity as a function of AB/2, obtained, from this experiment, is presented in the following figure (2.4.1).

Fig. 2.4.1. Apparent resistivity as a function of depth (AB/2), obtained, from the application of the electrical sounding methodology in Southern Africa by using the Cabora Bassa power line (Blohm, Worzyk and Scriba, 1977) as a power source.

   The transformation of the apparent resistivity depth function into the corresponding resistivity layers revealed the following resistivity layering, presented as model-1 and model-2, which models are both valid for the area of the study. The first model, model-1, is based on a four-layer composition of the lithosphere, as it is presented in the following figure (2.4.2).

Fig. 2.4.2. Crust resistivity model – 1 (after Blohm et al. 1977).

   A more detailed and equally valid model is the following model-2, which is composed by five layers (fig. 2.4.3).

Fig. 2.4.3. Crust resistivity model – 2 (after Blohm et al. 1977).

 
   By taking into account the equivalent principle on the interpretation of the resistivity graph, it is possible to accept more different “resistivity interpretations”, even composed with more electrical layers, but these two have the following, general characteristics:
  A resistive zone with a maximum thickness of 24Km is followed by an intermediate, conductive layer which has a thickness of about 18Km. The next layer exhibits a high resistivity and a thickness of 120Km. The resistivity of the last layer is not defined very well, but it is definitely lower.
   The resistivity of the upper crustal layer seems to be rather high (ρ>=50000 Ohmm), while the conductive layer, found, in the lower crust, coincides with the results, obtained, by the use of the magnetotelluric method (Keller, 1971; Adam, 1976). The resistivity of the upper mantle (8000 Ohmm) is in reasonable agreement with laboratory results (Brace, 1971), while the thickness of the lithosphere (160Km) agrees very well with the values 150 – 175Km, obtained, by teleseismic delay times (Fairhead and Reeves, 1976). The final layer has low resistivity (ρ<= 50 Ohmm) in accordance to the results, obtained, by the deep magnetic sounding method (Schmucker, 1974).
   A comparison between (P) and (S) velocity depth distribution and the resistivity in the crust indicates that, there is a good agreement of Moho depth, suggested by the P velocity distribution and the thickness of the resistive layers of the upper crust. Specifically, Moho depth is clearly indicated by the bottom of the conductive layer of the 3rd layer in model-2, while the (S) velocity depth distribution coincides, within the first 200Km, with the high resistant package of the crust.
   Summarizing these results, it is evident that the seismogenic zone in the crust behaves, particularly, as a highly resistant medium in terms of its electrical properties. The latter, combined with its mechanical properties, is very important in the overall development of the precursory phenomena which are used in earthquake prediction methodology.
   A final point which must be mentioned is the generalized electrical model of the crust. In terms of geological formations, the very few first kilometers of the crust from the ground surface, consist of mainly sedimentary formations. This implies that the top few first kilometers of the crust consist of highly conductive material, followed, by highly resistant layers up to a depth of about 160-170Km and finally, followed, by a highly conductive basement. This is presented in the following figure (2.4.4).

Fig. 2.4.4. Generalized crust conductivity model.

 
   This electrical conductivity sequence can be viewed as an electrical system which is composed on an electrical shield (top and bottom conductivity) which “shields” the inner resistant conductor from external, electrical, electromagnetic influences. The latter will be speculated about in more details in the sections to follow.
 

2.5. Physical models, used, in this methodology.

   The issue of the earthquake prediction has raised strong debates in the scientific community. One of the arguments, which are presented against it, is the “absence of a physical model” that will justify any methodology to this end. In other words, there is no physical mechanism, universally agreed, for the seismogenesis.
   The different physical models, which have been proposed up to date, such as the “Rebound Model” (Reid, 1911), the “Self Organized Criticality (SOC) model (Bak et al. 1988, Bak 1996), the Chaotic Non-linear Systems (Anderson 1990, Kagan 1994, Main 1996), highly depend on the initial conditions, thus suggest the unpredictability of EQs. Other studies have shown that a small EQ can grow into a strong one, depending on very small variations on elastic properties, fault parameters variations as friction, stored, strain energy (Otsuka 1972a, b; Vere-Jones1976, Ito et al. 1990).
   Moreover, the fault friction plays the most important role in seismogenesis, in the frame of the tectonic regime of the lithosphere, while fracturing is a secondary one (Scholz, 1998). The term “stick-slip frictional instability”, introduced, by Brace and Byerlee (1966), suggests that strain energy is accumulated during the “stick” period, while an EQ occurs at the “slip” period. Actually, the phenomenon of seismogenesis is a mixture of frictional slip failure and shear fracture (Ohnaka, 2003).
   The seismogenesis mechanisms, postulated to date, refer to the physical-mechanical procedures that take place in the focal area, before and during an EQ occurs. In a most indirect way they refer to EQ precursors and therefore, to the capability of prognostic parameters determination. An example of such a calculation is the “power law - time to failure model” that correlates the magnitude (M) of a future strong EQ to the remaining time (tr) towards its occurrence (Bufe and Varnes 1993, Bowman et al. 1998).
   The absence of a valid and robust relation between the different seismogenic mechanisms and the seismic precursors, required, for the determination of the prognostic parameters (location, magnitude, time), led the seismologists to apply statistical methods for the issue of the earthquake prediction. Methodologies, as the algorithms M8 (Keilis-Borok et al. 1990, Healey et al. 1992, Romashkova et al. 2002) and CN (Keilis-Borok et al. 1990a), have been applied with some success, mainly for long-term and medium-term prediction for rather large areas.
   These statistical and any other prognostic algorithms fail to satisfy the postulated logical equation (1.1). At this point it must be stressed out that, what is really needed for a successful earthquake prediction is a physical model upon which the appropriate calculations will be based in order to analyze the precursory available data, so that earthquake prediction can be implemented. The difference between these physical models and the physical, seismogenic mechanisms is essential. The seismogenic mechanisms refer to what actually happens in the seismogenic region and therefore, they are a physical “close-up” view in the seismogenic area, while the physical models, used, for the calculation of the prognostic parameters, are a more generalized, physical approach of the seismogenic region, and they are independent from the actual seismogenic mechanisms which take place in it.
   Moreover, it is anticipated that, a single physical mechanism cannot provide answers for all earthquake prognostic parameters. This is evident from the large number of publications, which are related to the topic of earthquake prediction. The majority of them refer only to time, or regional area in relation to magnitude. None of them deals, simultaneously, with all the prognostic parameters.
   In the methodology to be presented, a different approach was followed. During the course of this research (1981 – 2003) it was found that each prognostic parameter required, a different generalized physical model, to be used. Each one of them fulfilled the logical equation (1.1). Therefore, the simultaneous use of them, applied on the appropriate precursory data, facilitate the implementation of the sort-term earthquake prediction.
These models are the following:

   a. The lithospheric seismic energy flow model. This is a direct application of the energy conservation law of physics. The physical system, in concern, is the lithospheric, seismogenic region. This model is used for the calculation of the magnitude of an imminent, strong EQ. Furthermore it explains, in energy transfer terms through the lithosphere, the “accelerated deformation” and “seismic quiescence” methodology, used, by the seismologists.

   b. The oscillating lithospheric plate model. Lithosphere is treated as a plate, which is forced into oscillation by the tidal forces. This methodology is in use by geophysicists, for the correction of gravity measurements, due to errors, induced, by the tidal oscillation of lithosphere. The latter is used for the time of occurrence determination in conjunction with the homogeneous ground model.


   c. The homogeneous ground Earth model. This is the start-up model in the analysis of the electrical methods, used in Applied Geophysics. It is used for the determination of the epicentral area.

   It must be mentioned, here, that all these models have been long ago, well-known, physical models, applied on Earth, for the application of geophysical studies for different geological and tectonic targets.
   The detailed use of each one of these physical models will be presented in the prognostic parameters calculation.
 

2.5.1. The lithospheric seismic energy flow model.

2.5.1.1. Theoretical analysis.

   It is generally accepted that stress energy built-up, in a focal area, is a very slow process which, closely, follows the motion of the lithospheric plates. It takes a long period of time (probably a large number of years) to reach the point when an earthquake will occur, because of rock fracturing.
   Under normal conditions, the stored energy is discharged through the background, small magnitude, seismicity of the seismogenic area. The latter is demonstrated in figure (2.5.1.1.1).

Fig. 2.5.1.1.1. The postulated model of seismic energy flow, through a seismogenic area of lithosphere.

   When a strong earthquake is in preparation, before its occurrence, the normal seismicity either decreases for a certain period of time and therefore, the seismic “quiescence”, which is detected, is used as a precursory indicator or increases and therefore, the “accelerated  deformation” is observed.
The mathematical analysis of the postulated model is as follows:

   Let us denote as (E_in) the inflow energy over a short period of time (dt) in a seismogenic area of the lithosphere and as (E_out) the energy outflow of the same area, due to its seismicity for the same period of time. (E_in) and (E_out) are time functions, E_in = E_in(t) and E_out = E_out(t).
   The seismogenic area is charged with energy E_st = E_st(t) as follows:

                                                            E_st = E_in – E_out                                                             (2.5.1.1)

   For any very short period of time (dt), equation (2.5.1.1) takes the following form:

                                                          E_st(dt) = d(E_st(t))                                                             (2.5.1.2)

   For successive values of dt, a discrete valued function Y(t) is defined as follows:

                                                           Y(t) = E_st(t)                                                                       (2.5.1.3)

   Equation (2.5.1.3) can take the following forms:

     a. Y(t) = C = 0                                                                                                                         (2.5.1.4)

   In this case, the stored energy in the seismogenic area, equals to zero, therefore, it is at a state of stable zero-charge conditions.

     b. Y(t) = C < 0                                                                                                                         (2.5.1.5)

   In this case, the stored energy decreases, therefore, the seismogenic area continuously discharges, towards a very stable state of uncharged conditions.

     c. Y(t) = C >0                                                                                                                          (2.5.1.6)

   In this case the stored energy increases continuously, therefore, the seismogenic area is charged towards a state of highly unstable conditions, leading to the occurrence of an earthquake.

   The term (C) of equation (2.5.1.6) can generally be either time dependent or time independent.

     - C is time independent:

   Equations (2.5.1.2) and (2.5.1.6) can be combined in equation (2.5.1.7).

                                                                        d(E_st(t)) = C >0                                                      (2.5.1.7)

   The cumulative energy (E_cum), stored in the seismogenic area, can be calculated as a function of time (E_cum(t)) by integrating both sides of equation (2.5.1.7), in respect to time (t),

                                                                   E_cum (t) = òd(E_st(t))dt = òcdt                                (2.5.1.7a)

   Since term C is constant and time independent, the calculated function E_cum(t) has the form of:

                                                                 E_cum (t) = C*t + b                                                        (2.5.1.8)


   Where, b denotes the integration constant.

   The linear equation (2.5.1.8) was firstly introduced by Thanassoulas et al. (2001), along with the postulated, lithospheric seismic energy flow model, for the calculation of the maximum magnitude of an imminent, strong EQ.


- C is time depended:

   In this case, equation (2.5.1.7) can be represented by an nth-order polynomial:

                                             d(E_st(t)) = a_ntⁿ + a_n-1t^n-1 + ......... + a₀                                      (2.5.1.9)

The cumulative energy E_cum(t), stored, in the seismogenic area, can be calculated as a function of time by integrating in time, both sides of equation (2.5.1.9).

                                E_cum(t) = òd(E_st(t))dt = ò( a_ntⁿ + a_n-1t^n-1 + ......... + a₀ )dt           (2.5.1. 9a)

   As an n+1 order polynomial,

                                             E_cum(t) = k_n+1 tⁿ⁺¹ +k_n tⁿ +………k₀                                    (2.5.1.9b)

   Where, ki values represent the polynomial constants.

   Summarizing the forms the equation E_cum(t) (that represents the energy flow) takes, the following cases are possible:

   a. linear polynomial - constant energy flow

   b. higher order polynomial - real acceleration

   c. accelerated for a period of time long before the main, seismic event, which is followed by a constant energy flow , just before the occurrence of the strong EQ.

   Cases (b) and (c) are represented, by the well-known “time to failure” function, very often.
  Although, mathematically, it is possible, to transform any polynomial to any arbitrary function, i.e. time to failure function by calculating the appropriate parameters of the latter, by using LSQ techniques, still remains the parameter of arbitrarity, as far as it concerns the validity of physics behind this transformation.
   Moreover, the time to failure function depends on two variables. The first one is the magnitude of the imminent EQ and the second one is the time to failure, left. In order to overcome the problem of solving a two parametric equation (infinite number of solutions), the parameter C (Bowman et al. 1998) was introduced, that is the ratio of power low fit error over the linear fit error, as far as it concerns the cumulative seismic energy release. Still, the notion of this ratio is set completely, arbitrarily.
   Therefore, it is suggested that the magnitude and time to failure of a strong, imminent EQ, calculated, by these methodologies, are not supported by any validated, physical mechanism and should be rejected.
 

2.5.1.2. Application of the model on real EQ cases.

   The application of the lithospheric, seismic energy flow model requires two parameters to be known, in advance. The first one is an initial estimation of the area extent of the physical system itself (seismogenic area). The second one is its seismic history.
   The choice of the first parameter is utilized by a) knowing the epicentral area of the future EQ by another methodology, as it will be explained later on b) by taking into account the sometimes observed increased low level seismicity of an under study seismogenic area, which is some times observed and finally c) by taking into account the deep, lithospheric fracture zones which are mapped by the conversion of the regional gravity map of the seismic region in to a gradient one.
   This procedure is explained by the figures to follow.
   As long as the seismogenic area has been activated, a small magnitude seismic activity, associated with the main rupture of the rock formation, increases in general, in the different order fracturing branches (nth order ridel) and therefore, this seismic activity is observed along and very close to the trace of the main seismic to be activated fault. The latter is presented in the following figure (2.5.1.2.1).

Fig. 2.5.1.2.1. The main, seismic fault branches into smaller ones, as it approaches the ground surface (Mattauer 1973, Vialon et al. 1976). The group of fractures that advances up to the ground surface forms the observed “fracture zone”.

   Since small-scale seismicity has started to emerge, the corresponding foci will be located close to the vicinity of the main fault, which will be activated in the future. A sketch drawing indicating this main earthquake precursory activity is presented in the following figure (2.5.1.2.2).

Fig. 2.5.1.2.2. Sketch drawing of the theoretically small magnitude seismicity, expected, to be observed on ground surface and in the vicinity of a main, seismic fault, prior to a strong EQ. The solid dots represent the epicenters of the precursory, seismic activity, while the thin lines indicate the nth order “ridel” different faults. The thick, black line represents the main fault which is going to be activated.

   Consequently, the physical system which will be studied must be chosen, close, to the main fracture, in a way that it takes into account this specific distribution of the small magnitude EQs. The second parameter, the seismic history of the seismogenic area, is taken from the EQ files from the Seismological State Observatories. In case of Greece, this is the National Observatory of Athens (NOA). The application of the methodology is shown in the following figures. It is assumed that, the epicentral area has been approximately determined by other means. As a first step, is defined the seismogenic area which will be considered as the corresponding physical system of which the seismic energy release will be analyzed (fig. 2.5.1.2.3).

Fig. 2.5.1.2.3. The approximately, determined seismogenic area is indicated by the blue polygon. It has been drawn in such a way to center the associated deep, main, lithospheric fracture zone and is elongated along it.

   The next step to take is to calculate backwards the cumulative energy release of this area for a long period of time. The graph which results indicates whether this area has entered the acceleration phase or not. The latter is presented in the following figure (2.5.1.2.4).

Fig. 2.5.1.2.4. Cumulative seismic energy release calculated, for the seismogenic area, which is indicated in the previous figure (2.5.1.2.3) for the period 1998 – 2006.

   The study of this graph indicates that this example, seismogenic area has been set into cumulative, seismic energy release acceleration mode for the last 8 years. A slightly different presentation is shown in next figure (2.5.1.2.5) by fitting a 6th degree polynomial. This facilitates the analytical calculation of the time gradient of the cumulative seismic energy release.

Fig. 2.5.1.2.5. Cumulative, seismic energy release fits with a 6th degree polynomial. The black line indicates the cumulative, seismic energy release data, while the red one indicates the fitted, polynomial values.

   It is obvious that the acceleration has been initiated from the end of 2005, indicating that an earthquake will strike “soon”. Actually, a strong (M = 6.9R) earthquake took place at a short distance towards NW on this lithospheric fracture zone at the start of 2006. Its details will be presented in the examples to follow.
   Quite often, the start of the increase of the cumulative energy release is not defined very sharply, but there is a gradual change over a rather lengthy period. In such cases, instead of fitting only a polynomial, as an advanced step, the corresponding gradient of the polynomial is calculated analytically.
   The process of gradient calculation acts as a high-pass filter on the cumulative seismic energy release data, low-order terms of the fitted polynomial are excluded by this operation and therefore the resolving capability for detecting smaller level changes is larger.
This process is shown in the following figure (2.5.1.2.6).

Fig. 2.5.1.2.6. Cumulative, seismic, energy release time-gradient, calculated analytically, from the polynomial of figure (2.5.1.2.5).

2.5.1.2.1. Application of the theoretical model on the EQ in Zakynthos (02/12/2002,

                Ms = 5.8R).

   The first example, based, on the methodology which has been already presented, is that of Zakynthos EQ. A few months before the main, seismic event which took place on 02/12/2002, with a magnitude of Ms = 5.8R, an increased, small-scale, seismic activity had been observed at Zakynthos area. The area of interest is indicated by a red circle in the following figure (2.5.1.2.1.1).

Fig. 2.5.1.2.1.1. The red circle indicates the regional area of Zakynthos. The blue circles indicate the small magnitude seismic activity all over Greece for a time period of a couple of months before the EQ occurred.

   The seismicity which was observed to increase, dictated the application of the latter methodology. The area, to be investigated, was selected in such a way, to include the main, lithospheric fracture zones, which are present in the region. The latter is presented in the following figure (2.5.1.2.1.2). The blue frame indicates the area of interest.

Fig. 2.5.1.2.1.2. Area of interest (blue frame) for which the lithospheric seismic energy release methodology will be applied. The brown lines indicate the deep lithospheric fracture zones.

   Next step is to calculate backwards the cumulative seismic energy release for some years. The latter was utilized by using the EQ file, which is available online in the web by the NOA, Athens, Greece. This operation indicated that the area of Zakynthos was set in seismic acceleration mode from the start of the year 2000. This is presented in the following figure (2.5.1.2.1.3).

Fig. 2.5.1.2.1.3. Cumulative seismic energy released, from the seismogenic area of Zakynthos for the period 2000 – 2002.

   The 6th order polynomial, fitted in the cumulative, seismic energy release data for the area of Zakynthos, is presented in the following figure (2.5.1.2.1.4).

Fig. 2.5.1.2.1.4. The 6th order polynomial fitted, in the cumulative, seismic, energy release data for the area of Zakynthos. The brown line indicates the cumulative, energy data, the blue one indicates the fitted polynomial. The red arrow indicates the time when the Zakynthos EQ occurred.

   The analytical expression of the polynomial which is obtained through the LSQ fitting procedure, utilizes the calculation of the time gradient of the polynomial and therefore, the rate of change of the cumulative, seismic energy release. This is presented in the following figure (2.5.1.2.1.5).

Fig. 2.5.1.2.1.5. Rate of change in time of the cumulative, seismic energy release observed, prior to the EQ of Zakynthos. The red arrow indicates the time when this EQ occurred.

   It must be pointed out that all this analysis had been made before (one month) the EQ occurred and therefore, this EQ was expected to occur soon. The, expected, EQ did happened within a month and its location coincides with the location of the deep lithospheric fracturing which is located in this area (fig. 2.5.1.2.1.6),

Fig. 2.5.1.2.1.6. Location of the EQ in Zakynthos in relation to the deep, lithospheric fracture zones, which are identified in the same area, by analyzing the gravity field.

   thus verifying the validity of the corresponding, deep, lithospheric fracture zones / faults map (Thanassoulas, 1998).
 

2.5.1.2.2. Application of the theoretical model on the EQ in Kythira (08/01/2006,

                Ms = 6.9R).

   This example is an “a posteriori” one. The validity of the method was tested against Kythira EQ (08/01/2006, Ms = 6.9R), which is presented in the following figure (2.5.1.2.2.1).

Fig. 2.5.1.2.2.1. Location (blue concentric circles) of the EQ in Kythira (08/01/2006, Ms = 6.9R), in relation with the location of the deep, fracture zones / faults of the lithosphere.

   The seismogenic region, which is taken into account, is presented as a blue polygon frame in the figure (2.5.1.2.2.2) below.

Fig. 2.5.1.2.2.2. Area of interest (blue frame) for which will be applied the lithospheric, seismic energy release methodology. Brown lines indicate the deep, lithospheric fracture zones.

   The cumulative, seismic energy release, which is calculated for a period of 14 years (1992 – 2006), indicates that this seismogenic region was set in acceleration mode for the last 2 years (2004 – 2006). This is presented in figure (2.5.1.2.2.3).

Fig. 2.5.1.2.2.3. Cumulative, seismic energy release determined, for the seismogenic region of Kythira, for the period 1992 – 2006.

   The fitted, 6th order polynomial function of cumulative, released, seismic energy vs. time, indicates a rapid increase “graph knee”, almost 2 years before the corresponding earthquake occurred. This is presented in next figure (2.5.1.2.2.4).

Fig. 2.5.1.2.2.4. Cumulative, seismic energy release is fitted by a 6th degree polynomial. The red line indicates the cumulative, seismic energy release data, while the blue one indicates the fitted, polynomial values. The black arrow indicates the time when the EQ occurred.

   In this case, the time gradient, which is calculated from the analytical expression of the polynomial, resolves the time, when the seismogenic area entered the period of intense seismic acceleration in a much better way. The latter is shown in next figure (2.5.1.2.2.5).

Fig. 2.5.1.2.2.5. Rate of change, in time, of the cumulative seismic energy release observed, prior to the EQ in Kythira. The black arrow indicates when the EQ occurred.

   A point that must be clarified, more, is the way it is set the seismogenic area extent. What is accepted to date in the seismological community is that, the stronger the pending EQ is, the larger the area which is overcharged with strain deformation around the epicenter is. The latter is represented by the following formula:

                                                             log R = 0.42M – 0.68                                                  (2.5.1.10)

   where: (R) is the radius of the strain charged area and (M) is the magnitude of the pending earthquake (Papazachos et al. 2000, Papazachos et al. 2001).
   A radius of almost 200Km has been estimated, after having adopted this formula, for the assignment of the seismogenic region which is going to be used for the calculation of the seismic energy released for the case of the EQ in Kythira (Ms = 6.9R). As a first approach, this area is approximated by a surface, enclosed in a, closely, orthogonal frame of dimensions: 4 degrees (latitude) by 5 degrees (longitude). The latter is presented in the following figure (2.5.1.2.2.6).

Fig. 2.5.1.2.2.6. Seismogenic area (blue frame) theoretically charged with strain deformation, required for the generation of the EQ in Kythira (Ms = 6.9 R).

   It is clearly obvious that the seismogenic area, determined, by the use of equation (2.5.1.10) will be affected, not only by the main fault where this strong event took place, but also by the adjacent ones.
   The latter will contribute, with their seismic energy release, to the one, calculated for the entire frame. We must avoid it, when such calculations are made, since the seismic energy release of the fault which is expected to be activated, is studied. Moreover, calculations made for the identification of the accelerated deformation are masked by the presence of faults, which will not be activated at all. This is demonstrated in the following figures.
Let us assume a 2 by 2 degrees frame (fig. 2.5.1.2.2.7) that encloses some more main fracture zones, besides the one that was triggered.

Fig. 2.5.1.2.2.7. A 2 by 2 degrees seismogenic area (blue frame), initially, assumed, for the cumulative, seismic energy release determination.

   The cumulative seismic energy release is calculated as a function of time, for the same period (1992 – 2006) as the EQ of Kythira. This is presented in the next figure (fig. 2.5.1.2.2.8).

Fig. 2.5.1.2.2.8. The cumulative, seismic energy release calculated, for the same period (1992 – 2006) as for the EQ in Kythira.

   This graph fits quite well a straight line, except for the period 1998-2002, when a moderate magnitude seismic event took place.
   There is no evidence to show that any, accelerated, deformation mechanism has been initiated. The same is obvious from next figure (fig. 2.5.1.2.2.9).

Fig. 2.5.1.2.2.9. Polynomial fitting (6th order), as determined on the cumulative seismic energy release data.

   For the next example, was used a smaller frame (1 by 1 degree), which encloses only the activated area, but does not take into account the strike of the seismically, active fault. This frame is presented in next figure (2.5.1.2.2.10).

Fig. 2.5.1.2.2.10. A 1 by 1 degrees seismogenic area (blue frame), initially, assumed, for the cumulative, seismic energy release determination.

   The corresponding, cumulative, seismic energy release graph is shown in the following figure (2.5.1.2.2.11).

Fig. 2.5.1.2.2.11. Cumulative, seismic energy release, determined for the 1 by 1 degrees seismogenic region, for the period 1992 – 2006.

   The accelerating deformation characteristics of the seismic energy are initially revealed in this graph. This is more obvious in the following figure (2.5.1.2.2.12) where was used the polynomial, fitting, procedure.

Fig. 2.5.1.2.2.12. Polynomial fitting (6th order) as determined on the cumulative seismic energy release data.

   The same accelerating deformation character of the seismicity of the corresponding frame is indicated by the time gradient graph (fig. 2.5.1.2.2.13) of the cumulative seismic energy release data.

Fig. 2.5.1.2.2.13. Rate of change in time, for the corresponding 1 by 1 degrees seismogenic area, of the cumulative, seismic, energy release, observed.

   The comparison of the three cases presented (2 by 2, 1 by 1, and frame oriented along fault strike) indicates that, the most indicative case, for the detection of the accelerating deformation status of a seismogenic region, is the one, where the assumed regional seismogenic area (frame) conforms to the tectonic character of the seismogenic area and the strike of the seismogenic fault.
   The latter leads us to the following conclusion. In the wider regional area, where a strong EQ will take place, the equation:

                                                            logR = 0.42M – 0.68                                                       (2.5.1.11)

   is valid and defines the strain, deformed, area.

   However, in the vicinity of the fault, to be activated, the most drastic, elastic deformations take place and this explains why the EQ takes place in this fault, during the strain release. Consequently, this leads directly to the Rebound Theory of Reid. After having presented the above examples, what is proposed and probably holds in the seismogenic area is the following:

     - Up to a distance (R) from the epicenter location, calculated by the equation (2.5.1.11),

       during the preparation of a strong EQ, strain is accumulated.

     - Very close to the fault which will be activated, the elastic deformation prevails and

       the dominant model for the earthquake generation mechanism is the one of the

       Rebound Theory.

   It is obvious that, the accelerating deformation status of each deep, lithospheric fracture zone must be studied in detail. Further more, it is possible to compile maps which will indicate the acceleration deformation spatial distribution all over a wide area, and hence the corresponding, seismic risk will be evaluated.
   In conclusion, the combined use of the lithospheric fracture zones and the cumulative, seismic energy flow lithospheric model selected for detailed analysis can provide us with valuable information for the seismic energy charge of any area.
 

2.5.2. The lithospheric plate oscillating model.

   The Earth-tides mechanism was early recognized as a potential trigger for the occurrence of strong earthquakes. This approach was followed to study the time of occurrence of EQs and their correlation to Earth-tides. Knopoff (1964), Shlien (1972), Heaton (1982) and Shirley (1988) suggested the Earth-tides as a triggering mechanism of strong EQs, Yamazaki (1965, 1967), Rikitake et al. (1967) studied the oscillatory behavior of strained rocks, due to Earth-tides, while Ryabl et al. (1968), Mohler (1980) and Sounau et al. (1982) correlated Earth-tides to local micro-earthquakes and aftershock sequences.
   The following figure (2.5.2.1) demonstrates the mechanism which generates the Earth-tides. On the left are shown the forces which are applied on the Earth’s surface by the Moon or the Sun. Because of these forces, the lithospheric plate of the Earth deforms, as it is shown, on the right part of the figure.

Fig. 2.5.2.1. Forces applied on the Earth’s surface (Stacey, 1969) by the Moon or the Sun (left) and the resulted, deformation (Garland, 1971) of the lithospheric plate (right).

   Earth-tides exhibit an oscillating mode of behavior. The study of Earth-tide oscillations has shown that, basically, it consists of the following main components, presented, in the table below.

                                                              Tidal components

                                               Symbol    Name                  Period (hr)
                                              -----------    --- -----------------     ---------------
                                                    M2  Principal Lunar          12.42
                                                    S2  Principal Solar            12.00
                                                    N2  Lunar Ellipticity         12.66
                                                    K2  Lunisolar                     11.97
                                                    K1  Lunisolar                     23.93
                                                    O1  Lunar Declination     25.82
                                                    P1  Solar Declination      24.07
                                                    M1  Moon declination     14 days
                                                  Ssa  Moon declination       6 months

   Finally, an Earth-tide wave is generated with in a year’s period, because of the Earth’s motion in an ellipse with the Sun in the focus.
   The deformation of the lithosphere follows the oscillatory character of the Earth-tides. Garland (1971), Stacey (1969), Sazhina and Grushinsky (1971), study in detail this type of Earth’s oscillation, not to mention the majority of the Geophysical textbooks.
   The lithosphere, following the Earth-tide oscillatory forces, self-oscillates with maximum amplitude which varies as follows:
 

   According to Garland (1971), the maximum elevation of the lithosphere is of the order

    of 10cm

   According to Sazhina and Grushinsky (1971), the lithosphere oscillates with a maximum

   p – p value of 53.4cm for the Moon component, while for the Sun component is only 24.6cm. 

   That means that, the maximum oscillatory amplitude can reach the value of 78cm p – p, in cases, when the two components are added in phase.
   A sample of a day’s tidal oscillation is presented in the following fig. (2.5.2.2).

Fig. 2.5.2.2. Earth-tide calculated for the 18th November, 1997, Greece φ = 38⁰, λ = 24⁰.

   The vertical scale is in mgals, while the horizontal one is in minutes (1 day = 1440 minutes). In this figure it is possible to identify the K1 (23.93hr), K2 (11.97hr) and S2 (12.00hr) components.
   The next figure (2.5.2.3) represents the T = 14 days oscillation of the Earth-tides.

Fig. 2.5.2.3. Earth-tide oscillation with T = 14 days, for the time period of 05/08/2000 – 13/10/2000

   Figure (2.5.2.4) represents the yearly oscillation of the Earth-tide, along with the 14 days period oscillation, superimposed, on it.

Fig. 2.5.2.4. Yearly and 14 day’s period, Earth-tide oscillation is shown, for the year 2000.

   The next figure (2.5.2.5) represents the yearly and the 14 days period oscillation of the Earth-tides.

Fig. 2.5.2.5. Yearly and 14 day’s period Earth-tide oscillation calculated, for the period 01/01/2001 – 31/12/2002.

   Using the Rudman et al. (1977) methodology, tidal waves have been generated. As long as the oscillating, Earth-tides forces act upon the lithospheric plate, the latter is excited in the same mechanical mode of oscillation.
   According to this model, the lithosphere is considered as a rigid plate which can oscillate, because of applied oscillating, external forces. As a first approximated, simple model, is used the single beam, supported at point A and B (fig. 2.5.2.6).

Fig. 2.5.2.6. Load distribution (bottom graph) along a resting beam on two supports A, B (Papaioannou, 1952).

   Then the following equations hold:

                                           p=c (constant weight load along its axis)                                    (2.5.2.1)

                                           Hα = Hβ = pl/2                                                                                       (2.5.2.2)

                                          δχ(χ) = Dx = p(l/2 – x)                                                                          (2.5.2.3)

                                          Mx = (px/2)(l-x) (bend moment)                                                        (2.5.2.4)

   While:

                                          Mμ = max for Dμ = 0                                                                            (2.5.2.5)

   The load distribution of the beam (Papaioannou, 1952) along its length is presented by the bottom part of the figure (2.5.2.6).
   In the case of the lithospheric plate the maximum oscillating load is applied, locally, by the Earth-tides.
   Therefore, it is activated and deformed in an oscillatory mode. This is presented in the following figure (2.5.2.7). The solid line parallelogram represents the lithospheric plate at resting position, while the dashed lines show the two extreme positions of its maximum amplitude of oscillation.
   The lower part of the figure shows the stress load distribution of lithosphere, at the maximum or minimum of the Earth-tides.

Fig. 2.5.2.7. Lithospheric plate at rest (solid line) and maximum amplitude oscillation (dashed line). Bottom drawing indicates the stress-load, along the plate, during its oscillation.

   The oscillation of the lithospheric plate, due to Earth-tides, has two severe consequences:

     -The first one is that, at the peak amplitude of its oscillation, when the plate is at nearly critical stress load, it can reach the extreme conditions which are necessary for an earthquake to occur (Thanassoulas et al. 2001), provided that, the focal area itself, has been charged, enough, by the linear increase of stress, due to plate’s motion. This takes place, especially, when all the oscillating components of the Earth-tides are “in phase”.

     -The second one is that, at the very same stress load critical conditions, the focal area of an imminent earthquake generates electrical signals (Thanassoulas 1991, Clint 1999), due to its large crystal lattice deformation and to the piezoelectric properties of the quartzite content in it.
 

   An EQ will occur when the stress load that increases gradually with time, due to plate’s motion, will reach the failure stress level of the rock at the focal area. This is represented, schematically, in the following figure (2.5.2.8).

Fig. 2.5.2.8. The EQ occurs (red arrow) when the stress-load (in arbitral units) exceeds its critical level (red line).

   The Y-axis represents stress increase in arbitrary units, while the x-axis represents time in arbitrary units. The horizontal red line represents the rupture stress level of the rock, at the focal area, and the blue, inclined, line denotes the state of stress charge of the rock of the focal area. The red arrow indicates the occurrence of the earthquake when the stress load of the rock equals the critical rupture stress level.
   The following (fig. 2.5.2.9) stress charge load, increase form, takes place in the seismogenic region, by combining the two different stress increase mechanisms, the linear caused, by the plate’s motion and the oscillatory one, caused, by the Earth-tides.

Fig. 2.5.2.9. The Linear – oscillatory mode increase of stress is presented in sketch drawing.

   As long as the repetition period of strong EQs (red arrows) is of the order of some years, it is obvious that critical stress levels, for the occurrence of an EQ, will occur when daily, local, maximum peaks of the stress-load curve, reach the fracturing stress level.
In the case of a 30 years recharge period of a seismogenic area, the ratio of the recharge period to the oscillations one, is:

                                    T recharge / 1 day = 30*360/1 = 10800.                                                (2.5.2.6)

   That value of 10800 satisfies the criterion of:

                                    T of plate recharge >> T of oscillation                                                (2.5.2.7)

   Therefore, theoretically, strong EQs are expected to occur only at the maximum peaks of the Earth-tide waves, since these will, firstly, reach the fracturing stress level of the seismogenic area rock formation.
   The practical use of the lithospheric oscillating plate model will be presented at the chapter that analyzes the method of the determination of the time of occurrence of a strong earthquake.
 

2.5.3. The homogeneous ground Earth model.

   The upper part of the Earth (lithosphere – crust) consists of different lithological - geological formations which have been affected by the different tectonic processes, which took place in the past. This geological-tectonic complex, although it affects in a very general approach, the Earth’s physical properties spatial distribution, it can be considered as non-existent under certain physical assumptions. In other words Earth is assumed as a homogeneous ground. That means that its physical properties do not depend on the coordinates (x, y, z) which describe the Earth physical system.
   At a first glance, if we take into account the real geological, tectonic conditions which exist in a seismogenic region, this assumption appears to be false. However, in the specific case of the generation of an electric field in the focal area and, especially, when its wavelength is much more larger than the wavelength of the geological-tectonic elements, met, in the seismogenic area, then, for this specific case only, the seismogenic area can be considered as consisted of homogeneous ground, as far as it concerns its electrical properties.
   In physical terms, the X-ray spectroscopy is a similar example, to what was mentioned before. This methodology uses X-rays because their wavelength is comparable, in size, to the crystal lattice grid size and therefore the diffraction of the X-rays is utilized. On the other hand, the crystal lattice is non-existent for the incident large wavelength “light” and therefore, in general, it passes through it without any other consequence.
   The model of the homogeneous ground, in terms of its electrical properties, will be used in the analysis of the precursory, electrical signals, which will be presented later on. The basic, theoretical analysis of an electric field that is generated by some mechanism in a homogeneous material (in this case in the upper part of the Earth, Telford et al. 1976, Grant and West 1965, Bhattacharya and Patra 1968) is as follows:
   It is assumed that a current (I) flows in an isotropic, homogeneous medium. If (J) is the current density and (δA) is an element of surface, then the current (I) that passes through δA is: I = J* δA.
   The electric field (E) and the current density (J) are related through the Ohm’s law:

                                                                          J = σΕ                                                                      (2.5.3.1)

   Where (σ) is the conductivity of the medium.
   The electric field (E) is the gradient of a scalar potential,

          Since the current flows outwards radially, in all directions, from the point current source, the generated equipotentials which are orthogonal to the current flow lines, will be spherical surfaces, given by r = constant. This is presented in the following figure (2.5.3.1).

Fig. 2.5.3.1.  Equipotential surfaces generated, by a point current source in an isotropic, homogeneous medium (Telford et al. 1976).

   When the generated, equipotential, spherical surfaces are intersected by a plane of conductivity discontinuity, such as the ground surface in the case of the solid Earth, then equipotential circles are formed at the plane of conductivity discontinuity. This is illustrated in next figure (2.5.3.2).

Fig. 2.5.3.2. Plan view of equipotential circles, generated on the ground surface, by the intersection of the spherical, equipotential surfaces, generated by the Point Current Source (PCS) in the ground, and the ground surface plane.

   The latter is the main assumption, which is used throughout the methodology, and is applied, for the epicentral location determination of an imminent, strong earthquake.