# Comparison of FLAC3D to *SHAKE* for a Layered, Nonlinear-Elastic Soil Deposit

Note

To view this project in FLAC3D, use the menu command . Choose “Dynamic/ F3ShakeNonLinearElasticCase” and select “F3ShakeNonLinearElasticCase.prj” to load. The main data files used are shown at the end of this example. The remaining data files can be found in the project.

In this section, we compare a simulation of wave propagation through a geological profile
using FLAC3D with hysteretic damping to *SHAKE-91*.
The FLAC3D and *SHAKE* models simulate the following problem: a horizontally layered soil
deposit overlying a rigid bedrock is subjected to a horizontal acceleration base motion.
The soil deposit is 150 feet (45.72 m) deep and contains 10 different soil types. The soils
are treated as nonlinear elastic materials by assuming that shear moduli and damping are
strain-dependent. Table 1 summarizes the properties and
locations of the soil layers in the deposit. The dynamic characteristics of these soils are
governed by two sets of modulus reduction factor (\(G/G_{\rm max}\)) and damping ratio
(\(\lambda\)) curves: the first set (set 1) for clay, and the second set (set 2) for
sand. [1] These curves are shown in Figure 1 through
4, and are denoted by “SHAKE91” in the legend.

Soil | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 |

Shear Modulus (MPa) | 186 | 150 | 168 | 186 | 225 | 327 | 379 | 435 | 495 | 627 |

Density (Kg/m^{3}) |
2000 | 2000 | 2000 | 2000 | 2000 | 2082 | 2082 | 2082 | 2082 | 2082 |

Dynamic Property (set) | 2 | 2 | 2 | 1 | 1 | 2 | 2 | 2 | 2 | 2 |

Location (feet) | 1-5 | 5-20 | 20-30 | 30-50 | 50-70 | 70-90 | 90-110 | 110-130 | 130-140 | 140-150 |

This example is derived from the data file provided with the downloaded *SHAKE-91* code (see https://nisee.berkeley.edu/elibrary/software.html). The *SHAKE-91* data file is modified such that the deepest layer is given a large wave speed to correspond to rigid bedrock. Note that to simulate a rigid base, we increase the stiffness (shear modulus) of the bedrock from 3.33 GPa to 2080 GPa in the *SHAKE-91* data file for the purpose of comparison.

The base acceleration input is a set of seismic data recorded in the Loma Prieta Earthquake, which is also downloadable (see https://nisee.berkeley.edu/elibrary/getpkg?id=SHAKE91) with the *SHAKE-91* code (i.e., “DIAM.ACC”). The input accelerogram is shown in Figure 5. The record has a peak acceleration of approximately 0.11 g, duration of 40 seconds, and dominant frequency of 2.47 Hz.

*SHAKE model* — The free-field response of the layered soil deposit to the base-motion
acceleration is calculated using *SHAKE-91*. An iterative solution is performed until the secant shear moduli and the damping ratios throughout the column converge to within a given tolerance of the values corresponding to 50% of the maximum shear strain. Calculations are made for different input accelerations, with the record scaled to a maximum base acceleration. [2]

FLAC3D *model with hysteretic damping* — A 1D column FLAC3D model composed of 30 cubic
zones, each with a depth of 5 feet (1.524 m), is assigned the soil properties listed in
Table 1. Model conditions are prescribed to simulate the
free-field motion of a layered soil deposit with a rigid base. The horizontal acceleration
record is applied as an input boundary condition at the base of the model. FLAC3D takes the
scaled base motion from *SHAKE-91*. For example, note that the input table “Diam-flac-0001.acc” is scaled with a maximum value of 0.0001 g. Vertical movement is
prevented in the model.

Hysteretic damping is added to this model by fitting a secant-modulus function curve to the (\(G/G_{\rm max}\)) versus strain data, as shown in Figures 1 and 3. Different fitting functions are available in FLAC3D. In this example, the default (two-parameter) function is used. The parameters of the secant-modulus functions are determined by best-fitting (e.g., least squares regression analysis) the laboratory modulus-reduction curves, assuming that the modulus reduction factor versus shear strain follows the default relation. The best-fit values for dynamic property set 1 — sand are listed here, and for dynamic property set 2 — clay here CS: this is not cool, but understandable. Preceding x-refs want better labeling. A comparison between the secant-modulus function curves and the laboratory curves is shown in Figure 1 and 3.

Damping ratio curves that correspond to the secant-modulus function curves are then calculated. [3] Both shear modulus reduction factor (\(G/G_{\rm max}\)) versus strain and damping ratio (\(\lambda\)) versus strain should be monitored, and parameters adjusted, if necessary, to provide a reasonable fit to both curves over a specified strain range. As shown in Figures 1 through 4, fairly reasonable fits for both modulus-reduction factor and damping-ratio curves are obtained with the default damping function and the parameters listed in Numerical Fits to Seed and Idris over a strain range of 0.0001% to 0.1% for both the sand and the clay.

The fitting curves of the FLAC3D default model
(Figures 1 through
4) show that, for small strain (corresponding to small
input acceleration), the behavior is approximately linear (i.e., both shear modulus and
damping ratio are constants for both FLAC3D and *SHAKE-91*). Both codes are thus expected
to give similar results in this circumstance. Here we compare accelerograms and response
spectra at the top of the model for very low input acceleration.
Figure 6 shows the horizontal acceleration at the top of the
model (gridpoint 121 in FLAC3D and sub-layer 1 in *SHAKE-91*) as a function of time with
maximum input acceleration amplitude of 0.0001 g. Both records are very similar; maximum
accelerations calculated by FLAC3D and *SHAKE-91* have less than 0.1% difference.

Figures 7 and 8 provide
pseudo-acceleration and pseudo-velocity spectra, calibrated in *SHAKE-91* and FLAC3D when
the maximum input acceleration amplitude is small (0.0001 g). In *SHAKE-91*, the response
spectra are calculated using option 9. In FLAC3D, response spectra are computed using a
FISH function, spectra, implemented in file “spectra.dat”. Here the damping ratio,
minimum period and maximum period of interest are 5%, 0.01 and 10, respectively. From these
plots it can clearly be seen that the FLAC3D and *SHAKE-91* results correspond quite
closely for small-amplitude input motion.

For validation purposes, acceleration amplifications at the top of the model, defined as the
ratio of the maximum acceleration value at the top to the maximum base input acceleration
value, are compared between FLAC3D and *SHAKE-91*. We measure the ratios when applying
scaled maximum input acceleration value in “DIAM.ACC” of 0.0001 g, 0.0005 g, 0.001 g, 0.005
g, 0.01 g, 0.05 g, 0.1 g, 0.5 g, and 1 g, where g is the acceleration of gravity (9.81 m/s^{2}). FLAC3D simply takes the scaled base motion from *SHAKE-91*. According to the
calculation results shown in Figure 9, it is confirmed that
FLAC3D and *SHAKE-91* predict similar acceleration amplification.

References

Hardin, B. O., and V. P. Drnevich. “Shear Modulus and Damping in Soils: I. Measurement and Parameter Effects, II. Design Equations and Curves,” Technical Reports UKY 27-70-CE 2 and 3, College of Engineering, University of Kentucky, Lexington, Kentucky. [These reports were later published in the *Journal of Soil Mechanics and Foundation Division, ASCE,* Vol. 98, No. 6, pp. 603-624 and No. 7, pp. 667-691, in June and July 1972.]

Idriss, I. M., and Joseph I. Sun. User’s Manual for SHAKE-91. University of California, Davis, Center for Geotechnical Modeling, Department of Civil & Environmental Engineering (November 1992).

Seed, H. Bolton, and I. M. Idriss. “Soil Moduli and Damping Factors for Dynamic Response Analysis,” Earthquake Engineering Research Center, University of California, Berkeley, Report No. UCB/EERC-70/10, p. 48 (December 1970).

Seed, H. Bolton, et al. “Moduli and Damping Factors for Dynamic Analyses of Cohesionless Soils,” *Journal of the Geotechnical Engineering Division, ASCE*, Vol. 112, No. GT11, pp. 1016-1032 (November 1986).

Sun, J. I., R. Golesorkhi and H. Bolton Seed. “Dynamic Moduli and Damping Ratios for Cohesive Soils,” Earthquake Engineering Research Center, University of California, Berkeley, Report No. UCB/EERC-88/15, p. 42 (1988).

Vucetic, M., and R. Dobry. “Effect of Soil Plasticity on Cyclic Response,” in *J. Geotech. Eng. Div., ASCE*, Vol. 111, No. 1, pp. 89-107 (January 1991).

Data Files

**F3ShakeNonLinearElasticCase.dat**

```
;---------------------------------------------------------------------
; Comparison of Hysteretic Damping with Shake91
; This file uses SHAKE91-comparison-flac3d.f3dat
; to create the acceleration
; amplification ratio at the top of layer 1. It changes
; the acceleration amplifications for different runs.
; It also creates the graphs and checks for errors.
;---------------------------------------------------------------------
model new
fish automatic-create off
model title 'Comparison of Hysteretic Damping with Shake91'
model configure dynamic
;----------------------------
; Load applied acceleration
;----------------------------
table 'iacc' import 'Diam-flac-0001.acc'
; Run 9 different cases of acceleration multiplier,
; storing the acceleration amplification result for each
[global all_mult = map(1,1,2,5,3,10,4,50,5,100,6,500,7,1000,8,5000,9,10000)]
; Map of acceleration multiplier cases, 9 cases from 1 to 10000
fish define run_all_cases
global acc_case
loop foreach acc_case map.keys(all_mult) ; All case numbers
global acc_mult = all_mult(acc_case) ; Acceleration multiplier
; for that case
command ; program call the data file that runs it
program call 'comparison'
end_command
end_loop
end
[run_all_cases]
model save 'final'
```

**comparison.dat**

```
model large-strain off
;---------------------------------------------------------------------
; Comparison of Hysteretic Damping with Shake91
;---------------------------------------------------------------------
; First, reset model without clearing everything
zone delete ; Remove all zones
history delete ; Remove all histories
mode dynamic time-total 0 ; Reset dynamic time to 0
;-------------------------------------
;Grid generation and model properties
;-------------------------------------
zone create brick size 1 1 30
zone cmodel assign elastic
zone property bulk 300e6 shear 186e6 density 2000 ...
range position-z 29 30
zone property bulk 200e6 shear 150e6 density 2000 ...
range position-z 26 39 ; 5- 20 ft
zone property bulk 200e6 shear 168e6 density 2000 ...
range position-z 24 26 ; 20- 30 ft
zone property bulk 270e6 shear 186e6 density 2000 ...
range position-z 20 24 ; 30- 50 ft
zone property bulk 350e6 shear 225e6 density 2000 ...
range position-z 16 20 ; 50- 70 ft
zone property bulk 480e6 shear 327e6 density 2082 ...
range position-z 12 16 ; 70- 90 ft
zone property bulk 550e6 shear 379e6 density 2082 ...
range position-z 8 12 ; 90-110 ft
zone property bulk 600e6 shear 435e6 density 2082 ...
range position-z 4 8 ;110-130 ft
zone property bulk 750e6 shear 495e6 density 2082 ...
range position-z 2 4 ;130-140 ft
zone property bulk 900e6 shear 627e6 density 2082 ...
range position-z 0 2 ;140-150 ft
;------------
; Histories
;------------
model history name 'time' dynamic time-total ; 2
zone history name 'atop' acceleration-x position (0,0,30) ;top accel 0'
;--------------------
;Boundary Conditions
;--------------------
zone gridpoint fix velocity-y
zone gridpoint fix velocity-z
zone face apply acceleration-x [acc_mult*9.81] table 'iacc' ...
time dynamic range position-z = 0
;--------------------
; Hysteretic damping
;--------------------
zone dynamic damping hysteretic default -3.325 0.823 range position-z 0 16
zone dynamic damping hysteretic default -3.156 1.904 range position-z 16 24
zone dynamic damping hysteretic default -3.325 0.823 range position-z 24 30
;----------------------------------------------
; Increase size by 1.524 to match shake
;----------------------------------------------
zone gridpoint initialize position (1.524,1.524,1.524) multiply
;----------------------------------------------
; Solve to 40.48s, taking histories every 0.02s.
;----------------------------------------------
model dynamic timestep fix 0.0002
history interval 100
model solve time-total 40.48
;top accel hist to table in order to generate response spectra
history export 'atop' vs 'time' table ['tacc'+string(acc_case)]
```

**NonlinearElasticCase-SHAKE.dat**

```
option 1 -- dynamic soil properties -- (max is thirteen):
1
3
11 #1 modulus for clay (Seed & Sun 1989) upper range
0.0001 0.0003 0.001 0.003 0.01 0.03 0.1 0.3
1. 3. 10.
1.000 1.000 1.000 0.981 0.941 0.847 0.656 0.438
0.238 0.144 0.110
11 damping for clay (Idriss 1990) --
0.0001 0.0003 0.001 0.003 0.01 0.03 0.1 0.3
1. 3.16 10.
0.24 0.42 0.8 1.4 2.8 5.1 9.8 15.5
21. 25. 28.
11 #2 modulus for sand (Seed & Idriss 1970) -- upper range
0.0001 0.0003 0.001 0.003 0.01 0.03 0.1 0.3
1. 3. 10.
1.000 1.000 0.990 0.960 0.850 0.640 0.370 0.180
0.080 0.050 0.035
11 damping for sand (Idriss 1990) -- (about LRng from SI 1970)
0.0001 0.0003 0.001 0.003 0.01 0.03 0.1 0.3
1. 3. 10.
0.24 0.42 0.8 1.4 2.8 5.1 9.8 15.5
21. 25. 28.
8 #3 ATTENUATION OF ROCK AVERAGE
.0001 0.0003 0.001 0.003 0.01 0.03 0.1 1.0
1.000 1.000 0.9875 0.9525 0.900 0.810 0.725 0.550
5 DAMPING IN ROCK
.0001 0.001 0.01 0.1 1.
0.4 0.8 1.5 3.0 4.6
3 1 2 3
Option 2 -- Soil Profile
2
1 17 Example -- 150-ft layer; input:Diam [] .1g
1 2 5.00 .050 .125 1000.
2 2 5.00 .050 .125 900.
3 2 10.00 .050 .125 900.
4 2 10.00 .050 .125 950.
5 1 10.00 .050 .125 1000.
6 1 10.00 .050 .125 1000.
7 1 10.00 .050 .125 1100.
8 1 10.00 .050 .125 1100.
9 2 10.00 .050 .130 1300.
10 2 10.00 .050 .130 1300.
11 2 10.00 .050 .130 1400.
12 2 10.00 .050 .130 1400.
13 2 10.00 .050 .130 1500.
14 2 10.00 .050 .130 1500.
15 2 10.00 .050 .130 1600.
16 2 10.00 .050 .130 1800.
17 3 .010 .140 100000.
Option 3 -- input motion:
3
1900 4096 .02 diam.acc (8f10.6)
.0001 25. 3 8
Option 4 -- sublayer for input motion within (1) or outcropping (0):
4
17 0
Option 5 -- number of iterations & ratio of avg strain to max strain
5
0 8 0.50
Option 6 -- sublayers for which accn time histories are computed & saved:
6
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15
0 1 1 1 1 1 1 1 1 1 1 1 1 1 1
1 0 0 0 0 1 0 0 0 1 0 0 0 1 0
Option 6 -- sublayers for which accn time histories are computed & saved:
6
16 17 17
1 1 0
0 1 0
option 7 -- sublayer for which shear stress or strain are computed & saved:
7
4 1 1 0 1800 -- stress in level 4
4 0 1 0 1800 -- strain in level 4
option 7 -- sublayer for which shear stress or strain are computed & saved:
7
8 1 1 0 1800 -- stress in level 8
8 0 1 0 1800 -- strain in level 8
option 9 -- compute & save response spectrum:
9
1 0
1 0 981.0
0.05
option 10 -- compute & save amplification spectrum:
10
17 0 1 0 0.125 -- surface/rock outcrop
execution will stop when program encounters 0
0
```

Endnotes

[1] | More variations of \(G/G_{\rm max}\) and \(\lambda\) for soils are available in the literature mentioned in the SHAKE-91 manual (Idriss and Sun 1992) and duplicated here for easy reference (e.g., Hardin and Drnevich 1972, Seed and Idriss 1970, Seed et al. 1986, Sun et al. 1988, and Vucetic and Dobry 1991). |

[2] | The scaling is accomplished through the built-in functionality in SHAKE-91, using either of two scaling parameters in OPTION 3. |

[3] | The data file in the project Several Cyclic Strain Levels can be used to perform a series of simple shear tests on a single-zone FLAC3D model in order to develop damping ratio curves that correspond to the secant-modulus functions. |

⇐ Comparison of FLAC3D to

*SHAKE*for a Layered, Linear-Elastic Soil Deposit | Slip Induced by Harmonic Shear Wave ⇒

Was this helpful? ... | FLAC3D © 2019, Itasca | Updated: Mar 10, 2023 |