Pure appl. geophys. 165 (2008) 373–407
0033–4553/08/020373–35
DOI 10.1007/s00024-008-0309-9

Ó Birkhäuser Verlag, Basel, 2008

Pure and Applied Geophysics

Estimation of Nonlinear Soil Behavior During the 1999 Chi-Chi,
Taiwan, Earthquake
O.V. PAVLENKO1 and K.-L. WEN2

Abstract—The 1999 Chi-Chi, Taiwan, earthquake (Mw = 7.6) was one of the strongest earthquakes in recent
years recorded by a large number of strong-motion devices. Though only surface records are available, the obtained
strong-motion database indicates the variety of ground responses in the near-fault zones. In this study,
accelerograms of the Chi-Chi earthquake were simulated at rock and soil sites, and models of soil behavior wereconstructed at seven soil sites (TCU065, TCU072, TCU138, CHY026, CHY104, CHY074, and CHY015), forwhich parameters of the soil profiles are known down to depths of at least *70 m and at 24 other soil sites, forwhich parameters of the soil profiles are known down to 30–40 m; all the sites were located within *50 km fromthe fault. For reconstructing stresses and strains in the soil layers, we used a method similar to that developed for theestimation of soil behavior based on vertical array records. As input for the soil layers, acceleration time histories
simulated by stochastic finite-fault modelling with a prescribed slip distribution over the fault plane were taken. In
spite of the largeness of the earthquake’s magnitude and the proximity of the studied soil sites to the fault plane, thesoil behavior at these sites was relatively simple, i.e., a fairly good agreement between the spectra of the observed
and simulated accelerograms and between their waveforms was obtained even in cases where a single stress-strain
relation was used to describe the behavior of whole soil thickness down to *70–80 m during strong motion.
Obviously, this is due to homogeneity in the characteristics of soil layers in depth. At all the studied sites, resonant
phenomena in soil layers (down to *40–60 m) and nonlinearity of soil response were the main factors defining soil
behavior. At TCU065, TCU110, TCU115, CHY101, CHY036, and CHY039 liquefaction phenomena occurred in
the upper soil layers, estimated strains achieved *0.6–0.8%; at other stations, maximum strains in the soil layers
were as high as 0.1–0.4%, according to our estimates. Thus, valuable data on the in situ soil behavior during the
Chi-Chi earthquake was obtained. Similarity in the behavior of similar soils during the 1995 Kobe, 2000 Tottori(Japan), and Chi-Chi (Taiwan) earthquakes was found, indicating the possibility of forecasting soil behavior in
future earthquakes. In the near-fault zones of the three earthquakes, ‘‘hard-type’’ soil behavior and resonant
phenomena in the upper surface layers prevail, both leading to high acceleration amplitudes on the surface.
Key words: 1999 Chi-Chi earthquake, Nonlinear soil behavior, Liquefaction, Stress-strain relations.

1. Introduction
A destructive earthquake of magnitude Mw = 7.6 and depth of 7.5 km occurred incentral Taiwan with its epicenter near the town of Chi-Chi on September 21, 1999. The
1
Institute of Physics of the Earth, Russian Academy of Sciences, B. Gruzinskaya 10, Moscow 123995,
Russia. E-mail: olga@ifz.ru
2
National Center for Research in Earthquake Engineering, 200, Sec. 3, Hsin-hai Rd, Taipei 106, Taiwan,
R.O.C. E-mail: wenkl@earth.ncu.edu.tw

374

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,earthquake had a predominantly thrust focal mechanism, and the surface rupture wasmapped over a distance of *100 km in a SN direction. Field observations show thegreatest surface offset to be more than 8 meters vertically, occurring near the northernend of the fault. The earthquake triggered more than 400 strong motion devices acrossTaiwan. The strong-motion wavefield was captured by a dense network of stations, andthe obtained records represent one of the most complete strong-motion databases to date.
The majority of near-fault records were obtained at soil sites (virtually all stations on the
footwall of the fault are soil stations); these records show signs of nonlinearity in soil
response, such as, predominating low-frequency components in oscillations, smoothed
spectral shapes close to the form of E(f) * f -n in cases of strong nonlinearity (PAVLENKO
and IRIKURA, 2005), and characteristic spiky waveforms (ARCHULETA, 1998).
The possibilities of studying soil behavior during the 1999 Chi-Chi earthquake are
restricted by the fact that only surface records of the earthquake are available. A single
vertical array at Dahan, located on dense soils, did record this earthquake but at a rather
long distance (*80 km) from the fault plane (CHIU, 2001; PAVLENKO and LOH, 2003).
Nonlinearity of soil response during the Chi-Chi earthquake has been previously studied
by ROUMELIOTI and BERESNEV (2003) and by K.-L. Wen (Wen, pers. communic.). Since
only surface records of the Chi-Chi earthquake are available, in these studies, spectra ofacceleration time histories at soil sites are investigated in order to find characteristics ofsoil nonlinearity. Lower amplification of seismic waves in some frequency bands duringstrong ground motions compared to weaker ones (aftershocks) is interpreted as anindicator of nonlinearity in soil response, and respective quantitative characteristics are
introduced.
In ROUMELIOTI and BERESNEV (2003), nonlinearity of soil response during the Chi-Chi
earthquake is studied by comparing amplification of seismic waves at soil sites during themain shock and aftershocks of the Chi-Chi earthquake. The authors simulate records of
the Chi-Chi earthquake by the stochastic finite-fault method and calibrate these against
the data recorded at twenty-four rock sites. Then soil-site records are simulated using thelinear-response assumption whereby the simulated soil-site input motions are amplified
by weak-motion amplification functions, estimated by the spectral ratio technique fromthe aftershock records. Their comparison with actual observations allows the authors toreveal an average reduction in strong-motion amplification to about 0.5–0.6 of that inweak motions, with an acceleration ‘‘threshold’’ for detectable nonlinearity near 200–
300 cm/s2. The conclusions of this study confirm in the whole earlier obtained average
estimates of decreasing amplification of seismic waves during strong ground motion at
soil sites; however, the authors emphasize a large interstation response variability, which
did not allow them to derive a statistically significant difference in weak and strongmotion amplifications based on the responses available at the studied sixteen stations.
The present study is aimed at a more detailed investigation of soil response during the1999 Chi-Chi earthquake. Accounting for the observed variability of soil response
patterns, we concentrate on constructing models of soil behavior at separate stations,instead of studying average regularities, which seems to be more informative in this case.

Vol. 165, 2008

Nonlinear Soil Behavior

375

To do this, we simulate acceleration time histories at soil sites at the bottoms of soil
layers and then calculate their propagation up to the surface, selecting models of soil
behavior that show best fit to the observed acceleration records at the surface. Models ofsoil behavior during the Chi-Chi earthquake are constructed for seven soil sites, for whichparameters of the soil profiles were known down to depths of at least 70 m, and fortwenty-four soil stations, for which the profiling data were available down to depths of*30–40 m. Thus, models of soil behavior during the Chi-Chi earthquake are constructed
for thirty-one soil sites located within *50 km of the fault plane. In constructing models
of soil behavior during the Chi-Chi earthquake, we use our previous experience instudying the soil behavior during the 1995 Kobe and 2000 Tottori (Japan) earthquakes
based on vertical array records (PAVLENKO and IRIKURA, 2003, 2005, 2006). The
constructed models of soil behavior testify to the strong nonlinearity of the soil responsein the near-fault zones of the Chi-Chi earthquake.

2. Simulation of Acceleration Time Histories of the Chi-Chi Earthquake at Rock Sites
The epicenter of the 1999 Chi-Chi earthquake, the surface rupture and locations ofeighteen rock and thirty one soil stations studied in this work are shown in Figure 1. We
used strong-motion data of the Chi-Chi earthquake disseminated on CD-ROM by LEE
et al. (2001b), accelerograms provided by the website of the Central Weather Bureau(CWB) of Taiwan and profiling data provided by the National Center for Research in
Earthquake Engineering (NCREE) of Taiwan.
Firstly, we calculated accelerograms of the Chi-Chi earthquake at rock sites. Since thepurpose of our study is to simulate soil behavior, we tried to simulate most accurately
motion at the ‘‘input’’ of the soil layers. In the calculations, we used Boore’s stochastic
approach (BOORE, 2003); however, we modified his method to account for the dimensionsof the source and slip distribution over the fault plane.
The availability of such a complete near-source strong-motion dataset being
generated by the Chi-Chi earthquake, combined with broadband teleseismic displacement
waveforms, well-distributed Global Positioning System (GPS) data, and the considerable
volume of geological and geophysical data collected in Taiwan over the past decades has
provided an excellent opportunity to study the kinematic source process of this large
earthquake, and various versions of slip distribution over the fault plane have beenobtained (MA et al., 2001; ZENG and CHEN, 2001; OGLESBY and DAY, 2001; WU et al.,
2001; CHI et al., 2001; WANG et al., 2001). In our calculations we used the slip
distribution obtained by W.-C. Chi, D. Dreger and A. Kaverina (CHI et al., 2001).
Velocity waveforms recorded at 21 stations were inverted for spatial variation in slip on aplanar fault model composed of 416 subfaults, each with a dimension of 3.5 km by
3.5 km. To account for possible temporal source complexity, each subfault was allowed
to slip within 10 overlapping 3-second time windows. The kinematic rupture process wasstudied in detail in this work (CHI et al., 2001), and various inversion problems were376

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 1
Map showing the epicenter of the 1999 Chi-Chi earthquake (star), the surface rupture and locations of eighteenrock (black triangles) and thirty one soil stations (gray triangles in circles) used in this work. Dashed lines boundthe areas where Rc < 15 km and Rc < 50 km. Dotted lines mark liquefaction areas, and the area, where basininduced surface waves were identified, is marked by grey color.

discussed and justified, with error analysis being performed; the results of the inversion
are provided on CD-ROM (LEE et al., 2001b) and could be easily inserted in thecalculations.
In simulating accelerograms of the Chi-Chi earthquake, we represented its source(119 km by 35 km) as a set of subfaults with dimensions, locations, parameters ofradiation, as well as the method for summarizing their energy being selected inaccordance with recommendations by GUSEV and PAVLOV (2006), who studied all aspectsof earthquake ground motion simulation for many years and achieved success in this area
(BAZZURRO et al., 2001) (we represented the space-time structure of the slip rate function
through a grid of subsources, where the cell size of a subsources was defined by site-tofault distance; onset times and seismic moments of subsources were determined by theVol. 165, 2008

Nonlinear Soil Behavior

377

generalized Haskell source model; time histories of subsources were taken uncorrelated;
final slip and rupture velocity were defined as random functions and varied along x and y;each point subsource was randomly shifted from a node of the perfect regular grid(typically by 0.3–0.6 cell size)). Each of the subfaults was assigned an x2 spectrum.
Some other parameters used in the calculations, such as, spatial orientation of the fault
plane, the nucleation point, coordinates and depth of the hypocenter, slip velocity and slip
distribution over the fault plane, were derived from the paper by CHI et al. (2001). We
selected the subfault dimension of 7 km by 7 km, averaging slips over the adjacent areasof 3.5 km by 3.5 km of the ‘‘preferred model’’ by CHI et al. (2001); the estimatedtemporal variations of the rupture process were also averaged, and the final simplified andsmoothed representation of the rupture process was obtained (Fig. 2). This representation
was used for generating two spatial weight functions, ‘‘dip-slip’’ and ‘‘left-lateral
strike-slip’’, indicating the relative contributions of the subfaults in the two directions.
Time delays for subfault radiation were calculated based on rupture velocity, which wastaken as constant and equal to 2.6 km/s according to (CHI et al., 2001). Subfault radiation
was simulated as stochastic, empirically attenuated to the observation site, andsummarized with introducing slight random variations in subfault locations and rupture
velocity. Local rise times were taken as equal for all subfaults. Thus, acceleration time
histories in ‘‘dip-slip’’ and ‘‘left-lateral strike-slip’’ directions were calculated with spatial
weight functions (Fig. 2) applied to subfault radiation and rotated in NS and EW

Dip-slipSouth

Distance along dip, km

-28

-14

0

14

28

North
42

56

70

-7

Slip, cm

7

4800

21

4000

Distance along strike, km

3200

Left-lateral strike-slip
-28

-14

0

14

28

42

56

70

2400
1600

-7

Km

800

7
0

21

Km

Figure 2
Distribution of slip over the fault plane of the Chi-Chi earthquake used in our simulations: Averaged slip
distribution of the ‘‘preferred model’’ by CHI et al. (2001). Star indicates the hypocenter.

378

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Table 1
Locations of rock stations used in simulations
Station code

Latitude (°)

Longitude (°)

North (forward direction):
TCU04924.179120.690
TCU05324.194120.669
TCU04624.468120.854
TCU02524.707121.176
ILA04324.629121.735
ILA05224.609121.849
East and West (normal to the fault):
TCU07523.983120.678
TCU07623.908120.676
TCU08923.904120.857
HWA02324.080121.596
HWA02423.352
121.297
CHY07523.567119.555South (backward direction):
CHY10223.246
120.614
CHY11023.252120.530
KAU05023.163120.757
KAU04723.082120.583
TTN02422.973121.108
KAU07722.747120.723

Elevation
(m)

Rc
(km)

Distance to the closest
subfault (km)Weighted aver. dist.
to subfaults (km)

129
133
231
310
36
20

10
15
50
50
50
50

6.5
5.4
13.3
36.8
88.6
99.9

28.3
28.9
44.4
80.3
118.9
128.3

104
99
705
30
260
11

15
10
12
50
50
50

6.8
7.5
16.2
79.2
57.7
100.9

31.3
34.8
38.1
93.5
103.2
135.9

560
215
640
278
277
835

50
50
50
50
50
50

33.2
32.7
41.9
50.6
70.1
85.6

95.0
95.8
103.2
112.6
129.7
147.2

directions. Simulated acceleration time histories at rock sites were compared to theobserved ones, and calibration of the simulation procedure was performed, i.e.,
parameters characterizing source and path effects were found, which show the best fit
to the observations.
We calculated acceleration time histories at all near-fault rock sites, for which goodquality acceleration records of the Chi-Chi earthquake at two horizontal components areavailable. We found 18 such sites; their coordinates, elevations, characteristic distances
Rc (discussed below), the distances to the closest subfaults, and weighted (weights being
proportional to the subfault radiation intensity) average distances to subfaults are given in
Table 1. Some stations classified in (LEE et al., 2001a) as Class B Sites, have shown
amplification at low frequencies (the problem is discussed in (LEE et al., 2001a)), which
we interpreted as site effects. We decided to exclude such stations (CHY042, CHY052,
CHY074, CHY086, CHY109, KAU069, KAU078, ILA063, HWA022, HWA026,
HWA046, HWA056, TTN018, TTN026, TTN041, TAP034, and TAP035) from the listof ‘‘rock sites’’ in this study. At the same time, other stations were added to this list basedon the results of simulations, such as, TCU049, TCH053, TCU075 and TCU076 stations.
Figure 3 shows the results of simulations at rock sites. The observed and simulated(one of the series) accelerograms at two horizontal components and spectra of theobserved and simulated (averaged over the series of simulated accelerograms)

Vol. 165, 2008

Nonlinear Soil Behavior

379

Figure 3
The results of simulations of the acceleration time histories of the Chi-Chi earthquake at rock sites: Observed
and simulated accelerograms and their spectra. After the station codes weighted average distances to thesubfaults and distances to the closest subfaults are shown.

380

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,accelerograms are shown. Comparison of the observed and simulated accelerograms andtheir spectra gives a representation of the accuracy of simulations. As seen from thefigure, a rather good agreement was achieved between the observed and simulated
accelerograms, and some misfit in their spectra (noticeable at stations TCU089, TCU046,
CHY102, and at island station CHY075) is obviously due to local effects of wave
amplification in some frequency bands.
To trace effects of directivity (GUSEV and PAVLOV, 2006), the stations are arranged
according to their orientation with respect to the direction of the crack propagation (CHI
et al., 2001). Directivity effects are clearly seen as decreased durations of acceleration
time histories in forward (northern) directions and increased durations in backward(southern) directions; eastern and western directions represent intermediate cases (Fig. 3).
Table 2 shows parameters used in the simulations. The spatial orientation anddimensions of the fault plane, the coordinates and depth of the hypocenter, the slip
velocity and slip distribution over the fault plane were derived from the paper by CHI
et al. (2001). Frequency-dependent attenuation effects of the propagation path aremodeled through the function Q(f) as described in BOORE (2003). We assumed therelationship estimated for the Taiwan area from coda-waves by CHEN et al. (1989):
Q(f) = 117f 0.77. Kappa operator was used for additional attenuation of the spectra(ANDERSON and HOUGH, 1984; BOORE, 2003); we took j = 0.07 sec. Thus, theseparameters were taken as the same as those used in the work by ROUMELIOTI and BERESNEV(2003), because our calculations show that these values well fit the observed data andonly slightly define the simulations (i.e., their variations within some limits do not greatly
influence the obtained accelerograms).
At the same time our results show that for near-fault stations, geometrical spreading is
one of the most important factors defining the shapes of the acceleration time histories.

Table 2
Parameters used in simulations
Fault orientation (strike/dip)
Fault dimensions along strike and dip (km)Coordinates of the hypocenter (lat., lon.)
Depth of the hypocenter (km)Subfault dimensions (km)Stress drop Dr (bars)
Number of subsources summed
Q(f)Geometrical spreading
Kappa j, (sec)
Crustal amplification

Crustal shear-wave velocity (km/s)
Rupture velocity (km/sec)
Crustal density (g/cm3)

5°/30°
119 by 35
23.8689, 120.84
8.75
797
50
85
117f 0.77
1/r for r < Rc km 1/Rc for Rc B r < 150 km
0.07
BOORE and JOYNER (1997) western North America
generic rock site; at underlined stations – hardrock (BOORE, 2003)
3.2
2.6
2.6

Vol. 165, 2008

Nonlinear Soil Behavior

381

For example, simulating acceleration time histories at TCU075, TCU076, TCU089,
TCU072 and other near-fault stations located within the area marked by a dashed line in
Figure 1 with conventional geometrical spreading relationships:
1=r for r\50 km;

1=50 for 50  r\150 km

ð1Þ

(BOORE, 2003), we obtain too strong attenuation in the medium and final parts of thesimulated accelerograms; this does not satisfy the observations. As seen from Figure 3,
accelerograms at these stations represent successive arrivals of groups of seismic waves
of approximately equal intensity, which evidently correspond to radiation from different
parts of the fault plane (subfaults). Application of conventional relationships (1) produces
accelerograms with the intensity of these waves quickly decreasing with time. However,a good agreement between observations and simulations can be obtained if describing
geometrical spreading by the relationship:
1=r for r\Rc km;

1=Rc for Rc  r\150 km,

ð2Þ

where distance Rc varies from 10 km in the vicinity of the fault plane to 50 km beyond
the area marked by the dash line in Figure 1 (Table 1). For all the stations, Rc values were
found in calculations as showing the best fit to the observations. It is seen from Table 1
that Rc < 50 km in cases where the average removal of the station from the subfaults isless than 44–50 km. As shown below (listed in Table 3), Rc < 50 km also for stations
located to the south of the fault plane at average distances of 50–85 km; the area of
decreased Rc values is stretched in the southern direction (Fig. 1). Apparently, this effectshould be attributed to peculiarities in the geometry of the earthquake source to theeffects of ‘‘subsurface waveguide’’.
At the majority of the studied rock sites, the upper-crustal amplification was describedby parameters established by BOORE and JOYNER (1997) for the western North America
generic rock site. However, at some stations (underlined in Fig. 1 and Table 1) theseparameters gave excessively higher amplification not corresponding to the observations,so they were changed to ‘‘hard rock’’ amplification parameters (BOORE, 2003), which
better satisfy the observations. In simulating acceleration time histories at soil sites, we
did not account for the upper-crustal amplification.
In the whole, a fairly good agreement was achieved between the observed andsimulated acceleration time histories at the studied rock stations, especially in near-fault
zones, in many respects because of the use of slip distribution along the fault plane (CHI
et al., 2001). Some disagreements between observations and simulations seem to berelated to site effects and therefore are not significant within the framework of this study.

3. Simulation of Acceleration Time Histories of the Chi-Chi Earthquake at Soil Sites
At the second step of our study, we simulated acceleration time histories at soil sites
and constructed models of soil behavior during the Chi-Chi earthquake, i.e., estimated382

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Table 3
Information on soil stations used in simulation
Latitude(°)

Longitude(°)

Elevation Rc
(km)
(m)

Closestdist.
to subfault
(km)

Aver.
dist.
to subfaults(km)

Aver. Vs
in
upper
30 m
(m/s)

Site Class(LEE et al., 2001),
soil
composition

24.14824.16124.24924.68624.781120.676
120.675
120.721120.862
121.007

84
99
227
6
86

15
20
22
50
50

5.7
6.1
7.8
16.9
29.7

28.3
28.4
29.6
64.2
78.1

476
454
735
450
600

D
D
D
D
?

23.059120.691

56

12

6.9

28.9

300

TCU072 24.041

120.849

292

12

14.9

32.4

563

TCU110
TCU138
TCU116
TCU115South:
CHY025
CHY092
CHY101
CHY026
CHY002
CHY104
CHY094
CHY03623.96223.92223.85723.960120.570
120.596120.580
120.469

23
37
42
15

14
12
15
35

2.4
3.8
3.4
11.4

35.5
36.2
40.3
41.6

207
605
483
190

D sandy clays, silt,
gravelD colluvium,
breccia, shale,siltstone
E sand, clay, silt
D gravel
E silt, gravel
E silt23.78023.79123.68623.79923.71923.67023.79423.607120.514120.478120.562120.411120.413120.465
120.321120.479

42
34
75
25
26
33
10
45

16
25
10
30
25
30
50
20

5.9
9.5
1.9
16.1
15.0
9.5
24.9
7.5

48.2
49.5
53.0
53.1
58.1
58.7
59.7
63.1

277
260
260
226
225
225
227
282

E
E
D
E
E
E
E
D

CHY082 23.724
CHY027 23.752
CHY074 23.510

120.298
120.247120.805

13
8
2413

50
50
16

26.6
32.1
12.8

65.4
67.6
67.8

210
220
546

E
E
C

CHY032
CHY039
CHY093
CHY004
CHY033
CHY015
CHY044
CHY055
CHY012120.294120.344120.147120.172120.215
120.405
120.164120.271
120.152

9
16
4
6
3
34
2
8
1

50
30
50
35
50
18
50
35
50

25.6
20.8
41.2
38.2
33.6
25.9
43.0
41.5
46.6

75.6
77.3
81.0
82.5
83.7
89.4
99.4
103.4
104.2

202
198
199
279
194
234
193
240
202

E
E
E
E
E
D
E
E
E

Station
Code

North:
TCU082
TCU054
TCU102
TCU033
TCU017
Centre:TCU06523.58023.52123.65423.60123.54123.35523.38323.27023.333gravelsilt, gravel
gravel, siltstone
gravel, siltstone
gravel, mudstonesilt, gravel
silt
silt, sand, gravelsand, clay, silt
silt
sand,clay,silt
silt
silt, clay, sand,gravelsilt
silt
sand, clay, siltstone, breccia
silt, clay, sandsilt, clay, sandsilt, sandsilt, sandsilt, sandsilt, clay, sandsilt, sandsilt, clay, sandsilt, clay, sandstresses and strains induced in the soil layers during strong motion. The soil behavior wasstudied at sites, for which the profiling data were available. At seven soil stations, suchas, TCU065, TCU072, TCU138, CHY026, CHY104, CHY074, and CHY015, soil

Vol. 165, 2008

Nonlinear Soil Behavior

383

composition, P- and S-wave velocities and the results of SPT-N testing are available
down to depths of *70–140 m, and at other twenty four soil stations, these data are
known down to *30–40 m. The coordinates of these thirty one soil stations, theirelevations, distances Rc, distances to the closest subfaults, weighted average distances tothe subfaults, Site Class and soil composition are given in Table 3 (seven stations, forwhich profiling data were available down to *70–80 m, are marked in bold). Site
Classification of Taiwan free-field strong motion stations is described in LEE et al. (2001).
Based on the available geologic and geomorphologic data, free-field strong motion
station sites were classified using a scheme compatible with the 1997 Uniform Building
Code (UBC) provisions. The response spectral shapes (RSS) and the horizontal-tovertical spectral ratio (HVSR) methods and field inspections were used for checking
purposes. However, some problematic sites do exist (LEE et al., 2001).
Concentric lines in Figure 1 mark stations with elevations above sea-level of higher
than 200 m. Though topographic effects were not accounted for in simulating
acceleration time histories, they could contribute to amplification or de-amplificationof oscillations on the surface. The influence of topographic effects on strong motion
stations in Taiwan is discussed in LEE et al. (2001), where a number of stations obviously
influenced by surface topography are listed; none of these stations is studied in thisresearch. Also, as seen from Figure 1 and Table 3, all the studied soil stations, except
TCU072, TCU102, and CHY074, are placed in valleys, and therefore, surface topography
does not influence the records. However, the influence of subsurface topography is seenin the records at some stations and discussed below.
In reconstructing stresses and strains at these stations, our previous experience instudying soil behavior during the 1995 Kobe and the 2000 Tottori earthquakes based onvertical array records was used (PAVLENKO and IRIKURA, 2003, 2004, 2005, 2006). The
method for the estimation of nonlinear stress-strain relations in soil layers based onvertical array records is described in detail in PAVLENKO and IRIKURA (2003).
In this study, to estimate stresses and strains in soil layers during strong motion, we
calculated the propagation of vertically incident ‘‘input’’ seismic waves (simulated
accelerograms) up to the surface. In the calculations, we used the algorithm of nonlinear
analysis by JOYNER and CHEN (1975), modified by Pavlenko so that, instead of thenormalized stress-strain relations used by Joyner and Chen (similar to those obtained inlaboratory experiments and described by HARDIN and DRNEVICH (1972)), any other stressstrain relations of the ‘‘hard’’ (declining to the stress axis at large strains) or ‘‘soft’’(declining to the strain axis at large strains) types can be used to describe the soil
behavior. The relations can differ for different soil layers and at different time intervals.
In our calculations, these stress-strain relations were prescribed parametrically, anddiffered from each other by their shape. For calculations, we have to divide the studied
soil profiles into groups of layers, according to the soil composition, depth and saturation
with water, and assume certain ‘‘hard’’ or ‘‘soft’’ types of stress-strain relations todifferent groups of layers. The choice of ‘‘hard’’ or ‘‘soft’’ types of soil behavior is basedon the established regularities of the behavior of cohesive and non-cohesive soils in situ

384

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

in conditions of dynamic loading. It is discussed in detail in PAVLENKO and IRIKURA (2003,
2004, 2006). ‘‘Hard-type’’ stress-strain curves describe well the behavior of watersaturated sandy soils. Analysis of accelerograms of the 1995 Kobe and 2000 Tottori
earthquakes has shown that ‘‘hard-type’’ soil behavior in strong ground motion is
widespread and typical for sandy soils in cases when the level of underground water isabove *10 m (PAVLENKO and IRIKURA, 2003, 2004, 2006). Whereas soft soils, like dryloess loams and sands, possess stress-strain relations declining to the strain axis at large
strains, which are often called ‘‘soft’’ diagrams.
In LEE et al. (2001a), site classification of free-field strong-motion stations in
Taiwan was performed. Based on available geologic and geomorphologic data, RSS(response spectral shapes) and the HVSR (horizontal-to-vertical spectral ratio) methods
and field checks, the stations were classified using a scheme compatible with the 1997
Uniform Building Code (UBC) provisions. Ten stations, including TCU082, TCU054,
NCU102, TCU065, TCU072, TCU138, CHY015, and others, were categorized as class
D sites, station CHY074 as a class C site, and the rest twenty stations, such as,
TCU110, TCU116, TCU115, CHY025, CHY026, CHY104, etc., as class E sites. Data on
the composition and thickness of the soil layers, P- and S-wave velocities, and SPT-N
values are taken from the National Center for Research in Earthquake Engineering(NCREE) of Taiwan.
First we have studied soil behavior at seven stations possessing the most complete
information on the profiling data, such as, TCU065, TCU072, TCU138, CHY026,
CHY104, CHY074, and CHY015 stations. At TCU065, CHY026, CHY104, and CHY015
stations, soil profiles in the upper *50 m, *100 m, *90 m, and *60 m, respectively,
represent alternating layers of sandy clays with silt and silty soils with S-wave velocities
increasing with depth from *160–200 m/s at the surface to 420–600 m/s at depth(parameters of the soil profiles are shown below in Figs. 4a-g). These layers areunderlayed by denser colluvium soils. At TCU072, TCU138, and CHY074 stations
colluvium deposits are closer to the surface. At TCU072 and CHY074 stations, S-wave
velocities increase with depth from 250 m/s and 360 m/s, respectively, in the upper
8–10 m to *1120 m/s at depths of *70 m. At TCU138 station, the soil layers down to
depths of *150 m represent mostly dense soils with S-wave velocities varying within*400–800 m/s.
According to P-wave velocities, the underground water level is rather close to thesurface at all the studied sites, and therefore, ‘‘hard-type’’ stress-strain relations wereselected to describe the behavior of upper soil layers at all the sites. In this work, todescribe the behavior of the soil layers down to *70–80 m at each station, we usedone ‘‘hard-type’’ normalized stress-strain relation. That is, we treated all the soil
thickness at all the stations, except TCU065, as one group of layers. With these simple
models of soil behavior, we could obtain a good agreement between observations andsimulations.
At TCU065 to fit the observations, we had to apply a more complicated model: We
used one ‘‘hard-type’’ stress-strain relation to describe the behavior of the upper *37 m

Vol. 165, 2008

Nonlinear Soil Behavior

385

Figure 4
The results of simulations of the acceleration time histories of the Chi-Chi earthquake at soil sites: simulated
accelerograms at the bottoms of the soil layers (upper left corner); below — the profiling data (VS – solid line
and smax – dashed line) and spectra of the simulated accelerograms at the bottoms of the soil layers (low-level
grey line) and of the observed (dashed line) and simulated (solid line) surface accelerograms; the observed andsimulated accelerograms on the surface (upper right corner), and the corresponding vertical distributionsof stresses and strains in the soil layers, changing in time during the strong motion for NS and EW components:a –TCU065; b –TCU072; c –TCU138; d –CHY026; e –CHY104, f –CHY074; and g –CHY015. After the station
codes weighted average distances to the subfaults and distances to the closest subfaults are shown.

386

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 4(Contd.)and the stress-strain relation similar to that obtained in laboratory experiments by HARDIN
and DRNEVICH (1972) to describe the behavior of the underlying more dense soil layers.
It was shown in PAVLENKO and IRIKURA (2003, 2004, 2006) that these relations wellVol. 165, 2008

Nonlinear Soil Behavior

387

Figure 4(Contd.)describe the behavior of dense cohesive and non-cohesive soils in situ at depth. The stressand strain were normalized in the manner used by HARDIN and DRNEVICH (1972): Stress
was normalized by multiplying by 1/smax (where smax is shear stress in failure) and strain

388

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 4(Contd.)was normalized by multiplying by Gmax/smax (where Gmax/smax is the low-strain modulus)(JOYNER and CHEN, 1975). The differences in the behavior of different soil layers result
from differences in the values of smax and Gmax assigned to the different layers.

Vol. 165, 2008

Nonlinear Soil Behavior

389

Figure 4(Contd.)Two hundred fifty different curves of the selected type were generated, anditem-by-item examination was applied to identify the curves, showing the best-fit
approximations to the observed accelerograms on the surface. ‘‘Inputs’’ to soil layers, i.e.,

390

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 4(Contd.)simulated accelerograms, were used as prescribed motion at the base of the soil columns(the free surface factor was taken equal to 2, when calculating acceleration time histories
at rock sites, and equal to 1 in this case). That is, we assume infinite rigidity in theVol. 165, 2008

Nonlinear Soil Behavior

391

Figure 4(Contd.)underlying medium. This approximation, that allows no energy to be radiated back into
the underlying medium, is rather accurate in this case, because the energy dissipation
within the soil columns (*80–140 m) is rather high. Multiple reflections and resonances

392

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

occur in the upper tens of meters (down to *40–60 m), as will be shown below, andresult in substantial losses of energy in these layers. Our preliminary calculations havealso shown that the deviations from the case of the imposed motion (when motion at thebase of the soil columns represents the sum of incident and reflected waves, and stressstrain relations in the soil layers are selected to fit observations on the surface and at thebottom of soil) are rather small.
Parameters, which were not known from measurements, but used in calculations, suchas, smax and damping in the soil layers, were estimated based on available information on
the characteristics of soil layers. To achieve a better agreement between the observations
and simulations, sometimes we calculated the soil response to a number of simulated
acceleration time histories within one series, in order to find the ‘‘best-fit’’ simulated
accelerogram.
Analyzing vertical array records, we verified the validity of the assumed stress-strain
relations in the soil layers by calculating deviations of the simulated accelerograms fromthe observed ones at depths of locations of the array devices. In the case of the Chi-Chi
earthquake, only surface records are available, and they can be compared with thecalculated soil response to simulated accelerograms, essentially stochastic signals.
Therefore, to find the best fit to the observations, we compared the spectra, averageintensities (as mean squares of acceleration time histories), and peak ground accelerationsof the observed and simulated accelerograms at two horizontal components. These
criteria were used to calculate a combined mean-square-root value of the deviation of thesimulated spectra, average intensities, and peak ground accelerations from the observed
ones at two horizontal components. This value was used as an estimate of the agreement
between the observed and simulated acceleration time histories, and thus ‘‘the best-fit’’stress-strain relations in the soil layers were selected, which are shown below as models
of soil behavior. The observed and simulated accelerograms were analyzed on the whole,
without dividing them into time intervals, as was done in the case of vertical array record
analysis, because of the evident insufficiency of input information.
The results of simulation of acceleration time histories at soil sites are presented in
Figure 4a-g. For each station, simulated accelerograms at the bottom of soil layers,observed and simulated accelerograms on the surface, and corresponding vertical
distributions of stresses and strains in the soil layers, changing in time during the strong
motion are shown for two horizontal components. The profiling data (S-wave velocities
and shear stress in failure smax in the soil layers) and spectra of the observed andsimulated (averaged over the series of simulated accelerograms) accelerograms are also
shown. As explained above, to describe the behavior of each soil layer, only one stressstrain relation was used for all time intervals during the strong motion; however, in order
to trace changes in soil behavior, working intervals of these stress-strain relations withintwelve 4-second time intervals are presented for each soil layer (Fig. 4a-g).
At stations CHY026 and CHY104, surface accelerograms show the presence of intense
seismic waves with periods of 3–5 seconds in oscillations at the surface (Fig. 4d, e). The
stations are located to the southwest of the fault plane, in a valley, and these long-period

Vol. 165, 2008

Nonlinear Soil Behavior

393

oscillations apparently represent standing waves due to S-wave reverberation in the basin
sediments, i.e., basin-induced surface waves. These waves cannot be modeled in our 1-D
simulations, and therefore, only the beginning part of the strong motion (15–20 sec.) isproperly simulated at these stations. Spectra of the recorded accelerograms in Figs. 4d, ecorrespond to these beginning parts. As seen from the figures, in spite of the substantial
uncertainty of the input information (motion at the bottom of the soil layers is estimated;shear stress in failure, damping in the soil layers, and models of the soil behavior are
suggested), a fairly good agreement between the observed and simulated accelerograms isachieved at the studied stations.
We found that in the case of the Chi-Chi earthquake, even simple models, in which
the behavior of all soil layers is described by the same stress-strain relation, areapplicable to all the studied sites (i.e., show an agreement between observations andsimulations), except TCU065, where the model is more complicated. At the TCU065
station, agreement between observations and simulations can only be obtained whendescribing the behavior of the upper (*37 m) and lower soil layers by different stressstrain relations.
Since the level of the underground water was close to the surface at all the stations(VP > 1300 m/s in surface layers), ‘‘hard-type’’ stress-strain relations were applied todescribe the soil behavior at all the stations. However, the stress-strain curves show
substantial ‘‘hardening’’ (declinations to the stress axis) mostly at stations closest to thefault plane (TCU065 and TCU072), where the levels of the input motion were rather
high. At these stations vertical distributions of stresses and strains in soil layers andspectra of oscillations at the bottoms and on the surfaces of soil profiles indicate resonant
phenomena in the upper soil layers induced by the strong motion, which lead toamplification of seismic waves at *1–2 Hz at TCU065 and at *1–3 Hz at TCU072(Fig. 4a, b). The ‘‘hard type’’ of soil behavior at these stations and rather low S-wave
velocities (*200–250 m/s in the upper *10 m) stipulate a substantial amplification ofoscillations on the surface. At TCU065 station, characteristic spiky waveforms in theacceleration time history are observed, which are often related to the strong nonlinearityof soil response (ARCHULETA, 1998) and to cyclic mobility and liquefaction of sandy
surface soils (LEE et al., 2001a).
As known, three main mechanisms of seismic wave transformation in subsurface soils
are: (1) transition of seismic waves to upper soil layers with smaller values of VS anddensity, q, leading to their amplification according to energy conservation law, (2) resonant
phenomena in the upper soil layers also leading to amplification of seismic waves, (3)nonlinearity of soil response, often leading to de-amplification of seismic waves.
At TCU065 and TCU072 sites (and at many other sites, as shown below), we observe
resonant phenomena in upper soil layers in conditions of strong nonlinearity of soil
response.
Increase of the intensity of seismic oscillations due to resonant effects naturally leads
to increasing nonlinearity of soil response; at the same time, nonlinearity decreases
resonant frequencies of soil layers (if compared to ‘‘linear’’ soil response in weak

394

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

motion), which was noticed by many researchers. Thus, in cases of rather strong
nonlinearity of soil response, resonant frequencies are related not only with geometricaland mechanical characteristics of the structure, but also with the intensity of ‘‘input’’
motion. In conditions of rather strong nonlinearity of soil response, spectra of oscillations
on the surface tend to take the smoothed form E(f)*f -n (PAVLENKO and IRIKURA, 2005).
These transformations of spectra of seismic waves propagating in soil layers should
be studied in special numerical experiments; in this paper, we can only establish facts of
excitation of resonant oscillations in soil layers by the Chi-Chi earthquake. In this case,resonant oscillations (or resonant phenomena) mean trapping of seismic waves due to theimpedance contrast between soft sediments and underlying denser layers. The interferences between these trapped waves lead to resonance patterns. These phenomena areclearly illustrated by spectra and vertical distributions of stresses and strains in soil layers
shown in Figures 4 and 5. Comparison of spectra of oscillations at the bottoms and on thesurfaces of soil layers reveals substantial amplification of their low-frequency components (up to 2–3 Hz, as seen from the figures); whereas vertical distributions of stresses
and strains in soil layers indicate that maximum strains occur in the upper (resonating)layers, i.e., trapping of seismic waves in the upper layers.
At TCU138, amplification of seismic oscillations in the soil layers was not observed(Fig. 4c), because the soil profile represents mostly dense soils, even in the surface layers(VS * 400 m/s). At stations: CHY026, CHY104, CHY074, and CHY015 amplification ofseismic oscillations occurs at low frequencies, as seen from the spectra of oscillations on
the surface and at the bottoms of the soil layers (Fig. 4d-g).
As seen from the figures, resonant oscillations were induced in the upper tens ofmeters (40–60 m) of soil at stations: TCU065, TCU072, CHY026, and CHY104(actually, at all the stations except TCU138) during the earthquake. Maximum strains
occurred at depths of *5–10 m at TCU065 (*0.6%), *8–12 m at TCU072 (0.3–0.4%),*40–100 m at TCU138 (*0.1%), *20 m at CHY026 (0.3–0.4%), *10–15 m at
CHY104 (0.2–0.3%), *40 m at CHY074 (0.1–0.2%), and *30–70 m at CHY015
station (0.1–0.2%).
According to our estimates, at TCU065, resonant phenomena occurred in the upper*9 m of sandy clays with silt (VS * 200 m/s) with maximum strains of *0.6% (LEEet al., 2001a, report liquefaction in the upper layers at this station). At the same time,resonant oscillations were also induced in the upper *37 m (with Vs increasing from*200 m/s to *380 m/s) with maximum strains of *0.2–0.4% (Fig. 4a).
At TCU072, resonances occurred in the upper *11 m of the soft colluvium (Vs *250 m/s), below which more dense breccia and shale layers are deposited (Vs * 550–
700 m/s) (Fig. 4b). At CHY026 and CHY104, resonant oscillations were induced in theupper *20 m (at CHY026, Vs *160–240 m/s) and *12 m (at CHY104, Vs *180–
190 m/s) of sandy clays with silt (Figs. 4d and e). At stations TCU138, CHY074, andCHY015, resonant phenomena were less pronounced and deeper: Maximum strains
occurred at *40–100 m at TCU138 (*0.1%), at *40 m at CHY074 (*0.2%), and at*30–70 m at CHY015 (0.1–0.2%). The corresponding stresses were as high as *0.6–

Vol. 165, 2008

Nonlinear Soil Behavior

395

Figure 5
The results of simulations of acceleration time histories of the Chi-Chi earthquake at soil sites (EW component):
The observed and simulated accelerograms on the surface (upper left corners), the corresponding vertical
distributions of stresses and strains in the soil layers changing in time during the strong motion (below), spectraof the simulated accelerograms at the bottoms of the soil layers (low-level grey line) and of the observed (dashedline) and simulated (solid line) surface accelerograms (upper right corner), and the profiling data (below): Vs –solid line and smax – dashed line. a – TCU082, TCU054, TCU102, and TCU110, b – TCU116, TCU115,
CHY025, and CHY092, c – CHY101, CHY002, CHY094, and CHY036, d - TCU033, CHY082, CHY027, andCHY032, e – CHY039, TCU017, CHY093, and CHY004, f – CHY033, CHY044, CHY055, and CHY012. Afterthe station codes weighted average distances to the subfaults and distances to the closest subfaults are shown.

396

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 5(Contd.)1.2 bars, *0.4 bars, and 0.2–0.4 bars, respectively. At TCU138 and CHY074, resonating
layers correlate with depths of slightly decreased S-wave velocities (Fig. 4c,f), whereas at
CHY015 the whole upper part of the soil profile down to *60 m (Vs *170–340 m/s),
underlayed by more dense soils (Vs > 600 m/s), is engaged into the resonant excitations

Vol. 165, 2008

Nonlinear Soil Behavior

397

Figure 5(Contd.)(Fig. 4g). According to our estimates, resonant oscillations were induced by the Chi-Chi
earthquake in the upper soil layers down to depths of *40–60 m at stations locatedwithin *30 km from the fault plane.

398

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 5(Contd.)As seen from Figures 4a-g, the strongest manifestations of soil nonlinearity wereobserved at stations TCU065, TCU072, and CHY074. Thus, strong nonlinearity of thesoil behavior during the Chi-Chi earthquake was observed at Sites Class D (TCU065 andVol. 165, 2008

Nonlinear Soil Behavior

399

Figure 5(Contd.)TCU072) and Site Class C (CHY074); TCU065 did show liquefaction, and its possibleclassification as site class F is discussed in LEE et al. (2001a).
Since we have obtained realistic models of soil behavior at the seven stations
possessing most complete information on the profiling data, we continued this study with

400

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,

Figure 5(Contd.)additional twenty-four soil stations, for which the profiling data, P- and S-wave velocities
and SPT-N testing results, were restricted to depths of *30–40 m. It was shown abovethat resonant oscillations were induced by the Chi-Chi earthquake in the upper soil layers

Vol. 165, 2008

Nonlinear Soil Behavior

401(*40–60 m) in near-fault zones. Therefore, knowledge of soil parameters in the upper*30–40 m is not enough for calculations, and profiling data should be completed down
to depths of at least *70–80 m. We reconstructed profiling data at depths *30–80 mbased on the results by LIN et al. (2006), who conducted microtremor measurements atseven stations spaced at a distance of *120 km along the western coastal plain of
Taiwan. Shallow S-wave velocity structures under all the studied sites down to *3 km
were estimated by the inversion of dispersion curves obtained by frequency-wavenumber(f-k) analysis (LIN et al., 2006). The most reliable estimates of S-wave velocities wereobtained for the upper *100 m, and they were verified by P-S logging data in the upper
40–70 m at some nearby TSMIP (Taiwan Strong Motion Instrumentation Program)stations. It was also found that the obtained S-wave velocity models correctly reflect thedepths of major interfaces including the tops of Pliocene and Miocene (LIN et al., 2006).
For the upper *100 m, fairly simple velocity structures and gradual increases of S-wave
velocities with depth were found (LIN et al., 2006). This helped facilitate estimation ofsoil parameters at depths of *30–80 m. The validity of our assumptions is supported by
the obtained (and discussed below) rather good agreement between simulated andobserved accelerograms on the surface.
The results of simulation of acceleration time histories of the Chi-Chi earthquake at
the twenty-four soil stations are presented in Figure 5f-f (only EW components areshown). For each station, the observed and calculated (one of the series of simulated
accelerograms) acceleration records on the surface, their spectra and spectra of ‘‘input’’
motion at the bottoms of the soil layers (averaged over the series of simulated
accelerograms), the profiling data (S-wave velocities Vs and shear stress in failure smax inthe soil layers), as well as vertical distributions of stresses and strains in soil layers,changing in time during strong motion, are shown. Also, weighted average distances toall the subfaults of the fault plane and distances to the closest subfaults (considered as
distances to the fault plane) are indicated.
As seen from the figures, in spite of the shortage of initial information, in some cases,
a fairly good agreement was achieved between observations and simulations: For
example, at stations: TCU082, TCU054, TCU102, TCU110 (Fig. 5a), TCU116,
CHY025, TCU115, CHY092 (Fig. 5b), CHY101, CHY036 (Fig. 5c), CHY027,
CHY032 (Fig. 5d), CHY039 (Fig. 5t), CHY033, CHY044, CHY055, and CHY012(Fig. 5f). For other stations, the agreement between the observations and simulations is
satisfactory.
At stations TCU115, CHY092, CHY002, CHY094, and CHY093 (Fig. 5b, c, e) final
parts of records contain low-frequency intense seismic waves with periods of
3–5 seconds, which cannot be modelled in our 1-D simulations. These oscillations seemto represent effects of subsurface topography; they resemble basin-induced surface
waves, i.e., standing waves due to S-wave reverberation in the basin sediments. Thisconclusion was made for stations CHY026 and CHY104 (Fig. 4d,e). Such waves can also
be distinguished, though not so clearly, at stations CHY082, CHY027, and CHY032(Fig. 5d). All these stations are located in a valley (elevations 4–34 m), close to each

402

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,other, and this area is marked by grey color in Figure 1. At these stations, the soil
composition in the upper *30–40 m is similar: Predominantly silty soils with sands andclays and with average S-wave velocities of *190–260 m/s (Table 3). For these stations,spectra were calculated based on the beginning parts of records, not containing basinwaves.
Long-period seismic waves are also noticed in the final parts of records at stations ofthe ‘‘northern’’ group, TCU033 and TCU017 (Fig. 5d and e). LEE et al. (2001) attributestation TCU017 to a typical case of a tree root problem, where the concrete slab base of
the station is affected by tree roots; it is not clear if basin waves are affecting records atthese stations. An agreement between observations and simulations was obtained for the
beginning parts of records at both of these stations; and spectra were estimated based onthese beginning parts. At the same time, the duration of the simulated strong groundmotion at both stations is shorter than that of the observed ones. This is a problem at allrock and soil stations of the ‘‘northern’’ group. Apparently, this is due to someinconsistency in the applied model of slip distribution over the fault plane.
At stations CHY101, CHY036, and CHY039 (Fig. 5c,e), the recorded accelerograms
possess characteristic shapes, such as, intense beginning parts (first 15–20 sec. of thestrong motion) and final parts of substantially smaller intensity, whereas the intensities of‘‘input’’ motions at the three stations (simulated accelerograms not shown) remained at
almost the same levels during the earthquake. We found that simulation of the recorded
accelerograms at the three stations with one stress-strain relation applied at all timeintervals during the strong motion is not possible; therefore, the beginning and final partsof these records were simulated separately. It was found that for the final parts of records,
an agreement between observations and simulations can only be achieved, if describing
the soil behavior by very slight stress-strain relations, with slopes close to horizontal ones
(Fig. 5c and e). A similar case was the simulation of soil behavior at the Port Island
vertical array site during the Kobe earthquake. At Port Island, liquefaction in the upper
soil layers (*13 m) occurred, and the behavior of these layers was described by similarvery slight stress-strain relations (PAVLENKO and IRIKURA, 2003).
All three stations, CHY101, CHY036, and CHY039, are located close to each other,in a bajada stretched at the southwestern corner of the fault, and their elevations are 75 m,
45 m, and 16 m, respectively (Fig. 1). The upper parts of the soil profiles represent silty
soils with sands and clays with average S-wave velocities of *200–280 m/s (Table 3),
underlayed by gravels and dense clays at depths of *26–34 m. Resonant oscillations
were induced in the upper softer layers during the strong motion, and maximum strainsachieved *0.8% at CHY101 station (with stresses up to *0.4 bars), *0.6% at CHY036
station (with stresses up to *0.4 bars), and *0.3% at CHY039 station (with stresses of*0.1 bars), according to our estimates. Evidently, the obtained results may indicateliquefaction phenomena occurring at stations CHY101, CHY036 and CHY039 during the
Chi-Chi earthquake (Fig. 5c and e), though there is likely another possible explanation.
The highest strains were achieved in the upper *30–40 m at stations TCU065(*0.6%), TCU110 (*0.6–0.8%), and TCU115 (*0.6–0.8%) (Fig. 5a, b), possessing

Vol. 165, 2008

Nonlinear Soil Behavior

403

similar models of soil behavior, such as, expressed resonant oscillations induced bystrong motion: in the upper *9 m and *37 m at TCU065 station, in the upper *30–
40 m at stations TCU110 and TCU115. The soil composition of the upper layers issimilar at these stations: mostly sandy, clayey, and silty soils with average S-wave
velocities of *300 m/s at TCU065 and 190–210 m/s at TCU110 and TCU115 stations(Table 3).
LEE et al. (2001) reported liquefaction occurred at stations TCU065 and TCU110(Figs. 4a and 5a) during the Chi-Chi earthquake. The constructed models of soil behavior
at these stations confirm this conclusion. At TCU115 station (Fig. 5b) strains in thesurface layers (in the upper meters) are not as high as at TCU065 and TCU110 stations,because it is located at a larger distance from the fault plane than TCU065 and TCU110(Table 3). Liquefaction effects are not as clearly expressed at TCU115 as at TCU065 andTCU110, and they involve deeper layers. The constructed models of soil behavior atthese stations indicate that resonant oscillations induced by strong motion in the upper
soft layers are accompanied by large strains and pore pressure build-up, leading to cyclic
mobility and liquefaction of the water-saturated sandy soils (Fig. 5a and b). Similar
phenomena are also observed at TCU102 (strains in the upper layers of *0.2%)(Fig. 5f), CHY025 (strains *0.4%) (Fig. 5b), and at stations TCU072 (strains *0.3–0.4%) and CHY026 (strains *0.3–0.4%), as discussed above. At other stations,maximum strains in the upper soil layers did not exceed *0.2%, according to ourestimates.
Thus, we can conclude that resonant oscillations were induced by the Chi-Chi
earthquake in the soil layers down to *40–60 m in the near-fault zones. As seen from
Figures 4 and 5, resonant phenomena occurred in the upper soft soil layers (with average
S-wave velocities not exceeding *300 m/s), and the resonances were most pronounced
at stations located within *20–22 km from the fault plane, such as, at: TCU065,
TCU110, TCU115, CHY025, CHY092, CHY101, CHY026, CHY036, and CHY039.
Since at these soil sites resonant phenomena occurred in conditions of a rather strong
nonlinearity of soil behavior (as seen from Figures 4 and 5), substantial changes inspectral composition of seismic waves propagating in soil layers appeared: Spectra ofoscillations on the surface took forms close to E(f)*f -n, i.e., resonant amplification wasobserved in the low-frequency range, up to 2–3 Hz (Figs. 4 and 5).
At some of these stations closest to the fault plane (TCU065, TCU110, TCU115,
CHY025, CHY101, and CHY036), our simulations indicate more complicated soil
behavior, i.e., the behavior of different (in their properties) soil layers described bydifferent stress-strain relations. Different stress-strain relations were applied to describe
the behavior of soil layers in cases when the upper sandy-clayey-silty resonating layers
were underlayed by more dense layers (with SPT-N > 50): silty clays (as at TCU065
station), and gravels (as at stations CHY025, CHY101, and CHY036), etc.
On the whole, we found that soil behavior in near-fault zones during the Chi-Chi
earthquake was rather simple, i.e., at the majority of the studied soil stations, the behavior
of all the layers (upper tens of meters) could be satisfactorily described by a single

404

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,stress-strain relation. However, if we possessed more complete information (for example,
vertical array records down to *80–100 m were available), we could reconstruct the soil
behavior in more detail. Apparently, we would obtain changes in stress-strain relations
describing the soil behavior during the strong motion, as was obtained in analysis of theKobe and Tottori earthquakes (PAVLENKO and IRIKURA, 2003, 2004). Evidently, theapproximation of ‘‘one stress-strain relation’’ could be improved upon. Therefore, in thisstudy we estimate some averaged-in-time and smoothed characteristics of soil behavior.
These characteristics are also somewhat approximate, because of the insufficiency of theinitial information: Acceleration time histories at the bottoms of soil layers and the
profiling data down to *70–80 m at some stations are approximated in our study.
Thus, valuable data on in situ soil behavior during the Chi-Chi earthquake is obtained,
and it can be compared with models of soil behavior during the 1995 Kobe (Mw = 6.8)
and 2000 Tottori (Mw = 6.7) (Japan) earthquakes obtained in previous studies (PAVLENKO
and IRIKURA, 2003, 2004, 2005, 2006).
At most of the studied stations in the near-fault zones of the Kobe and Tottori
earthquakes, surface soil layers represented sandy water-saturated soils (the levels of
underground water were less than *10 m), and during strong ground motion, porepressure build-up occurred in the upper soft layers, and the soil behavior was of the‘‘hard-type’’, most expressed in near-fault zones. In the upper soil layers, resonant
oscillations were induced by the strong motion, and amplitudes of oscillations in theupper *10–15 m and on the surface were high. This kind of soil behavior was clearlyobserved at the SGK vertical array (*6 km from the fault plane) during the Kobe
earthquake (PAVLENKO and IRIKURA, 2003) and at the Kik-Net stations TTRH02 and
SMNH01 (*7 km and *8 km from the fault plane, respectively) during the Tottori
earthquake (PAVLENKO and IRIKURA, 2004, 2006).
During the Chi-Chi earthquake, resonant phenomena involved not only the upper*10–15 m, as in the near-fault zones of the Kobe and Tottori earthquakes, but muchthicker layers, up to *40–60 m. Resonance effects of increasing amplitudes ofoscillations on the surface combined with ‘‘hard-type’’ soil behavior, are mostpronounced at stations closest to the fault plane (distances < 10 km from the fault
plane) and are still observed at distances of *10–16 km from the fault. At other (remote)stations, such effects also occur but are less pronounced because they decrease withincreasing the distance from the fault plane. The ‘‘hard’’ character of soil behavior
virtually disappears at distances greater than *12 km from the fault plane, and at
distances of more than *30–40 km, oscillations in the soil layers become very close tolinear ones. A similarity in the behavior of similar soils during these three earthquakes isobserved, indicating a possibility of forecasting soil behavior in future strong
earthquakes.
The constructed models of soil behavior allow estimation of various importantparameters, characterizing soil response during the Chi-Chi earthquake, such as,
amplification of seismic waves in soil layers, average and maximum stresses and strains
induced by the strong motion in soil layers at different depths down to *70–80 m, shear

Vol. 165, 2008

Nonlinear Soil Behavior

405

moduli reduction, nonlinear components of soil response, and others. The distribution ofthese parameters around the fault plane and their dependencies on the distance from thefault can be obtained and analysed. All these problems will be discussed in our next
paper, which is now in preparation.

4. Conclusions
Thus, acceleration time histories at rock and soil stations in near-fault zones of the
Chi-Chi earthquake are simulated, and models of the soil behavior are constructed, i.e.,
vertical distributions of stresses and strains in the soil layers down to *70–140 m are
estimated at thirty one soil sites located within *50 km from the fault.
Because of the essentially stochastic character of the simulated acceleration time
histories, the constructed models of the nonlinear soil behavior can be treated as an
approximation to reality. However, a fairly good agreement was obtained between
spectra of the observed and simulated accelerograms and even between their waveforms,
which testifies to the validity of the obtained representations and allows us to make somedefinitive conclusions.
A decreased attenuation of seismic waves due to geometrical spreading, less than that
described by the relationships 1/r for r < 50 km, 1/50 for 50 B r < 150 km, was found at
stations in the near-fault zones of the Chi-Chi earthquake, which was apparently due tothe effects of the fault zone geometry.
A decreased (compared to other Sites Class B) crustal amplification was establishedat stations: TCU049, TCU053, TCU046, HWA023, HWA024, and CHY110, indicating
probable ‘‘hard rock’’ Sites Class A.
In spite of the large magnitude of the Chi-Chi earthquake and the proximity of thestudied soil stations to the fault plane (and therefore, strong ‘‘input’’ motion), the soil
behavior at the studied sites was found to be relatively simple, i.e., an agreement between
the observed and simulated accelerograms was obtained even in cases of using only one
stress-strain relation for describing the behavior of all the soil layers down to *70–80 m
at all time intervals during strong motion. Obviously, this is due to the homogeneity of
the characteristics of the soil layers in depth.
At all the studied stations, soil behavior during the Chi-Chi earthquake was defined byresonant oscillations of soil layers and nonlinearity of soil response. Resonant oscillations
were induced in the upper meters of soil profiles (down to *40–60 m); strains achieved*0.6% in the vicinity of the fault plane, at TCU065, and 0.1–0.4% at other stations
located within *30 km of the fault. At near-fault stations the relationships between the
stresses and the strains and the soil behavior were substantially nonlinear.
Liquefaction phenomena were identified at TCU065, TCU110, TCU115, CHY101,
CHY036, and CHY039 stations, where estimated strains achieved *0.6–0.8%; at otherstations, maximum strains in the soil layers did not exceed 0.1–0.4%, according to ourestimates.

406

O. V. Pavlenko and K.-L. Wen

Pure appl. geophys.,Thus, the constructed models of soil behavior in the near-fault zones of the Chi-Chi
earthquake add to our representations of soil behavior during strong ground motion, andallow for comparisons with soil behavior in the near-fault zones of the Kobe and Tottori
earthquakes, estimated based on vertical array data, with the obtained data of soil
behavior during the Chi-Chi earthquake. In the near-fault zones of the three earthquakes,‘‘hard-type’’ soil behavior is prevalent along with resonant oscillations in the uppersurface layers; both of which represent unfavourable factors that lead to high acceleration
amplitudes on the surface. Also, a similarity in the behavior of similar soils during theKobe, Tottori, and Chi-Chi earthquakes was found, indicating the possibility offorecasting soil behavior in future earthquakes.

Acknowledgements
This research was supported by the National Science Council (Taiwan) under Grant
No. NSC 95-2811-M-492-001. The profiling data at soil stations provided by NCREE are
greatly appreciated.

REFERENCES
ANDERSON, J. G. and HOUGH, S. E. (1984), A model for the shape of the Fourier amplitude spectrum of
acceleration at high frequencies, Bull. Seismol. Soc. Am. 74, 5, 1969–1993.
ARCHULETA, R. J. (1998), Direct observations of nonlinearity in accelerograms. In The Effects of Surface
Geology on Seismic Motion (eds. Irikura, Kudo, Okada and Sasatani), (Balkema, Rotterdam 1998) pp. 787–
792.
BOORE, D. M. (2003), Simulation of ground motion using the stochastic method, Pure Appl. Geophys. 160, 635–
676.
BOORE, D. M. and JOYNER, W. B. (1997), Site amplifications for generic rock sites, Bull. Seismol. Soc. Am. 87, 2,
327–341.
BAZZURRO, P., SJOBERG, B., and LUCO, N. (2001), Post-elastic response of structures to synthetic ground motions,
http://geohazards.cr.usgs.gov/staffweb/nluco/manuscripts/0409–Bazzurro+Sjoberg+Luco.pdf.
CHEN, K. C., SHIN, T. C., and WANG, J. H. (1989), Estimates of coda Q in Taiwan, Proc. Geol. Soc. China. 32,
339–353.
CHI, W.-C., DREGER, D., and KAVERINA, A. (2001), Finite-source modeling of the 1999 Taiwan (Chi-Chi)
earthquake derived from a dense strong-motion network, Bull. Seismol. Soc. Am. 91, 5, 1144–1157.
CHIU, H. C. (2001), Data Files from the SMART-2 Strong Motion Array for the Chi-Chi Earthquake, Bull.
Seismol. Soc. Am. 91, 5, 1391–1392.
GUSEV, A. A. and PAVLOV, V. M. (2008), Broad band Simulation of Earthquake Ground Motion by a SpectrumMatching, Multiple-Pulse Technique, Earthquake spectra.
HARDIN, B. O. and DRNEVICH, V. P. (1972), Shear modulus and damping in soils: Design equations and curves,
Proc. Am. Soc. Civil Eng., J. Soil Mech. Found. Div. 98, 667–692.
JOYNER, W. B. and CHEN, T. F. (1975), Calculation of nonlinear ground response in earthquakes, Bull. Seismol.
Soc. Am. 65, 1315–1336.
LEE, W. H. K., CHENG, C.-T., LIAO, C.-W., and TSAI, Y.-B. (2001a), Site classification of Taiwan free-fieldstrong-motion stations, Bull. Seismol. Soc. Am. 91, 5, 1283–1297.
LEE, W. H. K., SHIN, T. C., KUO, K. W., CHEN, K. C., and FU, C. F. (2001b), CWB free-field strong-motion datafrom the 21 September Chi-chi, Taiwan, earthquake, Bull. Seismol. Soc. Am. 91, 5, 1370–1376.

Vol. 165, 2008

Nonlinear Soil Behavior

407

LIN C.-M., WEN, K.-L., and CHANG, T.-M. (2006), Estimation of S-wave velocity model in the western coastal
plain of Taiwan, Third Intern. Conf. on Urban Earthquake Engineer, Center for Urban Earthquake
Engineering, Tokyo, 617–622.
MA, K.-F., MORI, J., LEE, S.-J., and YU, S. B. (2001), Spatial and temporal distribution of slip for the 1999 ChiChi, Taiwan, earthquake, Bull. Seismol. Soc. Am. 91, 5, 1069–1087.
OGLESBY, D. and DAY, S. (2001), Fault geometry and the dynamics of the 1999 Chi-Chi (Taiwan) earthquake,
Bull. Seismol. Soc. Am. 91, 5, 1099–1111.
PAVLENKO, O. V. and IRIKURA, K. (2003), Estimation of nonlinear time-dependent soil behavior in strong groundmotion based on vertical array data, Pure Appl. Geophys. 160, 2365–2379.
PAVLENKO, O. V. and IRIKURA, K. (2004), Identification of the nonlinear behavior of soils in near-fault areas
during the 2000 Tottori Earthquake, Proc. of the 11th Intern. Conf. on Soil Dynamics and Earthquake
Engineering and the 3d Intern. Conf. on Earthquake Geotechnical Engineering, 7–9 Jan. 2004, Univ. of
California, Berkeley, USA, pp. 201–207.
PAVLENKO, O. V. and IRIKURA, K. (2005), Identification of the nonlinear behavior of liquefied and non-liquefied
soils during the 1995 Kobe earthquake, Geophys. J. Intern. 160, 539–553.
PAVLENKO, O. V. and IRIKURA, K. (2006), Nonlinear behavior of soils revealed from the records of the 2000
Tottori, Japan, Earthquake at stations of the digital strong-motion network Kik-Net, Bull. Seismol. Soc. Am.
96, 6, 2131–2145.
PAVLENKO, O. V. and LOH, C.-H. (2005), Nonlinear identification of the soil response at Dahan downhole arraysite during the 1999 Chi-Chi earthquake, Soil Dynamics and Earthquake Engin. 25, 3, 241–250.
ROUMELIOTI, Z. and BERESNEV, I. A. (2003), Stochastic finite-fault modeling of ground motions from the 1999
Chi-Chi, Taiwan, Earthquake: Application to rock and soil sites with implications for nonlinear site response,
Bull. Seismol. Soc. Am. 93, 4, 1691–1702.
WANG, W.-H., CHANG, S-H., and CHEN, C.-H. (2001), Fault slip inverted from surface displacements during the1999 Chi-Chi, Taiwan, earthquake, Bull. Seismol. Soc. Am. 91, 5, 1167–1181.
WU, C., TAKEO, M., and IDE, S. (2001), Source process of the Chi-Chi earthquake: A joint inversion of strong
motion data and global positioning system data with a multifault model, Bull. Seismol. Soc. Am. 91, 5, 1128–
1143.
ZENG, Y. and CHEN, C.-H. (2001), Fault rupture of the 20 September 1999 Chi-Chi, Taiwan, earthquake, Bull.
Seismol. Soc. Am. 91, 5, 1088–1098.
(Received April 7, 2007; accepted October 29, 2007)

To access this journal online:www.birkhauser.ch/pageoph