Dynamic response analysis of concrete lining structure in high pressure diversion tunnel under seismic load

Deng, Jian; Xiao, Ming

doi:10.21595/jve.2015.16289

Journal of Vibroengineering

Browse Journal

Submit article

Published: 31 March 2016

Check for updates

Dynamic response analysis of concrete lining structure in high pressure diversion tunnel under seismic load

Jian Deng¹

Ming Xiao²

^{1, 2}State Key Laboratory of Water Resources and Hydropower Engineering Science, Wuhan University, Wuhan 430072, China

^{1, 2}Key Laboratory of Rock Mechanics in Hydraulic Structural Engineering of Ministry of Education, Wuhan University, Wuhan 430072, China

Corresponding Author:

Ming Xiao

Cite the article Download PDF

Downloads 1354

CrossRef Citations 3

Abstract

Dynamic response of concrete lining structure in high pressure diversion tunnel under seismic load is of great importance for engineering design and safety assessment. In light of the basic idea of static damage constitutive in strain space, a dynamic damage and cracking model of concrete suitable for programming is established. Based on the explicit finite element method for analyzing the dynamic response of single-phase solid and fluid medium, combining with the boundary conditions at coupling interface, a dynamic explicit finite element solving format considering the coupling interaction of lining and inner water is presented. This method can directly integrate, solve without simultaneous equations, and greatly simplify the calculation process. Moreover, it can be easily used to analyze wave propagation problems with a variety of mediums. According to the damage evolution properties of concrete lining under seismic load, a calculation method of permeability coefficient considering the effect of damage and cracking is put forward, and the inner water exosmosis process with lining cracks is analyzed. The calculation results of engineering example can reasonably reflect the seismic response characteristics of concrete lining structure in high pressure diversion tunnel, and an effective analysis method is provided for the aseismic design of hydraulic tunnel.

1. Introduction

The spatial distribution of water resources is extremely uneven worldwide, which highlights the contradiction between water supply and demand. The development and utilization of water resources is gradually turning to the upstream of each watershed. So the long distance basin water diversion and power generation become the important developing directions of water conservancy and hydropower projects in the long run. The long distance diversion tunnels usually cross high earthquake intensity region, facing serious seismic stability problem [1-4]. Moreover, with the vigorous development of economic construction, more and more diversion tunnels with high head and large diameter appear, and the liquid sloshing with high pressure and great mass brings huge challenge to the seismic stability. Therefore, researching the anti-seismic stability of long distance water diversion project is of great significance for guaranteeing the safe and steady operation of the project.

The key to seismic stability analysis of high pressure diversion tunnel lies in the seismic response of concrete lining structure, so it is particularly important to establish a scientific and practical dynamic constitutive model for concrete. The current dynamic constitutive models are primarily nonlinear elasticity [5], plasticity [6], viscoplasticity [7, 8] and damage model [9, 10], which are usually established on the basis of experimental studies, and cannot be widely applied on account of many uncertain parameters, ignoring the fatigue damage characteristics under dynamic cyclic loading, iterative calculation which adversely affects the speed and precision. Furthermore, it lacks an effective simulation for dynamic damage and cracking process of lining structure. When analyzing the dynamic response of concrete lining structure under earthquake excitation, the interaction of inner water and lining should be taken into consideration, which is also proposed in aseismic design codes at home and abroad. However, current researches aimed at the dynamic interaction between lining structure and water in high pressure diversion tunnel under seismic load are very few, and the calculation and evaluation of hydrodynamic pressure are not clear. The interaction between inner water and lining structure belongs to the typical fluid-structure interface coupling problem [11, 12]. Studies on fluid-structure coupling problem have quickly developed in hydraulic, marine and other fields since 1980s, such as the liquid sloshing in the liquid reservoir vessel [13, 14], dam reservoir water fluctuation [15, 16] and the structure vibration of deep seabed tunnels[17, 18]. However, fluid-structure coupling analysis method is rarely applied to the high head and large diameter water conveyance tunnel. The effect of inner water is simply considered with the added mass method. The interaction law of lining structure and fluid cannot be reflected, and the seismic response characteristics of lining structure are not clear as well.

On the basis of former researches, a new dynamic explicit finite element analysis method considering the coupling interaction of fluid and lining is presented. Taking the dynamic damage and cracking process of concrete lining structure under seismic load into account, the inner water exosmosis process with lining cracks is analyzed. The research results are of important reference value in solving aseismic safety problem and enhancing the anti-earthquake designing level of hydraulic tunnel.

2. Dynamic damage and cracking model for concrete lining

Since concrete is a heterogeneous material with micro cracks, the macroscopic defects gradually produce and accumulate under seismic cyclic loading, resulting in the deterioration of material nature. Based on the theoretical research, combining with abundant indoor and outdoor tests, a dynamic damage constitutive model for concrete is established according to the basic idea of static damage constitutive [19, 20]. In light of damage states of elements, the dynamic cracking process is analyzed.

2.1. Dynamic damage evolution law

The static damage constitutive relation based on strain space can be expressed as follows:

1

σ_{s} = E_{s} ε (I - D_{s}),

where $E_{s}$ is the static initial elastic modulus for concrete; $D_{s}$ is the second order static damage tensor; $I$ is the unit matrix.

According to abundant indoor and outdoor tests [21, 22], the static and dynamic constitutive curves of concrete are shown in Fig. 1, in which curve 1 is static condition and curve 2 refers to dynamic conditon. The two curves are proved to be highly similar on the basis of some experimental results, and the geometric relationship is satisfied:

2

\frac{A B}{C D} = \frac{O A}{O C} .

With reference to the static damage constitutive relation, the constitutive equation for concrete under seismic load can be assumed as follows:

3

σ_{d} = E_{d} ε (I - D_{d}),

where $E_{d}$ is the dynamic initial elastic modulus for concrete, which could be approximatively equal to $E_{s}$ with no special experiments; $D_{d}$ is the second order dynamic damage tensor.

Fig. 1Static and dynamic constitutive curves of concrete

According to the similarity relation, we have:

4

[I - D_{d} (ε_{d}^{u})] = \frac{f_{σ} (\dot{ε})}{f_{ε} (\dot{ε})} [I - D_{s} (ε_{s}^{u})],

where $f_{σ} (\dot{ε}) = σ_{d}^{u} / σ_{s}^{u}$ , $f_{ε} (\dot{ε}) = ε_{d}^{u} / ε_{s}^{u}$ ; $σ_{d}^{u}$ and $ε_{d}^{u}$ are the ultimate stress and strain of dynamic damage for concrete; $σ_{s}^{u}$ and $ε_{s}^{u}$ are the ultimate stress and strain of static damage for concrete.

In fact, for each pair of corresponding points on the static and dynamic curves, the above formula is tenable. But $f_{σ} (\dot{ε})$ and $f_{ε} (\dot{ε})$ are varying along with the change of strain and strain rate, which can be written in general form:

5

[I - D_{d} (ε^{d})] = \frac{f_{σ} (\dot{ε}, ε)}{f_{ε} (\dot{ε}, ε)} [I - D_{s} (\frac{ε^{d}}{f_{ε} (\dot{ε}, ε)})] .

2.2. Dynamic damage constitutive model

In order to simply and approximately describe the change law of $f_{σ} (\dot{ε}, ε)$ and $f_{ε} (\dot{ε}, ε)$ , two basic assumptions are introduced. For any corresponding point on the dynamic curve: (1) Stress amplification coefficient is the same; (2) Stain amplification coefficient is the same. We have:

6

f_{σ} (\dot{ε}, ε) = f_{σ} (\dot{ε}), f_{ε} (\dot{ε}, ε) = f_{ε} (\dot{ε}) .

Based on many statistical results of dynamic loading tests for concrete, the calculation formulas for dynamic amplification coefficients of tension, compression stress and strain are provided by Euro-International Committee for Concrete (CEB) [23]:

7

\{\begin{array}{l} f_{σ}^{l} (\dot{ε}) = {(\frac{{\dot{ε}}_{d}}{{\dot{ε}}_{s}})}^{1.016 α}, \\ f_{σ}^{y} (\dot{ε}) = {(\frac{{\dot{ε}}_{d}}{{\dot{ε}}_{s}})}^{1.026 β}, \\ f_{ε} (\dot{ε}) = {(\frac{{\dot{ε}}_{d}}{{\dot{ε}}_{s}})}^{0.020}, \end{array}

where ${\dot{ε}}_{s} =$ 3×10^-6 s^-1 is static strain rate; ${\dot{ε}}_{d}$ is dynamic strain rate, which is between 3×10^-6 s^-1and 30 s^-1; $α = {(10 + 3 f_{c}^{'} / 5)}^{- 1}$ , $β = {(5 + 3 f_{c u} / 4)}^{- 1}$ ; $f_{c}^{'}$ and $f_{c u}$ are static cylinder and cube compressive strength of concrete, respectively.

According to the basic assumption Eq. (2), the dynamic damage threshold strain for concrete can be obtained:

8

ε_{d}^{0} = f_{ε} (\dot{ε}) ε_{s}^{0} .

It indicates the dynamic damage threshold strain for concrete is proportional to the static condition, of which the static damage threshold strain $ε_{s}^{0}$ refers to ultimate tensile strain $ε_{t}$ and compression strain $ε_{c}$ according to the stress states.

The dynamic damage constitutive equation for concrete can be obtained:

9

\{\begin{array}{l} σ_{d} = E_{d} ε, & ε \leq ε_{d}^{0}, \\ σ_{d} = E_{d} ε \frac{f_{σ} (\dot{ε})}{f_{ε} (\dot{ε})} [I - D_{s} (\frac{ε}{f_{ε} (\dot{ε})})], & ε > ε_{d}^{0} . \end{array}

In order to display the damage state of lining structure in the post-processing, the second-order damage tensor is converted to a first-order damage scalar when outputting the damage results, which is applied to represent the whole damage state of element, namely $D = \sqrt{D_{1}^{2} + D_{2}^{2} + D_{3}^{2}}$ .

The tensile and compressive damage evolution equations for concrete under static condition can be expressed as [24]:

10

\{\begin{matrix} D_{s}^{l} = 1 - e x p (- R \sqrt{(ε - ε_{t}) (ε - ε_{t})}), & ε > ε_{t}, \\ D_{s}^{y} = 1 - e x p (- R \sqrt{(ε - ε_{c}) (ε - ε_{c})}), & ε > ε_{c}, \end{matrix}

where $R$ is the damage constant of concrete material.

2.3. Dynamic cracking process simulation

Without cracks, the concrete could be regarded as a linear-elastic isotropic material. For simplification, this paper assumes that the trend of crack growth is perpendicular to the direction of maximum principal strain. It should be noticed that since the compressive strength of concrete is far more than tensile strength, the damage and cracking mainly refer to tensile crack. In the calculation process, if the actual strain of concrete element exceeds the ultimate tensile strain, crack will appear in the lining. The schematic of lining cracking is shown in Fig. 2, in which $X'$ refers to the direction of maximum principal strain, and it is perpendicular to the crack interface. The maximum principal stress of the cracked element will be released and stress will be redistributed. Consequently, the stiffness coefficient along the direction of maximum principal strain will be equal to zero, and the constitutive relation of the cracked element will become anisotropic, which should be corrected by the damage and cracking state of element. The stress matrix of cracked element under local coordinate system can be obtained [25]:

11

[H_{e}^{'}] = [\begin{array}{l} 0 & 0 & 0 \\ d_{22} - \frac{d_{12}^{2}}{d_{11}} & d_{23} - \frac{d_{23}^{2}}{d_{11}} & 0 \\ d_{33} - \frac{d_{22}^{2}}{d_{11}} \\ η d_{44} & 0 & 0 \\ s y m m e t r y & η d_{55} & 0 \\ η d_{66} \end{array}],

where:

d_{11} = d_{22} = d_{33} = \frac{E_{d} (1 - ν)}{(1 + ν) (1 - 2 ν)}, d_{12} = d_{13} = d_{23} = \frac{ν E_{d}}{(1 + ν) (1 - 2 ν)},

d_{44} = d_{55} = d_{66} = \frac{E_{d}}{2 (1 + ν)}, η = 1 - D,

and $ν$ , $E_{d}$ are Poisson’s ratio and elastic modulus of concrete, respectively; $η$ is the residual shear coefficient, which is used to describe the residual shear capacity after concrete cracking; $D$ is the element damage variable. Through the coordinate transformation, the stress matrix of damage and cracked elements in the global coordinate system can be obtained.

Fig. 2Schematic of lining cracking

Fig. 3Combined bearing model of lining and steel

The combined bearing model of lining and steel in high pressure diversion tunnel is shown in Fig. 3. When the load is small, lining and steel work together with continuous deformation, otherwise, the bond-slip phenomenon occurs, resulting in cracks along the radial direction of tunnel. It is well known that the crack width has a close relationship with the thickness of protection layer, the spacing and strain of steel. The following formula is often used to estimate the crack width [26, 27]:

12

W_{m a x} = 4 ε_{s} t_{e},

where $W_{m a x}$ is the maximal probability crack width; $ε_{s}$ is the strain of steel; $t_{e} = 0.725 \sqrt[3]{b A}$ , $A = 2 b r$ , $b$ and $r$ are the thickness of protection layer and spacing of steel respectively.

The bond of lining and steel is a complex interaction to transfer stress and deformation. The damage variable $D$ is used to describe the change of material mechanical characteristic. When $D$ is more than zero, the lining will partly lose its force-bearing capacity, and the steel will share partial load. When cracks occur in lining, the load will be redistributed, and the strain of steel can be described by the maximum principal strain of lining element:

13

ε_{s} = \frac{(1 - ρ + n ρ) D ε}{(η - η ρ + n ρ D)},

where $n = E_{r} / E_{d}$ is the elasticity modulus ratio of steel to lining; $ρ$ is the reinforcement ratio; $ε$ is the maximum principal strain of lining element.

According to the Eqs. (12) and (13), the maximal probability crack width of lining elements can be estimated as follows:

14

W_{m a x} = \frac{2.9 \sqrt[3]{b A} (1 - ρ + n ρ) D ε}{(η - η ρ + n ρ D)} .

3. Coupling interaction analysis model of inner water and lining

3.1. Dynamic explicit finite element method of fluid medium

In view of the horizontal seismic vibration in this paper, when performing finite element analysis aimed at the interaction system of inner water and lining structure, it can be assumed that the fluid in the limited computing area is a steady uniform flow, and mainly considering the normal interaction between fluid and structure, which almost meet the basic assumptions of fluid element based on displacement. The dynamic equation of non-viscous compressible fluid (ideal fluid) is obtained:

15

ρ_{1} {\ddot{u}}_{1 i} = K_{1} u_{1 j, i j},

where the subscript 1 is the fluid medium; $ρ_{1}$ means the fluid density, $K_{1}$ means bulk modulus; $\ddot{u}$ and $u$ refer to acceleration and displacement, respectively; $i$ and $j$ are Cartesian coordinate components.

When the explicit element method is combined with multi-transmitting artificial boundary, a very little damping which is proportional to the stiffness matrix is introduced to eliminate the high frequency instability problem [28]. Meanwhile, the explicit finite element discretization equation of non-viscous compressible fluid can be written as:

16

m_{1} {\ddot{u}}_{1} + c_{1} {\dot{u}}_{1} + k_{1} u_{1} = f_{1},

where $m_{1}$ , $c_{1}$ and $k_{1}$ respectively refer to mass matrix, damping matrix and stiffness matrix of fluid medium; $f_{1}$ is external load matrix of fluid-lining coupling interface; ${\ddot{u}}_{1}$ , ${\dot{u}}_{1}$ and $u_{1}$ are acceleration, velocity and displacement matrix of fluid medium.

From the viewpoint of stability and calculation workload, the explicit integration scheme based on time domain weighted residuals is adopted. The time domain explicit integral format [29] and the node force expression of fluid-lining coupling interface for dynamic analysis of non-viscous compressible fluid are obtained:

17

u_{1 l i}^{n + 1} = u_{1 l i}^{n} + Δ t {\dot{u}}_{1 l i}^{n} - \frac{1}{2} Δ t m_{1 l}^{- 1} c_{1 l i} (2 Δ t {\dot{u}}_{1}^{n} - u_{1}^{n} + u_{1}^{n - 1}) - \frac{1}{2} Δ t^{2} m_{1 l}^{- 1} (k_{1 l i} u_{1}^{n} - f_{1 l i}^{n}),

18

{\dot{u}}_{1 l i}^{n + 1} = \frac{1}{Δ t} (u_{1 l i}^{n + 1} - u_{1 l i}^{n}) - \frac{1}{2} m_{1 l}^{- 1} c_{1 l i} (u_{1}^{n + 1} - u_{1}^{n}) - \frac{1}{2} Δ t m_{1 l}^{- 1} (k_{1 l i} u_{1}^{n} - f_{1 l i}^{n + 1}),

19

f_{1 l i}^{n + 1} = \frac{2}{Δ t} m_{1 l} [{\dot{u}}_{1 l i}^{n + 1} - \frac{1}{Δ t} (u_{1 l i}^{n + 1} - u_{1 l i}^{n})] + \frac{1}{Δ t} c_{1 l i} (u_{1}^{n + 1} - u_{1}^{n}) + k_{1 l i} u_{1}^{n + 1},

where the subscript $l$ refers to the node number of fluid medium, $i$ is the coordinate direction of node.

3.2. Dynamic explicit finite element method of solid medium

For lining structure and rocks, in the case of no volume force, a dynamic equation of linear elastic solid medium is obtained according to the classic linear elastic dynamics theory:

20

(λ_{2} + μ_{2}) {\ddot{u}}_{2 j, i j} + μ_{2} u_{2 j, i j} = ρ_{2} {\ddot{u}}_{2 i},

where the subscript 2 refers to lining or rock medium; $λ_{2}$ and $μ_{2}$ are Lame constants; $ρ_{2}$ is the density of lining or rock.

Discretizing the Eq. (20) with the explicit finite element method and introducing the Rayleigh damping terms can lead to the following explicit FEM equation of solid medium:

21

m_{2} {\ddot{u}}_{2} + c_{2} {\dot{u}}_{2} + k_{2} u_{2} = f_{2} + R_{2},

where $m_{2}$ , $c_{2}$ and $k_{2}$ respectively refer to mass matrix, damping matrix and stiffness matrix of solid medium; $f_{2}$ and $R_{2}$ respectively refer to the external load matrix of fluid-lining coupling interface and seismic load matrix of boundary.

Likewise, the explicit integration scheme based on time domain weighted residuals is adopted. The displacement and velocity responses of solid medium and node force of fluid-lining coupling interface can be expressed as:

22

u_{2 l i}^{n + 1} = u_{2 l i}^{n} + Δ t {\dot{u}}_{2 l i}^{n} - \frac{1}{2} Δ t m_{2 l}^{- 1} c_{2 l i} (2 Δ t {\dot{u}}_{2}^{n} - u_{2}^{n} + u_{2}^{n - 1})

- \frac{1}{2} Δ t^{2} m_{2 l}^{- 1} (k_{2 l i} u_{2}^{n} - f_{2 l i}^{n} - R_{2 l i}^{n}),

23

{\dot{u}}_{2 l i}^{n + 1} = \frac{1}{Δ t} (u_{2 l i}^{n + 1} - u_{2 l i}^{n}) - \frac{1}{2} m_{2 l}^{- 1} c_{2 l i} (u_{2}^{n + 1} - u_{2}^{n}) - \frac{1}{2} Δ t m_{2 l}^{- 1} (k_{2 l i} u_{2}^{n} - f_{2 l i}^{n} - R_{2 l i}^{n}),

24

f_{2 l i}^{n + 1} = \frac{2}{Δ t} m_{2 l} [{\dot{u}}_{2 l i}^{n + 1} - \frac{1}{Δ t} (u_{2 l i}^{n + 1} - u_{2 l i}^{n})] + \frac{1}{Δ t} c_{2 l i} (u_{2}^{n + 1} - u_{2}^{n}) + k_{2 l i} u_{2}^{n + 1} - R_{2 l i}^{n + 1} .

With regard to fluid and solid medium, the following initial conditions can be introduced:

25

u^{0} = {\dot{u}}^{0} = 0, {\ddot{u}}_{l i}^{0} = m_{l}^{- 1} (f_{l i}^{0} - c_{l i} {\dot{u}}^{0} - k_{l i} u^{0}),

26

u_{l i}^{- 1} = u_{l i}^{0} - Δ t {\dot{u}}_{l i}^{0} + \frac{1}{2} Δ t^{2} {\ddot{u}}_{l i}^{0} .

3.3. Boundary conditions of coupling interface and explicit finite element solution

According to the mechanical properties of structure, the seismic response of diversion tunnel can be regarded as the superposition of conditional seismic response process of inner water and structure, as shown in Fig. 4. In which Fig. 4(a) is the sketch of seismic response of diversion tunnel; Fig. 4(b), 4(c) are respectively sketches of seismic response of structure and inner water; $f_{1}$ and $f_{2}$ are interaction forces of the fluid-lining coupling interface; $F$ and $R_{2}$ refer to seismic input load.

The continuity conditions of interface stresses and displacements are satisfied for ideal fluid medium and lining medium. The boundary conditions can be expressed as follows.

1) The normal stress of interface is continuous:

27

f_{1 s z} + f_{2 s z} = 0 .

2) The tangential stress of interface is zero:

28

f_{1 s k} = f_{2 s k} = 0 .

3) The normal displacement and velocity of interface are continuous:

29

u_{1 s z} = u_{2 s z}, {\dot{u}}_{1 s z} = {\dot{u}}_{2 s z},

where $s$ is the node number of fluid-lining coupling interface; $z$ and $k$ are normal and tangential direction of interface node.

Substituting the boundary condition equations from Eqs. (27) to (29) into the displacement and velocity Eqs. (17), (18) and (22), (23), the normal displacement and velocity of interface node of fluid and lining medium at time step $n + 1$ can be acquired:

30

u_{1 s z}^{n + 1} = u_{1 s z}^{n} + Δ t {\dot{u}}_{1 s z}^{n} - \frac{1}{2} Δ t^{2} (m_{1 s} + m_{2 s})^{- 1} (k_{1 s} u_{1}^{n} + k_{2 s} u_{2}^{n})

- \frac{1}{2} Δ t (m_{1 s} + m_{2 s})^{- 1} [c_{1 s} (2 Δ t {\dot{u}}_{1}^{n} - u_{1}^{n} + u_{1}^{n - 1}) + c_{2 s} (2 Δ t {\dot{u}}_{2}^{n} - u_{2}^{n} + u_{2}^{n - 1})],

u_{2 s z}^{n + 1} = u_{1 s z}^{n + 1},

31

{\dot{u}}_{1 s z}^{n + 1} = \frac{1}{Δ t} (u_{1 s z}^{n + 1} - u_{1 s z}^{n}) - \frac{1}{2} Δ t (m_{1 s} + m_{2 s})^{- 1} (k_{1 s} u_{1}^{n + 1} + k_{2 s} u_{2}^{n + 1})

- \frac{1}{2} Δ t (m_{1 s} + m_{2 s})^{- 1} [c_{1 s} (u_{1}^{n + 1} - u_{1}^{n}) + c_{2 s} (u_{2}^{n + 1} - u_{2}^{n})],

{\dot{u}}_{2 s z}^{n + 1} = {\dot{u}}_{1 s z}^{n + 1} .

Substituting the normal displacement and velocity expressions of interface node into Eqs. (19) and (24), the nodal load of interface at next time step can be solved. Then substituting the nodal load expressions into Eqs. (18) and (23), the internal node velocity of fluid and lining medium at time step $n + 1$ can be obtained.

It should be noted that this paper is just focused on the normal load of fluid-lining coupling interaction rather than the tangential load, and the normal load of interface at the initial time is the internal water pressure under the static condition.

Fig. 4Seismic response analysis model of diversion tunnel

a) Diversion tunnel

b) Lining and rocks

c) Inner water

4. Inner water exosmosis analysis after lining cracks

Under seismic load, concrete lining structure would enter into the damage evolution stage and eventually produce macroscopic cracks, resulting in inner water exosmosis along the cracks [30]. Because the formation time of seepage field is far more than seismic duration, studying the unsteady seepage problem under seismic load lacks practical significance. This paper just accounts for the situation of inner water exosmosis after earthquake.

It could be known from the results of permeation tests [31] that the permeability of concrete material will sharply increase with the formation and development of cracks, and it has a close relationship with the crack width. Because of the limitation of experiment condition, the fissure seepage rule [32] is based on the precondition that the fissures are smooth, flat and straight. In view of the roughness and partial contact of concrete cracks, the permeability along the propagation direction of the main crack can be calculated by the modified cubic law [33, 34]:

32

k_{f} = \frac{γ}{12 μ c} W_{m a x}^{2},

33

c = 1 + 8.8 {(\frac{Δ}{d})}^{1.5},

where $γ$ and $μ$ are the bulk density and dynamic viscosity factor of water respectively; $c$ is the roughness modification coefficient of crack; $Δ$ is the average raised height of granule in crack, which can be determined by tests; $d$ is the hydraulic diameter of crack, which is twice of the crack width.

The concrete with cracks can be regarded as an anisotropic medium. According to Darcy’s law, the basic equation of three-dimension stable seepage can be obtained as follows:

34

\frac{\partial}{\partial X} (k_{x x} \frac{\partial H}{\partial X}) + \frac{\partial}{\partial Y} (k_{y y} \frac{\partial H}{\partial Y}) + \frac{\partial}{\partial Z} (k_{Z Z} \frac{\partial H}{\partial Z}) = Q,

where $k_{x x}$ , $k_{y y}$ , $k_{Z Z}$ are the main permeability coefficient in direction $x$ , $y$ , $z$ ; $Q$ is the seepage quantity of unit volume; $H$ is the hydraulic head.

Based on the principle of FEM, the node hydraulic heads are interpolated through the interpolation principle after the seepage field is discretized. The basic equation of three-dimension stable seepage can be described as follows:

35

[K] \{H\} = \{F\},

where $K$ is the seepage conductive matrix; $F$ is the nodal load matrix which is integrated at the seepage boundary.

Combining with two kinds of boundary conditions for three-dimension stable seepage, the inner water exosmosis considering lining cracks after earthquake can be solved.

5. Engineering application

5.1. Engineering background

Some hydropower station located in western Yunnan province of China is the fifth cascade hydropower plant of upper reaches of Lancang River. The diversion tunnel of the hydropower station spans 7126 m, which goes across different geological tectonic zones. The tunnel section from 2+500.00 to 2+680.00 located in the middle of water diversion system is chosen for analysis in this paper, of which the seismic intensity is up to VII. The tunnel section is designed in circular structure, with the excavation diameter 10.2 m, inner diameter 9.0 m after secondary lining, center elevation of tunnel 32.00 m, buried depth of tunnel 241.67 m and maximum inner head 135.24 m. The reinforced concrete structure is applied in the tunnel with C25 concrete, and the lining thickness is 0.6 m.

Fig. 5Complete finite element model

Fig. 6Local finite element model and monitoring points

Three-dimensional finite element model of diversion tunnel including fluid, lining structure and surrounding rocks is established, which consists of 63307 nodes and 59584 hexahedron elements, and the lining and fluid element numbers are 2688 and 5264 respectively. Besides that, there is FSI of inner water and lining structure, of which the node number is 696. The model ranges along $X$ , $Y$ , and $Z$ -axes are 180.0 m, 180.0 m, and 321.67 m respectively. The complete finite element model of diversion tunnel is shown in Fig. 5, and the local finite element model and layout of monitoring points are shown in Fig. 6. The plane $Y =$ 90.0 m is chosen as the monitoring section, and point A on the vault, point B on the haunch, and point C on the inverted arch of the lining structure of monitoring section are selected as monitoring points to acquire the quantitative dynamic responses under seismic load.

It should be noted that static calculation of three-dimensional finite element should be conducted firstly to simulate the excavation and lining support of diversion tunnel in order to provide the seismic analysis with initial stress and displacement conditions. The mechanical parameters of the rock mass and lining material are shown in Table 1.

Table 1Mechanical parameters of materials

Material	Deformation modulus (GPa)	Cohesion (MPa)	Friction angle (°)	Tensile strength (MPa)	Permeability coefficient (cm∙s^-1)
Rock	10	0.8	45	1.2	5×10^-6
Lining	28	1.8	40	1.27	1×10^-7
Fluid	2.1	–	–	–	–

5.2. Calculation conditions

Free field boundary condition is applied to absorb the reflection waves along the two vertical boundaries which are perpendicular to $X$ -axes, and Viscoelastic artificial boundary condition is used to absorb the incident waves from the bottom of the model[35].

In this study, the acceleration time history curve of the wave of the 1940 EI Centro earthquake (north-south, magnitude $M =$ 6.7, and maximum acceleration = 2.49 m/s²) is adopted. According to the existing codes, the peak value is adjusted to 0.17 $g$ ( $g$ is gravity acceleration), which is equivalent to the seven degree fortification criterion, and duration time is 20 s. The horizontal (direction of $X$ -axes) vibration is considered in this paper, and after pretreatment, the $X$ -acceleration time history curve which is inputted from the bottom of the limited area is shown in Fig. 7.

Fig. 7X-acceleration time history curve of input seismic wave

5.3. Calculation results and analysis

5.3.1. Displacement-time history analysis of lining structure

$X$ -displacement time curves of monitoring points of lining structure are shown in Fig. 8. The conclusions can be obtained: (1) The displacement time curve of the vault, haunch and inverted arch almost coincide; (2) $X$ -displacement reaches the maximum 7.12 cm at $t =$ 2.66 s, and reaches the minimum –6.78 cm at $t =$ 4.08 s; (3) The displacement reaches the wave peak or trough at the same time, indicating that various locations are in a synchronous vibration state.

By comparing with the displacement time curves of monitoring points, it can be known that the displacement of haunch is less than the vault and inverted arch when considering the hydrodynamic pressure, and the maximal differential displacement exceeds 10 mm, which seems small, but the relative deformation of 10 mm can make great influence on the mechanical properties of lining structure.

These conclusions indicate that the haunch of tunnel structure is significantly affected by the hydrodynamic pressure, which impedes the deformation of lining structure.

Fig. 8Displacement-time curves of monitoring points

Fig. 9Maximum principal stress-time curves of monitoring points

5.3.2. Principal stress-time history analysis of lining structure

Since the compressive strength of concrete is far more than the tensile strength, the damage and failure modes of lining structure under seismic load mainly refer to tensile failure. This paper is focused on the change laws of maximum principal stresses of monitoring points. The maximum principal stress-time curves of monitoring points of lining structure are shown in Fig. 9, which indicate the following: (1) The maximum principal stress of lining structure is up to 0.8-0.9 MPa under the initial hydrostatic pressure and seismic load; (2) The maximum principal stress quickly increases, and gradually exceeds the tensile strength of concrete, in which the haunch of lining structure is remarkably influenced, leading to great stress increment; (3) The hydrodynamic pressure has a significant influence on the stress state of lining structure, and it is of great importance in the anti-earthquake design of tunnel to consider the coupling interaction of inner water and lining structure under seismic load.

5.3.3. Damage and cracking analysis of lining structure

The dynamic damage and cracking model of concrete lining structure proposed in this paper is applied to simulate the damage and cracking process under seismic loading, and damage coefficient and maximal crack width distribution of lining structure at the monitoring section are shown in Fig. 10 and Fig. 11 (unit: mm).

It could be seen from Fig. 10 that under the earthquake and hydrodynamic pressure, the maximum principal stress of lining structure gradually exceeds the tensile strength of concrete, and damage and cracking of lining structure occur and develop gradually from the haunch to each side. Since lining structure mainly refers to elastic deformation under seismic load, hydrodynamic pressure becomes the key factor to cause the damage and cracking. Because the haunch of tunnel structure is remarkably influenced by earthquake and hydrodynamic pressure, causing serious damage, which gradually expands to both sides. In addition, the damage state of inner lining is much more serious than the outer layer, which indicates the parts with direct contact to inner water suffer more serious damage.

We can come to the conclusions from Fig. 11 that with the lining structure entering into damage evolution stage, cracking failure produces accordingly, and the maximal crack width distribution law is consistent with the damage coefficient. The maximal crack width at the monitoring section is up to 0.25-0.5 mm, in which the haunch has the maximal crack width, up to 0.5 mm, and expanding to the vault and inverted arch. By comparing the crack width of inner and outer lining structure, it can be seen that inner lining is significantly greater than the outer for the serious influence of hydrodynamic pressure, which reveals that the direct interaction with inner water would make great influence on lining structure, and it is the weak link of seismic design for diversion tunnel.

Fig. 10Damage coefficient distribution of lining

Fig. 11Maximal crack width distribution of lining

5.3.4. Inner water exosmosis analysis after lining cracks

According to the three-dimensional seepage field analysis, the distribution of pressure head isolines of lining structure at the monitoring section is shown in Fig. 12. The isolines are distributed densely in the lining, and the value decreases quickly from the inner layer to the outer. Meanwhile, the distribution of the isolines in the inner layer is sparser than the outer layer, which implies that the crack width of the inner layer is bigger. Consequently, the hydraulic gradient increases from the inner layer to the outer in the lining. Since the lining structure of haunch and inverted arch suffers serious damage, leading to great increment of permeability coefficient, and resulting in the decrease of hydraulic gradient and sparser distribution of pressure head isolines than any other parts.

Furthermore, because of the damage and cracking of lining structure, inner water exosmosis makes the surrounding rocks bear more water pressure. It can be seen from Fig. 12 that the pressure head borne directly by surrounding rocks is up to 63.54 m, which accounts for about 47 % of inner head, and would greatly affect the stability of surrounding rocks.

Fig. 12Distribution of pressure head isolines of lining structure at the monitoring section

6. Conclusions

Based on the seismic response mechanism of concrete lining structure in high pressure diversion tunnel, a dynamic explicit finite element numerical simulation method considering the coupling interaction of lining and inner water is presented. Taking the dynamic damage and cracking process of concrete lining structure under seismic load into account, the inner water exosmosis process with lining cracks is analyzed. A diversion tunnel of some hydropower station is chosen as an engineering project to analyze the seismic response. We can reach the following conclusions.

1) The displacement time curves of monitoring points are consistent with the input wave, and various parts of diversion tunnel are in a synchronous vibration state. The displacement of haunch is obviously less than the vault and inverted arch due to the influence of hydrodynamic pressure.

2) Under seismic load and hydrodynamic pressure, the maximum principal stress quickly increases, and gradually exceeds the tensile strength of concrete, in which the haunch of lining structure is remarkably influenced, leading to great stress increment and serious damage.

3) The lining structure mainly enters into the stage of damage and cracking, in which the haunch is greatly affected, resulting in greater damage coefficient and crack width, and gradually extension to each side. The inner lining directly interacting with fluid is the weak link of seismic design for diversion tunnel.

4) Inner water exosmosis occurs after lining cracks. Affected by the distribution of crack width, hydraulic gradient is smaller and the distribution of pressure head isolines is sparser on the haunch. Inner water exosmosis makes the surrounding rocks bear more water pressure, and it would affect the stability of surrounding rocks.

References

Kazuhide Y., Yoshiyuki K. Historical earthquake damage to tunnels in Japan and case studies of railway tunnels in the 2004 Niigataken-Chuetsu Earthquake. QR of RTRI, Vol. 48, Issue 4, 2007, p. 136-141.

Publisher
Wang W. L., Wang T. T., Su J. J., et al. Assessment of damage in mountain tunnels due to the Taiwan Chi-Chi Earthquake. Tunneling and Underground Space Technology, Vol. 16, 2001, p. 133-150.

Publisher
Li T. B. Damage to mountain tunnels related to the Wenchuan earthquake and some suggestions for a seismic tunnel construction. Bulletin of Engineering Geology and the Environment, Vol. 71, 2012, p. 297-308.

Publisher
Wang Z. Z., Gao B., Jiang Y. J., et al. Investigation and assessment on mountain tunnels and geotechnical damage after the Wenchuan earthquake. Technological Sciences, Vol. 25, Issue 2, 2009, p. 546-558.

Publisher
Tedesco J. W., Powell J. C., Allen R. C., et al. A strain-rate-dependent concrete material model for ADINA. Computers and Structures, Vol. 64, Issue 5, 1997, p. 1053-1067.

Publisher
Winnicki A., Cichon C. Plastic model for concrete in plane stress state. Journal of Engineering Mechanics, Vol. 124, Issue 6, 1998, p. 591-602.

Publisher
Xiao S. Y., Li H. N., Du R. Q., et al. Four-parameter dynamic constitutive model of concrete. Journal of Harbin Institute of Technology, Vol. 38, Issue 10, 2006, p. 1754-1757.

Search CrossRef
Pandey A. K., Kumar R., Paul D. K., et al. Strain rate model for dynamic analysis of reinforced concrete structures. Journal of Structural Engineering, Vol. 132, Issue 9, 2006, p. 1393-1401.

Publisher
Shkolnik I. K. Influence of high strain rates on stress-strain relationship, strength and elastic modulus of concrete. Cement and Concrete Composites, Vol. 30, 2008, p. 1000-1012.

Publisher
Li Q. B., Zhang C. H., Wang G. L. Dynamic damage constitutive model of concrete in uniaxial tension and compression. Journal of Wuhan Institute of Hydraulic and Electric Engineering, Vol. 12, 1994, p. 55-60.

Search CrossRef
Wang J. T., Du X. L., Zhao C. G., et al. Explicit finite element scheme for dynamic analyses of wave propagation in multi-layered media. Chinese Journal of Rock Mechanics and Engineering, Vol. 22, Issue 1, 2003, p. 75-79.

Search CrossRef
Wang J. T., Du X. L., Zhao C. G. Explicit finite element method for dynamic analyses of fluid-saturated porous solid. Chinese Journal of Rock Mechanics and Engineering, Vol. 21, Issue 8, 2002, p. 1119-1103.

Search CrossRef
Wang T. T., Yang C. H., Yan X. Z., et al. Dynamic response of underground gas storage salt cavern under seismic loads. Tunnelling and Underground Space Technology, Vol. 43, 2014, p. 241-252.

Publisher
Hashemi S., Saadatpour M. M., Kianoush M. R. Dynamic behavior of flexible rectangular fluid containers. Thin-Walled Structures, Vol. 66, 2013, p. 23-28.

Publisher
Du X. L., Wang J. T. Seismic response analysis of arch dam-water-rock foundation systems. Earthquake Engineering and Engineering Vibration, Vol. 3, Issue 2, 2004, p. 283-288.

Publisher
Bouaanani N., Renaud S. Effects of fluid–structure interaction modeling assumptions on seismic floor acceleration demands within gravity dams. Engineering Structure, Vol. 67, 2014, p. 1-8.

Publisher
Nateghi R., Kiany M., Gholipouri O. Control negative effects of blasting waves on concrete of the structures by analyzing of parameters of ground vibration. Tunnelling and Underground Space Technology, Vol. 24, Issue 6, 2009, p. 508-515.

Publisher
Cheng X. S., Xu W. W., Yue C. Q., et al. Seismic response of fluid-structure interaction of undersea tunnel during bidirectional earthquake. Ocean Engineering, Vol. 75, 2014, p. 64-70.

Publisher
Song Y. P. Dynamic Constitutive Relation and Failure Criterion for Concrete. Academic Press, Beijing, China, 2013, p. 240-288.

Search CrossRef
Zhou X. Q., Hao H. Modelling of compressive behavior of concrete-like materials at high strain rate. International Journal of Solids and Structures, Vol. 45, 2008, p. 4645-4661.

Publisher
Bischoff P. H., Perry S. H. Compressive behavior of concrete at high strain rates. Materials and Structures, Vol. 24, 1991, p. 425-450.

Publisher
Brooks J. J., Al-Samarie N. H. Influence of Rate of Strength on Tensile Stress-Strain Behavior Concrete. Fracture of Concrete and Rock, Elsevier Science Publisher, 1989.

Search CrossRef
Chen J. Y., Hu Z. Q., Lin G. Seismic analysis of arch dam with joints based on a new strain-rate-dependent plastic damage model. Chinese Journal of Computational Mechanics, Vol. 24, Issue 1, 2004, p. 45-49.

Search CrossRef
Zhang Z. G., Xiao M., Chen J. T. Simulation of earthquake disaster process of large-scale underground caverns using three-dimensional dynamic finite element method. Chinese Journal of Rock Mechanics and Engineering, Vol. 30, Issue 3, 2011, p. 509-522.

Search CrossRef
Su K., Wu H. G. Nonlinear finite element analysis of reinforced concrete lining of hydraulic tunnels. Rock and Soil Mechanics, Vol. 26, Issue 9, 2005, p. 1485-1490.

Search CrossRef
Yuan J. L., Xiao M., Fu Z. H. Numerical analysis on lining structure of restricted cylindric surge chamber in rock masses with steep obliquity. Chinese Journal of Rock Mechanics and Engineering, Vol. 23, Issue 8, 2004, p. 1296-1300.

Search CrossRef
Chen J. T., Yang Y., Ye C., et al. Three-dimensional numerical analysis of compound lining in complex underground surge-shaft structure. Mathematical Problems in Engineering, 2015, p. 387379.

Publisher
Liao Z. P. An Introduction to Wave Motion Theories in Engineering. Academic Press, Beijing, China, 2002, p. 65-322.

Search CrossRef
Wang J. T., Du X. L. An explicit difference method for dynamic analysis of a structure system with damping. Engineering Mechanics, Vol. 19, Issue 3, 2002, p. 109-111.

Search CrossRef
Zhao Y. L., Cao P., Wang Y. X., et al. Coupling model of seepage-damage-fracture in fractured rock masses and its application. Chinese Journal of Rock Mechanics and Engineering, Vol. 27, Issue 8, 2008, p. 1634-1643.

Search CrossRef
Gerard B., Breysse D., Ammouche A., et al. Cracking and permeability of concrete under tension. Materials and Structures, Vol. 29, Issue 3, 1996, p. 141-151.

Publisher
Xiong X. B., Zhang C. H., Wang E. Z. A review of steady state seepage in a single fracture of rock. Chinese Journal of Rock Mechanics and Engineering, Vol. 28, Issue 9, 2009, p. 1839-1847.

Search CrossRef
Wu Y. Q. Types of rock mass structure and hydromechanical models in rock mass. Chinese Journal of Rock Mechanics and Engineering, Vol. 19, Issue 6, 2000, p. 687-691.

Search CrossRef
Bian K., Xiao M., Chen J. T. Study on coupled seepage and stress fields in the concrete lining of the underground pipe with high water pressure. Tunnelling and Underground Space Technology, Vol. 24, 2009, p. 287-295.

Publisher
Zhang Z. G., Chen J. T., Xiao M. Artificial boundary setting for dynamic time-history analysis of deep buried underground caverns in earthquake disaster. Disaster Advances, Vol. 5, Issue 4, 2012, p. 1136-1142.

Search CrossRef

Cited by

Using an Improved Ant Lion Optimizer Algorithm and Wireless Network to Explore Nondestructive Testing Technology for Water Conservancy Engineering

Hui Shang

(2021)

Numerical Simulation Analysis of Combined Seismic Response for Rock‐Lining‐Water in Hydraulic Tunnel

Guoqing Liu | Yanhong Zhang | Ming Xiao | Guangjian Liu

(2021)

Performance Investigation of Tunnel Lining with Cavities around Surrounding Rocks

(2020)

About this article

Received

17 August 2015

Accepted

13 December 2015

Published

31 March 2016

SUBJECTS

Seismic engineering

DOI

https://doi.org/10.21595/jve.2015.16289

Keywords

diversion tunnel

dynamic damage and cracking

fluid-structure coupling

explicit finite element

permeability coefficient

inner water exosmosis

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.