Impacts of hydrogeological characteristics on groundwater-level changes induced by earthquakes

Changes in groundwater level during earthquakes have been reported worldwide. In this study, field observations of co-seismic groundwater-level changes in wells under different aquifer conditions and sampling intervals due to near-field earthquake events in Taiwan are presented. Sustained changes, usually observed immediately after earthquakes, are found in the confined aquifer. Oscillatory changes due to the dynamic strain triggered by passing earthquake waves can only be recorded by a high-frequency data logger. While co-seismic changes recover rapidly in an unconfined aquifer, they can sustain for months or longer in a confined aquifer. Three monitoring wells with long-term groundwater-level data were examined to understand the association of co-seismic changes with local hydrogeological conditions. The finite element software ABAQUS is used to simulate the pore-pressure changes induced by the displacements due to fault rupture. The calculated co-seismic change in pore pressure is related to the compressibility of the formation. The recovery rate of the change is rapid in the unconfined aquifer due to the hydrostatic condition at the water table, but slow in the confined aquifer due to the less permeable confining layer. Fracturing of the confining layer during earthquakes may enhance the dissipation of pore pressure and induce the discharge of the confined aquifer. The study results indicated that aquifer characteristics play an important role in determining groundwater-level changes during and after earthquakes.


Introduction
Changes in groundwater level triggered by earthquakes have been observed in many places of the world, particularly in the seismic active areas (Wakita 1975;Igarashi et al. 1992;Roeloffs et al. 1997;Grecksch et al. 1999;Chia et al. 2001;Sato et al. 2004;Shi et al. 2013). When earthquake waves pass through a well, the seismic shaking of the water column in the well and its connected aquifer results in an oscillatory change of water level (Thomas 1940;Cooper et al. 1965;Liu et al. 1989;Ohno et al. 1997;Kitagawa et al. 2006;Weingarten and Ge 2014). In the vicinity of earthquake epicenters, abrupt steplike groundwater level changes have been recorded during earthquakes (Waller 1966;Roeloffs 1988;Grecksch et al. 1999;King and Igarashi 2002;. Such a co-seismic change usually sustains from a few days to several months (Montgomery and Manga 2003), and thus can often be identified immediately after strong earthquakes.
Sustained groundwater level changes are attributed to the redistribution of stress and strain due to fault rupture during earthquakes (Wakita 1975;Roeloffs et al. 1989;Montgomery and Manga 2003). In the seismic area, groundwater level is often changed by earthquakes in addition to tide and barometric pressure. Large co-seismic changes in the vicinity of seismogenic faults, from a fall of 11.10 m to a rise of 7.42 m, have been reported . Hydrogeological characteristics are associated with coseismic changes (Rojstaczer 1988;Chia et al. 2001); however, their relationship has rarely been studied. Poroelastic theory, which couples soil deformation and pore pressure diffusion, has been proposed to describe the occurrence of co-seismic groundwater level changes (Biot 1941(Biot , 1956Rice and Cleary 1976;Palciauskas and Domenico 1989;Roeloffs 1996). The theoretical basis provides the chance to simulate the process of the generation and recovery of co-seismic changes in the aquifer under crustal deformation to understand the role of hydrogeological conditions in the earthquake-triggered groundwater level changes.
In this paper, the focus is on the field observations and numerical simulation of groundwater level changes induced by the displacement due to fault rupture during earthquakes under various hydrogeological conditions. In addition to the comparison between observation and simulation results, longterm analysis of co-seismic changes was performed to provide a better understanding of the mechanisms and factors controlling co-seismic responses of groundwater level.

Earthquakes and monitoring wells
As Taiwan is located at the convergent boundary between the Philippine Sea plate and the Eurasian plate, strong earthquakes occur frequently in the area. On average, from 1998 to 2015, approximately 3.8 earthquakes of M6.0 or greater occur every year in Taiwan. These large earthquakes occasionally triggered co-seismic groundwater level changes in monitoring wells. The detailed information of 15 earthquakes of M6.0 or greater used for this study is listed in Table 1. The epicenter and focal mechanism of these earthquakes is shown in Fig. 1. The focal mechanism describes the orientation of the fault plane and the direction of the slip. The red color represents the quadrant of compression and the white color represents the quadrant of dilation.
A dense network of monitoring well stations has been established in Taiwan for the management of groundwater resources since 1992. By 2008, there had been 747 wells placed at 348 monitoring stations in the plain area (WRA 2013). Each station consists of one to five 6-in-diameter (~15 cm) wells screened in the aquifers at different depths. Prior to the placement of monitoring wells, an exploratory borehole was drilled for hydrogeological investigation and hydraulic testing, providing comprehensive hydrogeological information for each station. The groundwater level at each well is recorded every hour, but a few wells have been selected to record water level every second for the monitoring of earthquakes.
Earthquake-triggered groundwater level changes recorded at nine monitoring wells were used in this study. The locations of these wells, where the elevations are between 10.8 and 297.7 m above sea level (asl), are shown in Fig. 1 and the well information is listed in Table 2. The stratigraphic columns of these wells, ranging in depth from 24 to 270 m, are shown in Fig. 2. Each well monitors only one aquifer, which is primarily composed of gravel or sand.
Field observations of earthquake-triggered groundwater level changes Field observations of groundwater level changes in the confined and unconfined aquifers during and after earthquakes are presented to show the impact of hydrogeological characters. Earthquakes that triggered groundwater level changes are listed in Table 1. Long-term groundwater level changes associated with earthquakes at LJ3, PD and HL are analyzed to explore the factors influencing the polarity and the magnitude of co-seismic changes.

High-frequency data
On December 26, 2006, two consecutive M W 7.0 and M W 6.9 Hengchun earthquakes occurred at 20:26 and 20:34 local time, respectively. Figure 3a shows the earthquake epicenters located in the southern offshore of Taiwan (CWB 2006). The focal mechanism suggests normal faulting for the first shock and strike-slip faulting for the second (USGS 2006). As shown in Fig. 3b, during the earthquake doublet, hourly water level data show an abrupt sustained fall of 18 cm at CS3 well tapping a confined sand and gravel aquifer. However, 1-Hz data show not only the oscillatory changes induced by passing seismic waves, but two gradual sustained co-seismic falls of 4.3 and 14.6 cm, respectively, induced by the earthquakes (Fig. 3c). The change process revealed that it took approximately 215 s after the arrival of surface waves to reach 90% of the first fall and 204 s to reach 90% of the second. The epicentral distance of CS3 is 100 km and 72 km for the first and the second earthquakes, respectively.
Figure 4 displays 1-Hz groundwater level data recorded at HL well between 12:00 November 5 and 12:00 November 6, 2009. HL is located at the Longitudinal Valley, right on the boundary between the Eurasian plate and the Philippine Sea plate where oscillatory changes in groundwater level are commonly observed due to frequently occurring earthquakes. The well screen was installed from a depth of 161 to 185 m in a silty-gravel aquifer (Fig. 2). The well site is only 300 m from the Pacific Ocean, and thus the groundwater level is influenced by the tidal cycle. Oscillatory changes in groundwater level, as shown by spikes in Fig. 4, are triggered by two swarms of earthquakes that occurred respectively in central Taiwan and the Longitudinal Valley in eastern Taiwan (Table 3). Oscillatory changes are attributed to the pressure changes of the aquifer and well water induced by the dynamic strain associated with seismic waves (Cooper et al. 1965;Liu et al. 1989). All oscillatory changes ceased within a few minutes and the water level returned nearly to its pre-earthquake level, except a sustained co-seismic fall during the earthquake occurred at 03:25 local time on November 6, 2009. The 0.6 cm sustained fall, that took about 90 min to complete the change, clearly altered the pattern of well water-level change due to the ocean tide.

Hydrogeological conditions
Sustained co-seismic changes with a sudden rise or fall of water level were often observed in wells tapping highly permeable confined aquifers. Figure 5 shows temporal changes of groundwater level at DH3 well during a M W 6.0 earthquake that occurred on March 27, 2013 in central Taiwan. The hourly water level data show the cyclic changes induced by the semi-diurnal earth tide, revealing the existence of vertical confinement of the aquifer tapped by DH3 (Fig. 5a). An abrupt sustained rise of 5 cm was recorded during the earthquake, disturbing the semi-diurnal water-level pattern. On the other hand, the 1-Hz water level data at DH3, as shown in Fig. 5b, display both the oscillatory and the sustained co-seismic changes. The sustained rise, approximately 4.1 cm, due to the earthquake took only about 5 min to reach a steady level. The rapid response of well water level to the stress change induced by the fault displacement during the earthquake is likely associated with the highly permeable aquifer tapped by DH3.  Table 1) selected for this study. The focal mechanism of an earthquake is displayed by quadrants of compression (red) and dilation (white) The co-seismic change and its recovery process was quite different at SL well during the same earthquake, as shown in Fig. 6. As SL is screened between the depths from 9 to 18 m in a shallow unconfined gravel aquifer (Fig. 2), any pore pressure change will quickly recover due to the atmospheric pressure at the water table. Obviously it is difficult to observe any step-like sustained change in response to the earthquake (Fig.  6a). The water level data recorded every 5 s show oscillatory changes, with the amplitude up to 1 m approximately, when earthquake waves pass by (Fig. 6b); however, it returned to its pre-earthquake level immediately after the passing of earthquake waves. Similar phenomena have been observed at many other wells tapping the unconfined aquifer. Generally oscillatory changes cease a few minutes after earthquakes, and thus can only be found in high-frequency water level data.
The 1999 M W 7.6 Chi-Chi earthquake is the largest land earthquake ever recorded in Taiwan. Figure 7 shows groundwater level changes in response to the earthquake at three wells. The sustained co-seismic rise was 3.72 m at HH3 well tapping a confined aquifer (Fig. 7a). The post-seismic water level, which declined slowly, did not return to the level before the earthquake even after 4 months. Figure 7b shows the temporal changes of groundwater level at TC2 tapping an unconfined aquifer. The sustained co-seismic drop was 2.47 m.  After the earthquake, the water level change took only about 12 h to recover to a steady level, much faster than that at HH3. The confinement of a rock aquifer could be breached by fracturing due to seismic shaking during earthquakes, particularly in the mountainous area (Rojstaczer and Wolf 1992;Rojstaczer et al. 1995;Sato et al. 2000;King and Igarashi 2002;Wang et al. 2004a;Charmoille et al. 2005;Elkhoury et al. 2006;King et al. 2006;Jang et al. 2008). Such a phenomenon was evidenced by the groundwater level change at LY2 well before and after the Chi-Chi earthquake. LY2 well is located on the bank of the Chingshui Creek which is approximately 5 km west of the seismogenic fault ( Fig. 1). The elevation of the well is 169.7 m and the well is screened at the depths between 90 m and 114 m in a sandstone aquifer (Fig. 2). Prior to the earthquake, LY2 was an artesian well where the hydraulic head was about 8 m above the ground. During the earthquake, a groundwater level decline of 5.94 m was recorded ( Fig. 7c). Following the step-like co-seismic fall, the water level did not show any recovery. Instead, it slowly declined with intermittent rapid falls during major aftershocks. One month after the earthquake, the groundwater level declined more than 8.5 m, approximately 50 cm below the ground, and LY2 was no longer an artesian well. The lack of water level recovery and the loss of artesian conditions after the earthquake suggested the occurrence of groundwater discharge, possibly through fractures generated by the earthquake, to the nearby creek valley.

Long-term observations
The groundwater level data recorded between 1997 and 2015 provide the opportunity to study long-term co-seismic changes at single wells. Figure 8 shows changes in groundwater    (Table 3) level recorded at LJ3 well in response to six earthquakes with magnitude greater than 6.0 between 2005 and 2015. The well is located in southwestern Taiwan (Fig. 1), tapping a confined sand aquifer (Fig. 2). The hydraulic conductivity of the aquifer is 4.45 × 10 −5 m/s. 1-Hz data provide a better temporal resolution for detecting the process of sustained groundwater level changes. Table 4 lists some features of these sustained coseismic changes at LJ3. It is noted that all of the sustained co-seismic changes are water level falls, regardless of the focal mechanism and the location of earthquake epicenter. These sustained co-seismic falls revealed that crustal extension prevailed in the area near LJ3 during earthquakes. Such a persistent extension is not consistent with the concept of a simple dislocation model for the earthquakes (Okada 1992). The differences are attributed to the subsurface geologic conditions being far more complex than those simple conditions assumed in the dislocation model, particularly in the mountainous area.
While all of co-seismic groundwater level readings fell at LJ3 well, most of the co-seismic changes rose at PD well ( Fig. 1). The PD well is located at the top of Douliu Hills, with an elevation of 297.7 m asl. The well, screened between the depths of 144 m and 198 m, is tapping a gravel-confined aquifer where the hydraulic conductivity is 2.83 × 10 −7 m/s (Fig. 2). Figure 9 shows co-seismic changes in groundwater level at PD during 14 earthquakes with magnitude greater than 6.0 between 1997 and 2013. Of those, the 2006 M W 7.0 and M W 6.9 offshore Hengchun earthquakes are doublet, which occurred at the same hour, and hence they are counted as one event. The recorded sustained changes in response to these  The largest co-seismic fall of 223 cm at PD was induced by the 1999 M W 7.6 Chi-Chi earthquake (Fig. 9). It is noted that the co-seismic fall was followed by a rapid rise, instead of a slow recovery in normal cases, immediately after the earthquake. Within 24 h after the fall, the hydraulic head rose to a level approximately 80 cm above that before the shock. As the seismogenic fault of the Chi-Chi earthquake is located only 6.5 km east of the PD well, the fault rupture is likely to trigger an extensive strain in the adjacent area (Chia et al. 2001). The abrupt co-seismic water level drop of 223 cm was attributed to the elastic response to the crustal extension during the earthquake; however, the rise of 303 cm in 24 h after the earthquake was not likely caused by the slow groundwater flow. Such a large and rapid change is possibly associated with the postseismic crustal deformation in response to the rupture of nearby Chelungpu Fault, although more information is needed to understand this phenomenon. Another co-seismic fall of 14 cm was induced by its aftershock of M W 6.5. Except for these two falls, all the other 11 co-seismic changes are sustained rises, ranging from 2 to 216 cm. These sustained rises suggested that, except for the nearby fault rupture, compressive strains dominated the area near the PD well during the earthquakes. Besides, most of the co-seismic changes took one to 4 h to reach a steady level, suggesting that the low hydraulic conductivity retarded the groundwater flow to reach a steady rate at PD well.
Due to the convergence of the Eurasian plate and the Philippine Sea plate, many large earthquakes occur in the vicinity of the Longitudinal Valley of eastern Taiwan where the HL well is located. One would expect to observe large and frequent sustained co-seismic groundwater level changes at HL. While oscillatory changes were often observed (Fig. 4), only three sustained co-seismic changes were recorded during  Fig. 7 Co-seismic and post-seismic groundwater level changes associated with the 1999 M W 7.6 earthquake at three monitored wells. a A 3.72-m co-seismic rise was recorded at HH3 well tapping a confined aquifer. The recovery of the sustained change lasted beyond four months. b A coseismic fall of 2.47 m was recorded at TC2 well tapping an unconfined aquifer, which was followed by a rise of 2.23 m in 12 h. c A co-seismic fall of 5.94 m was recorded at LY2 well, which was followed by a continuous decline of water level    (Table 4) earthquakes of magnitude greater than 6.0 between 2005 and 2015 (Fig. 10). These changes, ranging from 0.2 to 2 cm, were triggered by three local earthquakes with epicentral distances of 26, 63, and 51 km, respectively. In fact, during large earthquakes, small or negligible sustained changes were observed at all monitoring wells in the Longitudinal Valley.  Fig. 9 Co-seismic changes in groundwater level at PD well in response to 14 large earthquakes from 1999 to 2013. Changes are observed on hourly data over a period of 2 days centered at the time of earthquakes. Vertical scale is in the unit of 50 cm for the 1999 M W 7.6 earthquake and the 2013 M W 6.0 and 2013a M W 6.3 earthquake, and 10 cm for the other earthquakes. Co-seismic rises were observed during these earthquakes, except the 1999 M W 7.6 earthquake and one of its aftershocks Numerical simulation

Conceptual model
The variation of groundwater level in the monitoring well reveals the changes of pore pressure in the aquifer. Pore pressure in the aquifer may change in response to many natural factors, such as tide and barometric pressure. Without any artificial factors, the change in pore pressure is usually a slow process; however, during the occurrence of an earthquake, the stress redistribution induced by the fault displacement may cause abrupt changes in crustal displacement. A co-seismic compressive displacement may cause a sudden rise in pore pressure of the aquifer, while an extensive displacement may induce an abrupt fall. Generally, the pore pressure change of the aquifer induced by the crustal displacement will recover after the earthquake, but the recovery rate varies with local hydrogeological conditions.
In this study, three two-dimensional (2D) conceptual geological models were established to simulate pore pressure changes in response to the displacement during earthquakes. Model A was developed for an unconfined aquifer, including both homogeneous and heterogeneous conditions, while model B was established for a confined aquifer and model C was set up for an artesian aquifer fractured during the earthquake. The dimension of these models is 3,000 m wide and 300 m deep, presenting the situation of the shallow crust. The physical properties of the three models are listed in Table 5. Typical values for gravel, sand and mud are selected referring to Das (2008), Domenico and Schwartz (1997) and Obrzud and Truty (2012). The compressibility is calculated from the formula: where v is Poisson's ratio and E is Young's modulus. The Brock^property reflects the semi-consolidated sedimentary rock and a newly formed fracture that did not contain gouge is considered for the Bfracture^material.
There is no vertical displacement along the bottom and no horizontal displacement on the left and right boundaries of the model. The excess pore pressure is set to zero (hydrostatic) along the upper boundary, and no flow conditions are assumed along the lateral and lower boundaries. The porous medium of the geological model is assumed to be water saturated. The pore pressure is considered to be hydrostatic initially for models A and B, whereas for the model C, the pore pressure in the confined aquifer is 8 m above the hydrostatic pressure initially.

Theoretical development
Earthquake induced groundwater level change can be described by a coupled process of soil deformation and pore pressure diffusion. Terzaghi's theory of one-dimensional (1D) soil consolidation (Terzaghi 1926), which states that the total stress is equal to the sum of effective stress and pore water pressure, provided the basis for the coupled process. Biot (1941Biot ( , 1956 developed the general threedimensional concept of poroelastic theory. Rice and Cleary (1976) presented the coupled process of elastic deformation and pore fluid diffusion in compressible porous media. The pore pressure change in a deforming porous medium, which is primarily controlled by pore fluid flow and loading, may be expressed as (Palciauskas and Domenico 1989) Fig. 10 Co-seismic changes in groundwater level at HL well in response to three large earthquakes from 2009 to 2013. Based on 1-Hz data, both sustained rise and sustained fall were recorded. Compared with LJ3 and PD, the magnitude of co-seismic changes at HL is much smaller where P ex is excess pore pressure, t is time, K is hydraulic conductivity, ρ w is density of pore water, g is gravitational constant, σ is stress, 1/R is the change in fluid mass with pore pressure at constant stress, and R/H is the change in pore pressure with stress at constant fluid mass. Eq. (1) can be used to describe the excess pore pressure change under both drained and undrained conditions. The soil deformation during earthquakes is dominated by a sudden load due to the redistribution of stress field induced by fault displacement. As the duration of stress redistribution and co-seismic change in pore pressure is usually very short, the influence of pore fluid flow can be neglected. Hence, Eq. (1) can be simplified to As the stress redistribution ceased after earthquakes, the influence of loading can be neglected. Thus, the change in excess pore pressure is primarily controlled by the pore fluid flow and Eq. (1) can be simplified to Apparently the post-seismic change in excess pore pressure is related to the hydraulic conductivity and specific storage of the formation.

Numerical models
A 2D numerical model using the finite element software ABAQUS was developed to simulate time-dependent pore pressure change in the aquifer due to fault displacement during earthquakes. The coupled pore fluid diffusion and stress analysis capability of ABAQUS can provide pore pressure changes in the porous media induced by the displacement. Three transient models, models A, B, and C, were developed to simulate earthquake-triggered groundwater level changes in the unconfined and the confined aquifers under various conditions. The CPE8RP element, an eight-noded linear strain quadrilateral element with eight displacement nodes and four pore pressure nodes, was used to perform the analysis. A three-step method was adopted for the simulation. The first is the geostatic step, which calculates the initial geostatic stress in equilibrium with initial conditions and boundary conditions. For models A and B, the second step simulates the generation of excess pressure in response to the compressive displacement during the earthquake, whereas for model C, the second step simulates the change of pore pressure in response to the extensive displacement. The third step simulates pore pressure changes due to groundwater flow after the earthquake.
In this study, model A1 is used to simulate the pore pressure change in a 300-m thick homogeneous unconfined gravel aquifer in response to the compressive displacement during earthquakes. The vertical profile of calculated excess pore pressure in the aquifer after the compressive displacement is shown in Fig. 11a. Here the instantaneous excess pressure head is about 30 cm everywhere. Evidently, the excess pressure is the elastic response of the porous medium to the displacement. The dissipation of excess pressure after the displacement occurs due to the upward pore water flow driven by the hydraulic gradient between the water table and the porous medium. Generally the dissipation rate of excess pressure is faster at the shallow depth, but it becomes slower at the deep depth. It takes approximately 250 s to reduce 90% of the excess pressure at 20 m depth, but 1,650 s at 120 m depth (Fig. 11c). If the aquifer thickness is reduced to 30 m, it takes only 22 s to dissipate 90% of the excess pressure. If the hydraulic conductivity of the aquifer is reduced from 10 −2 m/s to 10 −3 m/s, the dissipation time will increase one order approximately.
If the unconfined aquifer is heterogeneous, consisting of alternating 50-m thick layers of gravel and sand (model A2), the generated excess pressure head is 29.9 to 33.5 cm in the gravel layer and 20.2 to 23.7 cm in the sand layer due to the stress concentration in the less compressible gravel layer (Fig. 11b). It takes 63 s to reduce 90% of the excess pressure at 20 m depth in the shallowest gravel layer, but 7,300 s at 120 m depth (Fig. 11c). The different results between models A1 and A2 suggest that the excess pressure may dissipate faster if the formation is underlain by a less permeable layer, but slower if it is overlain by a less permeable layer. In both homogeneous and heterogeneous conditions, 99% excess pressure is dissipated in less than 1 day, indicating the strong influence of the atmospheric pressure condition at the water table of the unconfined aquifer.
Model B is used to simulate the excess pressure change in a confined sand aquifer, which is overlain by a mud layer and an unconfined sand aquifer, in response to compressive displacement. The excess pressure heads generated by the compressive displacement at 50, 150 and 250 m are 19.3, 3.5 and 19.1 cm, respectively, at the beginning. Apparently, the higher excess pressure generated in the aquifer results from the stress concentration in the less compressible sand formations. The dissipation rate of excess pressure in the confined aquifer is somewhat different from that in the unconfined aquifer (Figs. 11 and 12). The rapid recovery process of pore pressure in the unconfined aquifer is similar to those of models A1 and A2. In the confining layer and the confined aquifer, however, the recovery rate is much slower (Fig. 12a). As no flow condition is assumed along the lateral boundaries of model B, it takes more than 1 year to dissipate the excess pressure completely. The dissipation process in the unconfined aquifer, the confining layer, and the confined aquifer can be shown by the changes of excess pressure at the depths of 50, 150, and 250 m, respectively, in Fig. 12b. While the dissipation of the excess pressure at 50 m depth in the unconfined aquifer took less than 1 day, the dissipation of 50% excess pressure in the confined aquifer took more than 8 days. In the middle of the confining layer (150 m depth), the excess pressure declined during the first 2 days and then increased in the following 13 days before a longterm decline. The increase is attributed to the upward dissipation of excess pressure of the underlying confined aquifer. The simulated results suggest that the dissipation of excess pressure in the confined aquifer generated during earthquakes is a very slow process due to the retardation of the less permeable confining layer. Therefore, co-seismic groundwater level changes in the confined aquifer are likely to sustain for a long time.
Streamflow increases, due to fracturing of rock aquifers in the hilly area induced by earthquakes, have often been reported (Rojstaczer and Wolf 1992;Quilty et al. 1995;). Model C is used to simulate the pore pressure change in an artesian aquifer in response to the extensive displacement and rock fracturing during earthquakes. The model is composed of a 50-m thick unconfined sand aquifer, a 50-m thick confining mud layer, and a 200-m thick artesian aquifer where the hydraulic head is 8 m above the hydrostatic level. In addition, there is a 50-cm wide 300-m deep vertical fracture created by the earthquake. Figure 13a shows the vertical distribution of excess pressures in model C due to extensive displacement during earthquakes. The transient response of excess pressure at the depths of 75 and 150 m is shown in Fig. 13b. The pressure head in the artesian aquifer (at 150 m depth) fell 300.6 cm abruptly in response to the extensive displacement during the earthquake. This is followed by a slow decline, instead of a slow recovery. After 1 month, the excess pore pressure reduced to the hydrostatic level, and the confined aquifer was no longer artesian. The simulation results suggest that, after a strong earthquake, the state of an artesian aquifer becomes hydrostatic due to groundwater discharge through a highly permeable conduit, such as a fracture zone, induced by seismic shaking. The head in the confining layer between the unconfined aquifer and the confined aquifer dropped slightly during the extensive displacement, and declined slowly afterwards.

Discussion
The simulation of the pore-pressure change process due to fault displacement may provide proper explanation of groundwater level changes during and after earthquakes. The results of models A1 and A2 show a rapid recovery of co-seismic change in groundwater level in the unconfined aquifer (Fig. 11). The recovery rate is related to the hydraulic conductivity and the thickness of the aquifer. The thinner the aquifer, the faster the recovery of co-seismic change. Such a rapid recovery is evidenced by the field observations at the SL well (Fig. 6). In a thick unconfined aquifer, the calculated coseismic change at depths deeper than 100 m, particularly in a heterogeneous aquifer, may take several hours to complete the recovery process. Similarly, at the 270-m deep TC2 well tapping an unconfined aquifer where the hydraulic conductivity is 1.24 × 10 −3 m/s, it took approximately 12 h to complete the recovery process (Fig. 7b).
In model B, the excess pressure generated by the compressive displacement during earthquakes varies with the physical properties of formations. The simulation results indicate that the recovery of the excess pore pressure in the confined aquifer may take several months or longer. Such a slow recovery rate results in a Bsustained^change. In the field observation, more than 4 months for the recovery of co-seismic change of 3.72 m was observed in HH3  Fig. 13 Dissipation of excess pressure in an artesian aquifer overlain by a confining layer and unconfined aquifer in response to the extensive displacement during the earthquake. a Spatial and temporal changes in excess pressure in three formations fractured during the earthquake. b Temporal changes in excess pressure at 75 and 150 m depth well after the 1999 M W 7.6 earthquake (Fig. 7a). HH3 was installed in a confined aquifer where the hydraulic conductivity is 2.92 × 10 −4 m/s. According to the simulation results of model B, the recovery rate of co-seismic changes depends primarily upon the hydraulic conductivity of the confining layer rather than that of the confined aquifer itself.
After extensive displacement and fracturing induced by the earthquake in model C, the artesian aquifer becomes hydrostatic in 30 days as a result of fracture flow (Fig. 13). The simulation results are similar to the observed water-level changes at the artesian well LY2 triggered by the 1999 M W 7.6 earthquake (Fig. 7c). After the co-seismic drop at LY2, a continuous decline, instead of a recovery, of groundwater level was observed. The phenomenon suggested that after the artesian aquifer was fractured by the seismic shaking, the groundwater outflow through fractures would further reduce the pore pressure. Consequently the aquifer can no longer return to the original state. The earthquake-triggered fractures may create conduits through the confining layer to dissipate the pore pressure of the aquifer. When the fracture is connected to the ground surface such as at riverbeds or creek valleys, the fracture flow may result in streamflow increase immediately after the earthquake. The phenomenon has been observed after many strong earthquakes (Rojstaczer and Wolf 1992;Wang et al. 2004b;Cucci and Tertulliani 2015).
The simulation results suggest that the confining layer plays a major role in governing the recovery rate of coseismic groundwater level changes after earthquakes. The field observations indicated that while the recovery process of co-seismic changes takes minutes to several hours in the unconfined aquifer, it usually lasts for weeks to several months in the confined aquifer. Such a different pattern of water-level changes induced by large earthquakes has the potential to be used as a criterion to evaluate whether an aquifer is confined or unconfined. As the aquifer at a few kilometers depth is usually overlain by several confining layers, the dissipation of co-seismic rises is considerably slow. The excess pressure may build up over a long period of time as a result of frequent crustal compression during large earthquakes. Therefore, in the orogenic belt under compression, an abnormally high pore fluid pressure is likely to occur in the deep formation.
The polarity and magnitude of co-seismic changes was thought to be affected by the focal mechanism and the epicentral distance of earthquakes (Roeloffs 1998;Jónsson et al. 2003;Montgomery and Manga 2003). Modeling of volumetric strain due to fault movement during earthquakes has been conducted in previous studies Rojstaczer et al. 1995;Koizumi et al. 1996;Quilty and Roeloffs 1997;Grecksch et al. 1999;Lee et al. 2002;Matsumoto et al. 2003;Koizumi et al. 2004;Itaba et al. 2008;Shi et al. 2014;Cucci and Tertulliani 2015). The polarities of some calculated strains are moderately consistent with the observed changes; however, the magnitude of most calculated strains did not match the observed changes Grecksch et al. 1999;Matsumoto et al. 2003;Koizumi et al. 2004;Itaba et al. 2008;Shi et al. 2014). Long-term data of groundwater level in the near-field revealed that only co-seismic falls were recorded at LJ3 and 11 out of 13 changes were co-seismic rises at PD. The polarity of co-seismic changes at LJ3 and PD is possibly affected more by local geological conditions and less by focal mechanism-for instance, as the PD well is located at the axis of an anticline, the crust at the well-site is likely under compression during most earthquakes, resulting in the coseismic rises in groundwater level. Long-term data also indicated that the co-seismic changes at HL were consistently small, even during the earthquakes occurring nearby. According to the simulation results, the stress distribution and the co-seismic pore pressure change during earthquakes were affected by the formation compressibility. Such a phenomenon suggests that the stress change is small in the sediments of the Longitudinal Valley during earthquakes, while the stress change is likely to concentrate in the rock formation of the nearby hills.

Conclusions
The study indicated that field observations of groundwater level changes during the earthquake are consistent with the simulation results of pore pressure changes induced by the fault movement. Co-seismic groundwater level in the unconfined aquifer can quickly return to its original level, but at the depth of a few hundred meters, the recovery process may take several hours. In a confined aquifer, however, the co-seismic change may sustain for several months due to the confinement of less permeable layers. Thus, a co-seismic groundwaterlevel change has the potential to be used as a geophysical indicator of whether the aquifer is confined or not.
In the confined aquifer, the dissipation of excess pore pressure after the earthquake is slow due to the less permeable confining layers. The thickness and the hydraulic conductivity of the confining layers play an important role in determining the dissipation of the layers under it. The accumulation of sustained co-seismic rises due to historical large earthquakes may lead to abnormally high pore pressure in deep formations in the orogenic belt. Additionally, fracturing of the confining layer by seismic shaking in the hilly area may change a confined artesian aquifer to an unconfined aquifer. Fracturing may create a highly permeable conduit and allow the excess pore pressure to return to hydrostatic condition in a short period of time. The dewatering of such an aquifer may result in a prolonged increase of streamflow in adjacent streams.
The analysis of long-term data revealed that the rise or fall of co-seismic groundwater level is associated primarily with the tectonic conditions of the well site. The polarity and magnitude of co-seismic changes are related to the focal mechanism and epicentral distance as well as local geological setting. The difference in compressibility between the unconsolidated deposits and rock formations affects the stress distribution and pore pressure changes during earthquakes.