Underground Cavern Drift due to Creep Effect (FLAC3D)
Note
The project file for this example may be viewed/run in FLAC3D.[1] The main data file used is shown at the end of this example.
The implemented enhanced Munson-Dawson model is applied to a typical multi-zone boundary-value example, specifically the drift simulation of an underground cavern excavation in rock salt. This model is a 2.5D representation with dimensions of 32 m × 80 m in the vertical plane. Due to symmetry, it effectively represents a calculation domain of 64 m × 80 m. The model’s top is situated underground at a depth of 200 meters, and the salt rock is characterized by the material properties of Set I listed in Table 1, serving as a benchmark scenario. The temperature is maintained constant at 313.15 K (40 °C).
Initial in-situ stress is assumed that the model top is under 200 meters from the ground surface and the vertical stress is determined through gravity calculations based on a uniform density of 2300 \(kg/m^3\). The horizontal stresses are assumed based on \(K_0\) = 1.0.
The material properties are summarized in the table below.
shear
\(G\)
12.4 (GPa)
bulk
\(K\)
20.7 (GPa)
hosford
\(a\)
4
f-exponent
\(\chi\)
2
A-0
\(A_0\)
5.617e1 (\(s^{-1}\))
activation-ratio-0
\(Q_0/R\)
5123 (K)
n-0
\(n_0\)
1.595
A-1
\(A_1\)
8.386e22 (\(s^{-1}\))
activation-ratio-1
\(Q_1/R\)
12580.5 (K)
n-1
\(n_1\)
5.5
A-2
\(A_2\)
4.415e16 (\(s^{-1}\))
activation-ratio-2
\(Q_2/R\)
5123 (K)
n-2
\(n_2\)
6.279
B-1
\(B_1\)
6.086e6 (\(s^{-1}\))
B-2
\(B_2\)
3.034e-2 (\(s^{-1}\))
stress-limit
\(\sigma_0\)
20.57e6 (MPa)
stress-constant
\(q\)
5335
k-0
\(k_0\)
5.277e-2
c-0
\(c_0\)
8.882e-3 (\(K^{-1}\))
m-0
\(m_0\)
0.9201
k-1
\(k_1\)
3.052e12
c-1
\(c_1\)
8.882e-3 (\(K^{-1}\))
m-1
\(m_1\)
5.282
alpha-hardening
\(\alpha_h\)
3.367
beta-hardening
\(\beta_h\)
-0.6838
alpha-recovery
\(\alpha_h\)
0.58
beta-recovery
\(\beta_r\)
0.0
A horizontal circular cavern with a radius of 2 meters, located at the center right of the model domain, is excavated at t = 0 days, and the model is cycled to static equilibrium with the creep option disabled. Subsequently, the creep option is activated, and the model is run for 10 years. An automatic timestep scheme is employed, starting with a timestep of 1e−20 s and a maximum timestep of 1000 s. To ensure a gradual adjustment of the timestep, a lower multiplier of 1.002 and an upper multiplier of 0.998 are utilized. The horizontal drift at the utmost side point of the cavern wall is monitored throughout the simulation. Figure 1 plots the horizontal drift at the utmost side point of the cavern wall throughout 10 years. Figure 2 presents the horizontal displacement contour at t = 5 years.
Figure 1: Horizontal drift horizontal drift at the utmost side point of the cavern wall.
Figure 2: Horizontal displacement contour at t = 5 years.
References
Cheng, Z. Implementation and Verification of Enhanced Munson-Dawson Creep Model for Rock Salt. In US Rock Mechanics/Geomechanics Symposium (D041S047R003), ARMA (2024).
Reedlunn, B. Enhancements to the Munson-Dawson model for rock salt (No. SAND-2018-12601). Sandia National Laboratories, Albuquerque, NM, USA (2018).
Data File
ini.dat
model new
model large-strain off
model configure creep
model gravity 9.81
;
program call "grid.dat"
zone face skin
;
zone face apply velocity-x 0.0 range group "west" or "east"
zone face apply velocity-y 0.0 range group "north" or "south"
zone face apply velocity-z 0.0 range group "bottom"
;
model creep active off
;
zone cmodel assign elastic
zone property density 2300 shear 1.24e10 bulk 2.07e10
;
zone face apply stress-zz [-2300*9.81*200] range group "top"
zone initialize-stresses ratio 1.0 overburden [-2300*9.81*200]
model solve convergence 1
;
zone null range group "excavation"
model solve convergence 1
;
model save "exc"
undergroundDrift.dat
model restore "exc"
;
zone gridpoint initialize displacement (0,0,0)
zone gridpoint initialize velocity (0,0,0)
zone initialize state 0
;
model creep active on
zone cmodel assign munson-dawson range group "rock"
zone property density 2300 temperature 313.15
;; Set I
zone property bulk 20.7e9 shear 12.4e9 f-exponent 2 hosford 4
zone property A-0 5.617e1 activation-ratio-0 5123 n-0 1.595
zone property A-1 8.386e22 activation-ratio-1 12580.5 n-1 5.5
zone property A-2 4.415e16 activation-ratio-2 5123 n-2 6.279
zone property B-1 6.086e6 B-2 3.034e-2
zone property stress-limit 20.57e6 stress-constant 5335
zone property k-0 5.277e-2 c-0 8.882e-3 m-0 0.9201
zone property k-1 3.052e12 c-1 8.882e-3 m-1 5.282
zone property alpha-hardening 3.367 beta-hardening -0.6838
zone property alpha-recovery 0.58 beta-recovery 0.0
;
model creep timestep starting 1e-20
model creep timestep automatic
model creep timestep minimum 1e-20
model creep timestep maximum 5e3
model creep timestep upper-multiplier 0.998
model creep timestep lower-multiplier 1.002
;
history interval 1000
model history name 'time' creep time-total
zone history name 'disp-x' displacement-x position 30 0 0
;
model solve creep time-total [3600*24*365*10]
history export 'time' 'disp-x' file 'drift-md.txt' truncate
model save "drift-md"
⇐ Munson-Dawson Model: Biaxial Test (FLAC3D) | WIPP-Drucker Model: Comparison with Drucker-Prager Model (FLAC2D) ⇒
| Was this helpful? ... | Itasca Software © 2024, Itasca | Updated: Dec 09, 2025 |