# Unified Description of the Mechanical Properties of Typical Marine Soil and Its Application

Authors: Yongqiang Li; Gang Wang; Liping Jing; Lei Zhang; Xinjun Cheng
Journal: Shock and Vibration
 (2017)

---

## Abstract

This study employed a modified elastoplastic constitutive model that can systematically describe the monotonic and cyclic mechanical behaviors of typical marine soils combining the subloading, normal, and superloading yield surfaces, in the seismic response analysis of three-dimensional (3D) marine site. New evolution equations for stress-induced anisotropy development and the change in the overconsolidation of soils were proposed. This model can describe the unified behaviour of unsaturated soil and saturated soil using independent state variables and can uniquely describe the multiple mechanical properties of soils under general stress states, without changing the parameter values using the transform stress method. An effective stress-based, fully coupled, explicit finite element–finite difference method was established based on this model and three-phase field theory. A finite deformation analysis was presented by introducing the Green-Naghdi rate tensor. The simulation and analysis indicated that the proposed method was sufficient for simulating the seismic disaster process of 3D marine sites. The results suggested that the ground motion intensity would increase due to the local uneven complex topography and site effect and also provided the temporal and spatial distribution of landslide and collapse at the specific location of the marine site.


---

## Body

## 1. Introduction

The seismic safety of marine sites in strong earthquake is receiving increased attention because of the potential serious hazard of earthquakes to offshore engineering. Three-dimensional (3D) elastoplastic seismic response analysis of marine sites is a focus of and frontier topic in geotechnical earthquake engineering. The investigation and analysis of the typical earthquake damage phenomena and failure characteristics of 3D marine sites in large earthquake are of great significance.In the seismic response analysis of 3D marine sites, one of the key points is to build a constitutive model that can suitably describe the mechanical properties of marine soils; in particular, a reasonable description of the dynamic characteristics of marine soils is the most important task. Considering the coexistence of unsaturated and saturated soils on an island, developing a unified constitutive model that can describe unsaturated and saturated soil continuously is difficult.The mechanical behaviours of soil are dependent on the state conditions of the soil (overconsolidation (influence of density), structure (whether the soil is disturbed or not, Asaoka et al. [1]), anisotropy (mainly stress induced), and degree of saturation) as well as the loading conditions (drained or undrained condition, monotonic or cyclic loading, and other loading paths aside from the conventional triaxial condition). Hashiguchi and Ueno [2] proposed the concept of “subloading,” which made the description of soil overconsolidation possible. The concept of “superloading,” proposed by Asaoka et al. [1], together with the concept of subloading, not only allows the description of overconsolidation but also explains the effect of the soil structure. Zhang et al. [3] introduced a new approach to describe the stress-induced anisotropy. They noted that the change in density is influenced by not only plastic stretching and elastic unloading but also stress-induced anisotropy. Based on their model, the mechanical behaviour of soil subjected to different loadings under different drainage conditions was simulated. The above constitutive models are all based on the Cam-Clay model (Roscoe et al. [4]) or the modified Cam-Clay model (Schofield and Wroth [5]). Therefore, the physical meaning of the model is easy to understand, and relatively few parameters need to be employed. Based on these models, this study establishes a more widely applicable constitutive model.Considerable research has been performed experimentally, empirically, and theoretically on the mechanical behaviours of unsaturated soils, which are widely found on islands. Several constitutive models for unsaturated soils have been developed since the work of Alonso et al. [6]. These researchers proposed the Barcelona basic model (BBM) using LC conception, and this model is typically regarded as one of the basic models for unsaturated soils. Constitutive models based on the framework of the BBM have been presented by Kohgo et al. [7, 8], Cui and Delage [9], and Sun et al. [10], among others. In these models, the stress-suction-strain relations of unsaturated soil were considered explicitly, whereas the degree of saturation was not considered directly. Constitutive models considering the influence of the degree of saturation were also proposed in the same period by Kato et al. [11], Gallipoli et al. [12], Muraleetharan et al. [13], and Sheng et al. [14].The majority of the models in the literature use stress and suction as independent state variables. Ohno et al. [15] made an assumption that the void ratio-logarithmic effective stress relation e-lnp is dependent on the degree of saturation based on various experimental results and modelled the dependency with a simple relation. In the model proposed by Ohno et al. [15], only the effective stress, which was defined by a combination of net stress, suction, and “effective degree of saturation,” was used as a state variable, and a Cam-Clay-type model for unsaturated soil was proposed based on the effective stress. The advantage of using effective stress or skeleton stress instead of net stress and suction as independent state variables in the modelling of unsaturated soils is that it is considerably easier and smoother to describe the behaviours of soil from an unsaturated state to a saturated state and vice versa if the moisture characteristics of the soil are properly described and incorporated into the constitutive model. Based on this type of model, a numerical analysis related to boundary value problems of saturated and unsaturated grounds can be performed in a unique manner. Pioneering research on this topic was conducted by Kanazawa et al. [16].Zhang and Ikariya [17] proposed a constitutive model for unsaturated soil that used skeleton stress and degree of saturation as independent state variables and can accurately describe the influence of the degree of saturation. The constitutive model described the behaviour of unsaturated soil as well as saturated soil, as the skeleton stress can smoothly shift to effective stress if the degree of saturation changes from an unsaturated condition to a saturated condition. Overconsolidation, one of the main features of soils that were discussed in the models for saturated soils (Zhang et al. [3]), was also considered together with the degree of saturation. Research on unsaturated soil considering the overconsolidation effect for saturated soils includes the work by Kohgo et al. [8] which employs the subloading concept and the work by Russell and Khalili [18] regarding the concept of the bounding surface.Other mechanical features, such as the soil structure formed during the sedimentary process of the soil (Asaoka et al. [1]) and stress-induced anisotropy (Zhang et al. [3]), are often discussed for saturated soil but are rarely incorporated into the existing models. These factors can be incorporated into the proposed model within the framework of Zhang and Ikariya [17]. The primary objectives of this study are to (1) modify the constitutive model considering the structure and stress-induced anisotropy of marine soils and propose new evolution equations for the development of stress-induced anisotropy and the change in overconsolidation; (2) make the above model uniquely describe the multiple mechanical properties of soils under general 3D stress states without changing the parameter values; (3) introduce the improved model into soil-water-air three-phase coupled field theory, conduct an effective stress-based, fully coupled, explicit finite element–finite difference method, (4) perform an elastoplastic seismic response analysis of a 3D typical marine site, and discuss the local uneven complex topography and site effect (Borcherdt [19] and Anderson et al. [20]). These results will help bridge information gaps regarding the seismic safety of marine sites.


## 2. Unified Description of Soil Behaviour

### 2.1. State Line

Based on the experimental results, a quantitative relation for the void ratio-logarithmic mean skeleton stress e-lnp is established using the degree of saturation as a state variable. It is assumed that a normally consolidated line in the unsaturated state (N.C.L.S.) is parallel to the normally consolidated line in the saturated state (N.C.L.) but in a higher position than the N.C.L. [17], which means that unsaturated soil can maintain higher void ratio than saturated soil under the same mean skeleton stress. It is also assumed that the critical state line in the unsaturated state (C.S.L.S.) is parallel to critical state line in the saturated state (C.S.L.) (Roscoe et al. [4]) but in a higher position than the C.S.L.The N.C.L.S. and C.S.L.S. are expressed as follows:(1)N.C.L.S.:e=NSr-λln⁡pprη=qp=0,C.S.L.S.:e=ΓSr-λln⁡pprη=qp=M,where N(Sr) and Γ(Sr) are the void ratios of the N.C.L.S. and C.S.L.S.,with a reference mean skeleton stress pr and a certain degree of saturation Sr. λ is the compression index and has the same value for the saturated and unsaturated states. p=σii/3 and q=3(σij-pδij)(σij-pδij)/2 are the mean skeleton stress and second invariant of the deviatory skeleton stress tensor, respectively. η is the stress ratio, and M is the stress ratio at the critical state and has the same value for the saturated and unsaturated states (Cui and Delage [9]).


### 2.2. Yield Surface

Based on the subloading (Hashiguchi and Ueno [2]), normal, and superloading (Asaoka et al. [1]) yield surfaces in the p′-q plane (shown in Figure 1) and the cyclic mobility model (Zhang et al. [3], Zhang et al. [21], and Ye et al. [22]), a new yielding function can be written as shown in (2) (Li [23]). Using the transform stress method (Yao et al. [24]), this model can uniquely describe the multiple mechanical properties of soils under general stress states, without changing the parameter values. The model is suitable for 3D static-dynamic analysis for complex 3D marine sites.(2)f=ln⁡σ~mσ~m0+ln⁡M2-ς~2+η~∗2M2-ς~2-ρ~s1+e01Cp-ln⁡R~+ln⁡R~∗-εvpCp=0,where σ~m=p~m∗′, σ~m0 is the reference mean skeleton stress; ρ~s is the state variable of the degree of saturation (Zhang and Ikariya [17]); e0 is the void ratio at the reference mean skeleton stress; εvp is the plastic volume strain. The similarity ratio of the normal yield surface to the superloading surface R~∗ and the similarity ratio of the subloading surface to the superloading surface R~ in the transform stress space are given as(3)R~∗=p~′p-′=q~q-,0<R∗≤1,q~p~′=q-p-′,R~=p∗′p-′=q∗q-,0<R≤1,R=1OCR,q-p-′=q~p~′=q∗p∗′,where (p∗′,q∗), (p~′,q~), and (p-′,q-) represent the present stress state, the corresponding normally consolidated stress state, and the structured stress state in the transform stress p~′-q~ plane, respectively. The variables involved in (2) and (3) are defined as(4)η~∗=32η^ijη^ij,η^ij=η~ij-βij,η~ij=s~ijσ~m,ς~=32βijβij,Cp=λ-κ1+e0,where s~ij is the deviatory stress tensor, βij is the anisotropic stress tensor, and κ is the swelling index.Figure 1
Yield surfaces in the transform stress space in the present model.


### 2.3. Consistency Equation for the Yield Surface

Considering the advantages of the associated flow rule (Hashiguchi and Chen [25] and Zhang et al. [3]), the consistency equation for the new yield surface can be given as(5)df=∂f∂σijdσij+∂f∂βijdβij-dρ~s1+e01Cp+1R~∗dR~∗-1R~dR~-1Cpdεvp=0.In (5), evolution equations must be provided for the development of the state variables, namely, the anisotropic stress tensor, degree of structure, degree of overconsolidation, and degree of saturation.


#### 2.3.1. Evolution Rule for the Anisotropic Stress Tensor

Different from the work by Hashiguchi and Chen [25] and Asaoka et al. [26], the following evolution rule is used for the anisotropic stress tensor here:(6)dβij=Λ6MbrblM-ς~η^ijCpM2-ς~2+η~∗2σ~m,where br is the evolution parameter of anisotropy and bl is a limitation parameter, the value of which varies from 0.90 to 0.99 according to the definition of “completely liquefaction” in the engineering sense. A value of bl that is closer to 1.0 indicates that the mean effective stress will be closer to zero in the cyclic mobility region. Without losing physical meaning, we fixed the value of bl at 0.95 throughout this paper. In (6), the artificial limitation on the development of anisotropy originally proposed by Hashiguchi and Chen [25] is no longer necessary. The increment in anisotropy, which plays a critical role in the evolution rule for the overconsolidation, will be discussed in detail later.The increment of the evolution rule for the anisotropic stress tensor can be written as(7)∂f∂βijdβij=ΛM6brblM-ζ~η~∗2-2M2+3η~ijβijCpM2-ς~2+η~∗22M2-ς~2σ~m,where ζ<M provides a natural physical limitation on the development of anisotropy automatically; according to (6), the development of anisotropy will stop at the state when η=β.


#### 2.3.2. Evolution Rule for the Degree of the Structure

The following evolution rule for the degree of the structureR∗, which is the same as that proposed by Asaoka et al. [1, 26, 27], is adopted:(8)dR~∗=Λ2aMR~∗1-R~∗η^CpM2-ζ~2+η~∗2σ~m,where a is a parameter that controls the rate of the collapse of the structure during shearing.


#### 2.3.3. Evolution Rule for the Degree of Overconsolidation

In the present model, the changing rate of overconsolidation is assumed to be controlled by two factors, namely, the plastic component of stretching, which was employed as the only factor in the work by Asaoka et al. [26], and the increment in anisotropy. For degree of overconsolidation R~, the new definition is as follows:(9)dR~=ΛmMln⁡R~τ6η~∗2+1/3M2-η~22CpM2-ς~2+η~∗2σ~mσ~m/σ~m02σ~m/σ~m02+12+R~η~M∂f∂βijdβij·γ,where m is the degradation parameter of overconsolidation and τ and γ are two fitting parameters. In general, τ=1.0 and γ=1.0. dβij is proportional to the norm of the plastic component of stretching dεijp.The term(σ~m/σ~m0)2/(σ~m/σ~m0)2+12 is used because, in the region close to the zero point, the effective mean stress is extremely small and the change in the effective stress is slow compared to the rapid development of plastic strain. Therefore, the changing rate of overconsolidation must be reduced. When the effective mean stress is extremely large, that is, (σ~m/σ~m0)2/(σ~m/σ~m0)2+12=1.0, the change of the effective stress is also slow.


#### 2.3.4. Evolution Rule for the Degree of Saturation

The mechanical behaviours of unsaturated soil are strongly influenced by the hydraulic characteristics, which are expressed as the relationship between the degree of saturation and suction. A proper model of the moisture characteristic curve (MCC) should be proposed in advance to describe the mechanical behaviour of unsaturated soil, particularly the transition from the unsaturated state to the saturated state. Zhang and Ikariya [17] proposed a new MCC model that can properly consider hydraulic hysteresis using the boundary surface concept. For convenience, we will employ this MCC model.The evolution rule for degree of saturation is(10)NSr=N+Nr-NSrs-SrrSrs-Sr;Nr=NSrr,ρ~s=NSr-N=QSrs-Sr;Q=Nr-NSrs-Srr,dρ~s=-QdSr,where Srr and Srs are the degrees of saturation of the residual and saturated conditions, respectively. Equation (10) indicates that N(Sr) changes linearly with the degree of saturation. Nr is the void ratios at the N.C.L.S. under the reference mean skeleton stress Pr when the degrees of saturation is in the residual state; that is, Nr=N(Srr).Due to the associated flow rule(11)dεijp=Λ∂f∂σ~ij=0,the plastic volumetric strain increment and plastic shear strain increment can be expressed as(12)dεvp=Λ∂f∂σ~m=ΛM2-η~2M2-ς~2+η~∗2σ~m,dεdp=Λ∂f∂s~ij=Λ3η^ijM2-ς~2+η~∗2σ~m.The following new relation can be obtained by substituting (7), (8), (9), (10), and (12) into (5):(13)∂f∂σijdσij+1-η~M·γ6MbrblM-ζ~η∗2-2M2+3η~ijβijCpM2-ζ~2+η~∗22M2-ζ~2σ~m·Λ+2aM1-R~∗η~∗CpM2-ζ~2+η~∗2σ~m·Λ-1R~mMln⁡R~τ6η~∗2+1/3M2-η~22CpM2-ζ~2+η~∗2σ~mσ~m/σ~m02σ~m/σ~m02+12·Λ+Q1+e0dSr1Cp-1CpM2-η~2M2-ζ2+η~∗2σ~m·Λ=0,which results in the positive valuable Λ that can be written as(14)Λ=∂f/∂σijdσij+Q/1+e0dSr1/Cp1/CpM2-ζ~2+η~∗2σ~mMs2-η~2=∂f/∂σijdσij+Q/1+e0dSr1/Cphp,where(15)hp=1CpM2-ζ~2+η~∗2σ~mMs2-η~2,Ms2=M2+mMln⁡R~τR~σ~m/σ~m02σ~m/σ~m02+126η~∗2+13M2-η~22-2aM1-R~∗η~∗-1-η~M·γ6MbrblM-ζ~η~∗2-2M2+3η~ijβijM2-ζ~2+η~∗22M2-ζ~2.The strain increment can be divided into elastic and plastic parts as follows:(16)dεij=dεije+dεijp.Using Hooke’s theory with stiffness tensorEijkl, the incremental stress tensor can be expressed as(17)dσij=Eijkldεkl-dεklp=Eijkldεkl-EijklΛ∂f∂σkl.Then, we obtain(18)Λ=∂f/∂σijEijkldεkl+Q/1+e0dSr1/Cp∂f/∂σijEijkl∂f/∂σ~kl+hp.We define the loading criteria as(19)Λ>0loadingΛ=0neutralΛ<0unloading.The denominator of (18) is typically positive (Asaoka et al. [28] and Asaoka et al. [26]); therefore, Λ>0 is equivalent to the following relation:(20)∂f∂σijEijkldεkl+Q1+e0dSr1Cp>0.Substituting (18) into (11) yields(21)dεijp=∂f/∂σijEijkldεkl+Q/1+e0dSr1/Cp∂f/∂σijEijkl∂f/∂σ~kl+hp∂f∂σ~ij.Thus, we obtain(22)dσij=Eijkldεkl-dεklp=Eijkldεkl-EijqrEmnkl∂f∂σmn∂f∂σ~qr1D~dεkl-1CpQ1+e0dSr1D~Eijqr∂f∂σ~qr=Eijkl-Eijklpdεij-1CpQ1+e0dSr1D~Eijqr∂f∂σ~qr=Eijklepdεij-ERFSdSr,where(23)Eijklp=EijqrEmnkl∂f/∂σmn∂f/∂σ~qrD~,D~=∂f∂σmnEmnkl∂f∂σ~kl+hp,(24)Eijklep=Eijkl-Eijklp,ERFS=1CpQ1+e01D~Eijqr∂f∂σ~qr.Equation (24) will be used in the next soil-water-air coupled three-phase field theory.The above improved model can describe the loose to dense soils subjected to monotonic loading under undrained and drained condition. It can also describe the loose to dense soils subjected to cyclic loading under undrained and drained conditions.


## 2.1. State Line

Based on the experimental results, a quantitative relation for the void ratio-logarithmic mean skeleton stress e-lnp is established using the degree of saturation as a state variable. It is assumed that a normally consolidated line in the unsaturated state (N.C.L.S.) is parallel to the normally consolidated line in the saturated state (N.C.L.) but in a higher position than the N.C.L. [17], which means that unsaturated soil can maintain higher void ratio than saturated soil under the same mean skeleton stress. It is also assumed that the critical state line in the unsaturated state (C.S.L.S.) is parallel to critical state line in the saturated state (C.S.L.) (Roscoe et al. [4]) but in a higher position than the C.S.L.The N.C.L.S. and C.S.L.S. are expressed as follows:(1)N.C.L.S.:e=NSr-λln⁡pprη=qp=0,C.S.L.S.:e=ΓSr-λln⁡pprη=qp=M,where N(Sr) and Γ(Sr) are the void ratios of the N.C.L.S. and C.S.L.S.,with a reference mean skeleton stress pr and a certain degree of saturation Sr. λ is the compression index and has the same value for the saturated and unsaturated states. p=σii/3 and q=3(σij-pδij)(σij-pδij)/2 are the mean skeleton stress and second invariant of the deviatory skeleton stress tensor, respectively. η is the stress ratio, and M is the stress ratio at the critical state and has the same value for the saturated and unsaturated states (Cui and Delage [9]).


## 2.2. Yield Surface

Based on the subloading (Hashiguchi and Ueno [2]), normal, and superloading (Asaoka et al. [1]) yield surfaces in the p′-q plane (shown in Figure 1) and the cyclic mobility model (Zhang et al. [3], Zhang et al. [21], and Ye et al. [22]), a new yielding function can be written as shown in (2) (Li [23]). Using the transform stress method (Yao et al. [24]), this model can uniquely describe the multiple mechanical properties of soils under general stress states, without changing the parameter values. The model is suitable for 3D static-dynamic analysis for complex 3D marine sites.(2)f=ln⁡σ~mσ~m0+ln⁡M2-ς~2+η~∗2M2-ς~2-ρ~s1+e01Cp-ln⁡R~+ln⁡R~∗-εvpCp=0,where σ~m=p~m∗′, σ~m0 is the reference mean skeleton stress; ρ~s is the state variable of the degree of saturation (Zhang and Ikariya [17]); e0 is the void ratio at the reference mean skeleton stress; εvp is the plastic volume strain. The similarity ratio of the normal yield surface to the superloading surface R~∗ and the similarity ratio of the subloading surface to the superloading surface R~ in the transform stress space are given as(3)R~∗=p~′p-′=q~q-,0<R∗≤1,q~p~′=q-p-′,R~=p∗′p-′=q∗q-,0<R≤1,R=1OCR,q-p-′=q~p~′=q∗p∗′,where (p∗′,q∗), (p~′,q~), and (p-′,q-) represent the present stress state, the corresponding normally consolidated stress state, and the structured stress state in the transform stress p~′-q~ plane, respectively. The variables involved in (2) and (3) are defined as(4)η~∗=32η^ijη^ij,η^ij=η~ij-βij,η~ij=s~ijσ~m,ς~=32βijβij,Cp=λ-κ1+e0,where s~ij is the deviatory stress tensor, βij is the anisotropic stress tensor, and κ is the swelling index.Figure 1
Yield surfaces in the transform stress space in the present model.


## 2.3. Consistency Equation for the Yield Surface

Considering the advantages of the associated flow rule (Hashiguchi and Chen [25] and Zhang et al. [3]), the consistency equation for the new yield surface can be given as(5)df=∂f∂σijdσij+∂f∂βijdβij-dρ~s1+e01Cp+1R~∗dR~∗-1R~dR~-1Cpdεvp=0.In (5), evolution equations must be provided for the development of the state variables, namely, the anisotropic stress tensor, degree of structure, degree of overconsolidation, and degree of saturation.


### 2.3.1. Evolution Rule for the Anisotropic Stress Tensor

Different from the work by Hashiguchi and Chen [25] and Asaoka et al. [26], the following evolution rule is used for the anisotropic stress tensor here:(6)dβij=Λ6MbrblM-ς~η^ijCpM2-ς~2+η~∗2σ~m,where br is the evolution parameter of anisotropy and bl is a limitation parameter, the value of which varies from 0.90 to 0.99 according to the definition of “completely liquefaction” in the engineering sense. A value of bl that is closer to 1.0 indicates that the mean effective stress will be closer to zero in the cyclic mobility region. Without losing physical meaning, we fixed the value of bl at 0.95 throughout this paper. In (6), the artificial limitation on the development of anisotropy originally proposed by Hashiguchi and Chen [25] is no longer necessary. The increment in anisotropy, which plays a critical role in the evolution rule for the overconsolidation, will be discussed in detail later.The increment of the evolution rule for the anisotropic stress tensor can be written as(7)∂f∂βijdβij=ΛM6brblM-ζ~η~∗2-2M2+3η~ijβijCpM2-ς~2+η~∗22M2-ς~2σ~m,where ζ<M provides a natural physical limitation on the development of anisotropy automatically; according to (6), the development of anisotropy will stop at the state when η=β.


### 2.3.2. Evolution Rule for the Degree of the Structure

The following evolution rule for the degree of the structureR∗, which is the same as that proposed by Asaoka et al. [1, 26, 27], is adopted:(8)dR~∗=Λ2aMR~∗1-R~∗η^CpM2-ζ~2+η~∗2σ~m,where a is a parameter that controls the rate of the collapse of the structure during shearing.


### 2.3.3. Evolution Rule for the Degree of Overconsolidation

In the present model, the changing rate of overconsolidation is assumed to be controlled by two factors, namely, the plastic component of stretching, which was employed as the only factor in the work by Asaoka et al. [26], and the increment in anisotropy. For degree of overconsolidation R~, the new definition is as follows:(9)dR~=ΛmMln⁡R~τ6η~∗2+1/3M2-η~22CpM2-ς~2+η~∗2σ~mσ~m/σ~m02σ~m/σ~m02+12+R~η~M∂f∂βijdβij·γ,where m is the degradation parameter of overconsolidation and τ and γ are two fitting parameters. In general, τ=1.0 and γ=1.0. dβij is proportional to the norm of the plastic component of stretching dεijp.The term(σ~m/σ~m0)2/(σ~m/σ~m0)2+12 is used because, in the region close to the zero point, the effective mean stress is extremely small and the change in the effective stress is slow compared to the rapid development of plastic strain. Therefore, the changing rate of overconsolidation must be reduced. When the effective mean stress is extremely large, that is, (σ~m/σ~m0)2/(σ~m/σ~m0)2+12=1.0, the change of the effective stress is also slow.


### 2.3.4. Evolution Rule for the Degree of Saturation

The mechanical behaviours of unsaturated soil are strongly influenced by the hydraulic characteristics, which are expressed as the relationship between the degree of saturation and suction. A proper model of the moisture characteristic curve (MCC) should be proposed in advance to describe the mechanical behaviour of unsaturated soil, particularly the transition from the unsaturated state to the saturated state. Zhang and Ikariya [17] proposed a new MCC model that can properly consider hydraulic hysteresis using the boundary surface concept. For convenience, we will employ this MCC model.The evolution rule for degree of saturation is(10)NSr=N+Nr-NSrs-SrrSrs-Sr;Nr=NSrr,ρ~s=NSr-N=QSrs-Sr;Q=Nr-NSrs-Srr,dρ~s=-QdSr,where Srr and Srs are the degrees of saturation of the residual and saturated conditions, respectively. Equation (10) indicates that N(Sr) changes linearly with the degree of saturation. Nr is the void ratios at the N.C.L.S. under the reference mean skeleton stress Pr when the degrees of saturation is in the residual state; that is, Nr=N(Srr).Due to the associated flow rule(11)dεijp=Λ∂f∂σ~ij=0,the plastic volumetric strain increment and plastic shear strain increment can be expressed as(12)dεvp=Λ∂f∂σ~m=ΛM2-η~2M2-ς~2+η~∗2σ~m,dεdp=Λ∂f∂s~ij=Λ3η^ijM2-ς~2+η~∗2σ~m.The following new relation can be obtained by substituting (7), (8), (9), (10), and (12) into (5):(13)∂f∂σijdσij+1-η~M·γ6MbrblM-ζ~η∗2-2M2+3η~ijβijCpM2-ζ~2+η~∗22M2-ζ~2σ~m·Λ+2aM1-R~∗η~∗CpM2-ζ~2+η~∗2σ~m·Λ-1R~mMln⁡R~τ6η~∗2+1/3M2-η~22CpM2-ζ~2+η~∗2σ~mσ~m/σ~m02σ~m/σ~m02+12·Λ+Q1+e0dSr1Cp-1CpM2-η~2M2-ζ2+η~∗2σ~m·Λ=0,which results in the positive valuable Λ that can be written as(14)Λ=∂f/∂σijdσij+Q/1+e0dSr1/Cp1/CpM2-ζ~2+η~∗2σ~mMs2-η~2=∂f/∂σijdσij+Q/1+e0dSr1/Cphp,where(15)hp=1CpM2-ζ~2+η~∗2σ~mMs2-η~2,Ms2=M2+mMln⁡R~τR~σ~m/σ~m02σ~m/σ~m02+126η~∗2+13M2-η~22-2aM1-R~∗η~∗-1-η~M·γ6MbrblM-ζ~η~∗2-2M2+3η~ijβijM2-ζ~2+η~∗22M2-ζ~2.The strain increment can be divided into elastic and plastic parts as follows:(16)dεij=dεije+dεijp.Using Hooke’s theory with stiffness tensorEijkl, the incremental stress tensor can be expressed as(17)dσij=Eijkldεkl-dεklp=Eijkldεkl-EijklΛ∂f∂σkl.Then, we obtain(18)Λ=∂f/∂σijEijkldεkl+Q/1+e0dSr1/Cp∂f/∂σijEijkl∂f/∂σ~kl+hp.We define the loading criteria as(19)Λ>0loadingΛ=0neutralΛ<0unloading.The denominator of (18) is typically positive (Asaoka et al. [28] and Asaoka et al. [26]); therefore, Λ>0 is equivalent to the following relation:(20)∂f∂σijEijkldεkl+Q1+e0dSr1Cp>0.Substituting (18) into (11) yields(21)dεijp=∂f/∂σijEijkldεkl+Q/1+e0dSr1/Cp∂f/∂σijEijkl∂f/∂σ~kl+hp∂f∂σ~ij.Thus, we obtain(22)dσij=Eijkldεkl-dεklp=Eijkldεkl-EijqrEmnkl∂f∂σmn∂f∂σ~qr1D~dεkl-1CpQ1+e0dSr1D~Eijqr∂f∂σ~qr=Eijkl-Eijklpdεij-1CpQ1+e0dSr1D~Eijqr∂f∂σ~qr=Eijklepdεij-ERFSdSr,where(23)Eijklp=EijqrEmnkl∂f/∂σmn∂f/∂σ~qrD~,D~=∂f∂σmnEmnkl∂f∂σ~kl+hp,(24)Eijklep=Eijkl-Eijklp,ERFS=1CpQ1+e01D~Eijqr∂f∂σ~qr.Equation (24) will be used in the next soil-water-air coupled three-phase field theory.The above improved model can describe the loose to dense soils subjected to monotonic loading under undrained and drained condition. It can also describe the loose to dense soils subjected to cyclic loading under undrained and drained conditions.


## 2.3.1. Evolution Rule for the Anisotropic Stress Tensor

Different from the work by Hashiguchi and Chen [25] and Asaoka et al. [26], the following evolution rule is used for the anisotropic stress tensor here:(6)dβij=Λ6MbrblM-ς~η^ijCpM2-ς~2+η~∗2σ~m,where br is the evolution parameter of anisotropy and bl is a limitation parameter, the value of which varies from 0.90 to 0.99 according to the definition of “completely liquefaction” in the engineering sense. A value of bl that is closer to 1.0 indicates that the mean effective stress will be closer to zero in the cyclic mobility region. Without losing physical meaning, we fixed the value of bl at 0.95 throughout this paper. In (6), the artificial limitation on the development of anisotropy originally proposed by Hashiguchi and Chen [25] is no longer necessary. The increment in anisotropy, which plays a critical role in the evolution rule for the overconsolidation, will be discussed in detail later.The increment of the evolution rule for the anisotropic stress tensor can be written as(7)∂f∂βijdβij=ΛM6brblM-ζ~η~∗2-2M2+3η~ijβijCpM2-ς~2+η~∗22M2-ς~2σ~m,where ζ<M provides a natural physical limitation on the development of anisotropy automatically; according to (6), the development of anisotropy will stop at the state when η=β.


## 2.3.2. Evolution Rule for the Degree of the Structure

The following evolution rule for the degree of the structureR∗, which is the same as that proposed by Asaoka et al. [1, 26, 27], is adopted:(8)dR~∗=Λ2aMR~∗1-R~∗η^CpM2-ζ~2+η~∗2σ~m,where a is a parameter that controls the rate of the collapse of the structure during shearing.


## 2.3.3. Evolution Rule for the Degree of Overconsolidation

In the present model, the changing rate of overconsolidation is assumed to be controlled by two factors, namely, the plastic component of stretching, which was employed as the only factor in the work by Asaoka et al. [26], and the increment in anisotropy. For degree of overconsolidation R~, the new definition is as follows:(9)dR~=ΛmMln⁡R~τ6η~∗2+1/3M2-η~22CpM2-ς~2+η~∗2σ~mσ~m/σ~m02σ~m/σ~m02+12+R~η~M∂f∂βijdβij·γ,where m is the degradation parameter of overconsolidation and τ and γ are two fitting parameters. In general, τ=1.0 and γ=1.0. dβij is proportional to the norm of the plastic component of stretching dεijp.The term(σ~m/σ~m0)2/(σ~m/σ~m0)2+12 is used because, in the region close to the zero point, the effective mean stress is extremely small and the change in the effective stress is slow compared to the rapid development of plastic strain. Therefore, the changing rate of overconsolidation must be reduced. When the effective mean stress is extremely large, that is, (σ~m/σ~m0)2/(σ~m/σ~m0)2+12=1.0, the change of the effective stress is also slow.


## 2.3.4. Evolution Rule for the Degree of Saturation

The mechanical behaviours of unsaturated soil are strongly influenced by the hydraulic characteristics, which are expressed as the relationship between the degree of saturation and suction. A proper model of the moisture characteristic curve (MCC) should be proposed in advance to describe the mechanical behaviour of unsaturated soil, particularly the transition from the unsaturated state to the saturated state. Zhang and Ikariya [17] proposed a new MCC model that can properly consider hydraulic hysteresis using the boundary surface concept. For convenience, we will employ this MCC model.The evolution rule for degree of saturation is(10)NSr=N+Nr-NSrs-SrrSrs-Sr;Nr=NSrr,ρ~s=NSr-N=QSrs-Sr;Q=Nr-NSrs-Srr,dρ~s=-QdSr,where Srr and Srs are the degrees of saturation of the residual and saturated conditions, respectively. Equation (10) indicates that N(Sr) changes linearly with the degree of saturation. Nr is the void ratios at the N.C.L.S. under the reference mean skeleton stress Pr when the degrees of saturation is in the residual state; that is, Nr=N(Srr).Due to the associated flow rule(11)dεijp=Λ∂f∂σ~ij=0,the plastic volumetric strain increment and plastic shear strain increment can be expressed as(12)dεvp=Λ∂f∂σ~m=ΛM2-η~2M2-ς~2+η~∗2σ~m,dεdp=Λ∂f∂s~ij=Λ3η^ijM2-ς~2+η~∗2σ~m.The following new relation can be obtained by substituting (7), (8), (9), (10), and (12) into (5):(13)∂f∂σijdσij+1-η~M·γ6MbrblM-ζ~η∗2-2M2+3η~ijβijCpM2-ζ~2+η~∗22M2-ζ~2σ~m·Λ+2aM1-R~∗η~∗CpM2-ζ~2+η~∗2σ~m·Λ-1R~mMln⁡R~τ6η~∗2+1/3M2-η~22CpM2-ζ~2+η~∗2σ~mσ~m/σ~m02σ~m/σ~m02+12·Λ+Q1+e0dSr1Cp-1CpM2-η~2M2-ζ2+η~∗2σ~m·Λ=0,which results in the positive valuable Λ that can be written as(14)Λ=∂f/∂σijdσij+Q/1+e0dSr1/Cp1/CpM2-ζ~2+η~∗2σ~mMs2-η~2=∂f/∂σijdσij+Q/1+e0dSr1/Cphp,where(15)hp=1CpM2-ζ~2+η~∗2σ~mMs2-η~2,Ms2=M2+mMln⁡R~τR~σ~m/σ~m02σ~m/σ~m02+126η~∗2+13M2-η~22-2aM1-R~∗η~∗-1-η~M·γ6MbrblM-ζ~η~∗2-2M2+3η~ijβijM2-ζ~2+η~∗22M2-ζ~2.The strain increment can be divided into elastic and plastic parts as follows:(16)dεij=dεije+dεijp.Using Hooke’s theory with stiffness tensorEijkl, the incremental stress tensor can be expressed as(17)dσij=Eijkldεkl-dεklp=Eijkldεkl-EijklΛ∂f∂σkl.Then, we obtain(18)Λ=∂f/∂σijEijkldεkl+Q/1+e0dSr1/Cp∂f/∂σijEijkl∂f/∂σ~kl+hp.We define the loading criteria as(19)Λ>0loadingΛ=0neutralΛ<0unloading.The denominator of (18) is typically positive (Asaoka et al. [28] and Asaoka et al. [26]); therefore, Λ>0 is equivalent to the following relation:(20)∂f∂σijEijkldεkl+Q1+e0dSr1Cp>0.Substituting (18) into (11) yields(21)dεijp=∂f/∂σijEijkldεkl+Q/1+e0dSr1/Cp∂f/∂σijEijkl∂f/∂σ~kl+hp∂f∂σ~ij.Thus, we obtain(22)dσij=Eijkldεkl-dεklp=Eijkldεkl-EijqrEmnkl∂f∂σmn∂f∂σ~qr1D~dεkl-1CpQ1+e0dSr1D~Eijqr∂f∂σ~qr=Eijkl-Eijklpdεij-1CpQ1+e0dSr1D~Eijqr∂f∂σ~qr=Eijklepdεij-ERFSdSr,where(23)Eijklp=EijqrEmnkl∂f/∂σmn∂f/∂σ~qrD~,D~=∂f∂σmnEmnkl∂f∂σ~kl+hp,(24)Eijklep=Eijkl-Eijklp,ERFS=1CpQ1+e01D~Eijqr∂f∂σ~qr.Equation (24) will be used in the next soil-water-air coupled three-phase field theory.The above improved model can describe the loose to dense soils subjected to monotonic loading under undrained and drained condition. It can also describe the loose to dense soils subjected to cyclic loading under undrained and drained conditions.


## 3. Modified Soil-Water-Air Coupled Three-Phase Field Theory

Based on Zienkiewicz’s theory, Oka et al. [29] proposed a coupled field theory:(I)
Displacements of the solid and pore water pressure are chosen as unknown variables.(II)
The finite element method (FEM) is used to discretize the motion equation for the mixture.(III)
The finite difference method (FDM) is used to discretize the continuity equation for the pore fluid.We modified the above theory and established the equations of three-phase field theory based on a finite deformation algorithm.Based on mixture theory (Darcy [30], Truesdell and Toupin [31], Borja [32], and Fredlund [33]), the balance equation for soil-water-air coupled three-phase is expressed as (Li [23])(25)ρu¨is=∂σij∂xj+ρbi,where ρ is the appearance density of the mixture, σij is the total stress tensor, and bi is body force. ρ can be expressed as(26)ρ=ρ-s+ρ-w+ρ-a=1-nρs+nSrρw+1-Srρa,in which the appearance densities of the solid phase ρ-s, water phase ρ-w, and air phase ρ-a are defined as(27)ρ-s=1-nρs,ρ-w=nSrρw,ρ-a=n1-Srρa,where ρs,ρw, and ρa are the densities of the solid phase, water phase, and air phase, respectively; n is the porosity; Sr is the degree of saturation of soil.Based on Bishop theory, total stress tensor can be written as(28)σij=σij′+pFδij,pF=Srpw+1-Srpa,where σij′ is the effective stress and pF,pw, and pa are the mean pore pressure, pore water pressure, and pore air pressure, respectively.The continuum equation for water is expressed as(29)ρwε¨iis-∂2pw∂xi∂xi+γwkwε˙iisn-1Kwp˙dw+γwkw1SrS˙r=0.The continuum equation for air is expressed as(30)-kaγa∂2pda∂xi∂xi-ρaε¨iis+ε˙iisn-1Kap˙da-11-SrS˙r=0.In (29) and (30), ε˙iis=∂u˙is/∂xi, uis is displacement vector of solid phase; kw and ka are coefficients of permeability for water and air, respectively; Kw is the volume elastic coefficient of the water phase; Ka is the volume elastic coefficient of the air phase; pdw and pda are the excessive pore water pressure and excessive pore air pressure, respectively (Li [23]).After the spatial discretization, (25), (29), and (30) can be represented as (31), (32), and (33), respectively:(31)Mu¨→N+Cu˙→N+KΔu→N-K→SatΔSr+K→vSrpEw+K→v1-SrpEa=F→d-R→d+Q→d,(32)ρwK→vTu¨→N+γwnkwK→vTu˙→N-αwpdEw+∑i=1mαiwpdEiw-A+Fsrwp˙dEw+Fsrwp˙dEa=0,(33)ρaK→vTu¨→N+γankaK→vTu˙→N-αapdEa+∑i=1mαiapdEia-B+Fsrap˙dEa+Fsrap˙dEw=0,where different shape functions are used for displacement and pore pressure (water and air) in the FEM discretization. The variables of displacement are given at the nodes, and the variables of pore pressure (water and air) are given at the gravitational centre and are denoted by pdEw and pdEa. u→N is the nodal displacement vector of all nodes in one element. K→Sat,K→v,pEw,pEa,F→d,R→d,Q→d,K→vT,αw,A,Fsrw,αa,B,m, and Fsra are the calculated parameters (Li [23]).Based on the Newmark-β method, the new weak form formulations for the FEM are as follows:(34)M+γΔtC+βΔt2KSr∣tK→v+K→Satks-1K→v1-Sr∣t-K→Satks-1K→vT-A′+αw+F-srwF-srwK→vTF-sra-B′+αa+F-srau¨→N∣t+ΔtpdE∣t+ΔtwpdE∣t+Δta+0∑i=1mαiwpdEi∣t+Δtw∑i=1mαiapdE∣t+Δta=F→d∣t+Δt-R→d∣t+Δt+Q→d∣t+Δt-K→Satks-1pdE∣ta-pdE∣tw-Cu˙→N∣t+1-γΔtu¨→N∣t-KΔtu˙→N∣t+12-βΔt2u¨→N∣t-EwK→vTu˙→N∣t+1-γΔtu¨→N∣t-A′+F-srwpdE∣tw+F-srwpdE∣ta-EaK→vTu˙→N∣t+1-γΔtu¨→N∣t+F-srapdE∣tw-B′+F-srapdE∣ta,where γ,Δt,β,ks-1,A′,αw,F-srw,F-sra,B′,αa,Ew, and Ea are the calculated parameters (Li [23]).Based on (34), the program is capable of solving repeated soil-water-air coupled static and dynamic boundary value problem using an infinitesimal deformation algorithm or finite deformation algorithm. The constitutive model and computational methods have been verified strictly (Zhang et al. [21], Zhang and Ikariya [17], Ye et al. [22], and Li [23]).


## 4. Application

### 4.1. Numerical Model

A 3D numerical model of a typical submarine site is shown in Figure2. The horizontal dimension of the model is 300×300 m and is divided equally into 100×100 elements in the X and Y directions. The sea level elevation is ±0.000 m, the elevation of the maximum height of the soil at the monitoring point (MP) 15 is −5.000 m, and the elevation of the model bottom is −59.226 m. The vertical dimension of the soil and bedrock (54.226 m) is divided into 16 elements. Therefore, the total number of elements is 160,000. The numbers in circles shown in Figure 2 are the monitoring points for acceleration, velocity, and displacement.Figure 2
3D numerical model and monitoring points for a complex submarine site.As shown in Figure2, three monitoring surfaces (MS) and three monitoring lines (ML) are selected to present the results. MS1, including MPs 0, 1, 2, 3, and 4, is located at the bottom of the model; MS2, including MPs 5, 6, 7, 8, and 9, is located at the top of the bedrock; MS3, including MPs 10, 11, 12, 13, and 14, is located at the top of the model. ML1 is composed of MPs 0, 5, and 10 in a vertical line; ML2 is composed of MPs 15–19, 10, and 20–24; ML3, which is vertical to ML2, is composed of MPs 25–27, 10, and 28–30. MPs 5 and 0 are the projection points of point 10 on the top of the bedrock and the input surface. MPs 1, 6, and 11 are on the same vertical line, and the same situation applies to points (2, 7, and 12), (3, 8, and 13), and (4, 9, and 14).


### 4.2. Material Parameters

Experimental data obtained from static and cyclic triaxial compression tests were used to determine the parameters of the marine soil [1–5, 17]. The material parameters of marine soil involved in our model are listed in Tables 1, 2, and 3. The comparison between the simulated stress-strain-dilatancy relations of soil sample under monotonic/cyclic loading in conventional triaxial tests is shown in Figures 3 and 4, respectively. The bedrock is assumed to be elastic, and its parameters are shown in Table 4.Table 1
Material parameters of the soil.
Parameter
Soil
Compression indexλ
0.026
Swelling indexκ
0.005
Critical state constantRf
2.890
Void ratioN (p′=98 kPa on the N.C.L.)
0.980
Poisson’s ratioν
0.320
Degradation parameter of overconsolidationm
2.200
Degradation parameter of structurem∗
0.100
Evolution parameter of anisotropybr
0.500
Parameter of suctionb
0.500
Void ratioNr (p′=98 kPa on the N.C.L.S.)
1.120Table 2
Initial conditions of the soil.
Parameter
Soil
Initial overconsolidation ratioOCR
2.200
Initial degree of structureR0∗
0.600
Initial anisotropyζ0
0.000
Initial void ratioe0 (p′=98 kPa)
0.900
Initial mean effective stressp′ (kPa)
98.000
Gravity densityγ (kN/m3)
1.940
Volume elasticity coefficient of waterKw (kPa)
2.23
0
e
6
Permeability coefficientk (m/s)
1.000e-7Table 3
Parameters of the moisture characteristics curve of the soil.
Parameter
Soil
Saturated degrees of saturationSrs
0.960
Residual degrees of saturationSrr
0.120
Parameter corresponding to drying AEV (kPa)Sd
320.200
Parameter corresponding to wetting AEV (kPa)Sw
13.600
Initial stiffness of scanning curve (kPa)ks0
77600.000
Parameter of shape functionc1
0.013
Parameter of shape functionc2
0.009
Parameter of shape functionc3
22.000Table 4
Material parameters of the elastic model for bedrock.
Parameter
Bedrock
Linear elastic modulusEc (kPa)
1.375
e
6
Poisson’s ratioν
0.250
Gravity densityγ (kN/m3)
2.200
Permeability coefficientk (m/s)
1.000e-10Figure 3
Stress-strain-dilatancy relations of soil sample in strain-controlled monotonic loading CD tests (confining stress: 98 kPa).Figure 4
Comparison between the calculated and experimental stress-strain relations of soil sample under undrained cyclic loading.


### 4.3. Input Ground Motion

The amplitude, spectral characteristics, and duration of the input ground motion have considerable influences on the seismic response of a marine site. The real three-component ground motions, shown in in Figure5(a), are selected as the input ground motions in the analysis. The original records were filtered and baseline-corrected; the velocity time history before and that after adjustment are shown in Figure 5. As shown in Figure 5, the portion in the blue dashed rectangular box was the selected input time histories. Before inputting the model, the velocity time history should be preliminarily reduced by half and transformed into stress time history (Jing et al. [34]). EW is in the east-west direction (positive: east), NS is in the north-south direction (positive: north), and UD is in the up-down direction (positive: up).Figure 5
Input ground motions before and after adjustment.
(a)
Original records of acceleration
(b)
Acceleration after adjustment
(c)
Velocity
(d)
Displacement


### 4.4. Boundary Condition

For static analysis, the four vertical boundaries were simply supported, and the bottom boundary was fixed. However, in the dynamic analysis, the four vertical boundary conditions were free-field boundaries (shown in Figure6), and the bottom boundary condition was a viscous damping boundary (Lysmer and Kuhlemeyer [35] and Jing et al. [34]).Figure 6
Artificial boundaries for dynamic analysis.


### 4.5. Results of the Static Analysis

The effects of water pressure and the pore water pressure cannot be neglected in the seismic response analysis of marine site. As shown in Figure7(a), the pore water pressure varies with the depth. The initial stress, which was generated by the self-weight stress, is shown in Figure 7(b).Figure 7
Initial pore water pressure and stress states (WLE: water depth above monitoring point 15; IGMSR: input ground motion scaling ratio).
(a)
Pore water pressure (Pa), WLE: +5.000 m, IGMSR: 1
(b)
Z-direction initial stress (Pa), WLE: +5.000 m, IGMSR: 1


### 4.6. Results of the Seismic Response Analysis

#### 4.6.1. Spatial and Temporal Distribution of Displacement

Figure8 shows the contours of displacement at the following dynamic times: 5, 10, 20, 30, 40, 50, 60, and 67 s. The deformation increases from the bottom to top of the model. The residual deformations of the MPs on the top are considerably larger than the corresponding residual deformations at the foot of the hill. The numerical simulation results correspond well with the field survey.Figure 8
Contour of displacement at 5, 10, 20, 30, 40, 50, 60, and 67 s.
(a)
Time: 5 s, WLE: +5.000 m, IGMSR: 1
(b)
Time: 10 s, WLE: +5.000 m, IGMSR: 1
(c)
Time: 20 s, WLE: +5.000 m, IGMSR: 1
(d)
Time: 30 s, WLE: +5.000 m, IGMSR: 1
(e)
Time: 40 s, WLE: +5.000 m, IGMSR: 1
(f)
Time: 50 s, WLE: +5.000 m, IGMSR: 1
(g)
Time: 60 s, WLE: +5.000 m, IGMSR: 1
(h)
Time: 67 s, WLE: +5.000 m, IGMSR: 1The spatial and temporal distribution of the displacements shown in Figure8 illustrate that the slope has a large displacement. In the 3D simulation, we can identify the range of sliding, which is critical to our seismic design and reinforcement of marine slope sites.


#### 4.6.2. Dynamic Time History

Typical acceleration, velocity, and displacement time histories at MPs 0, 5, 10, and 23 are shown in Figures9–11.Figure 9
Acceleration time histories of MPs 0, 5, 10, and 23.Figure 10
Velocity time histories of MPs 0, 5, 10, and 23.Figure 11
Displacement time histories of MPs 0, 5, 10, and 23.The PGD, PGV, and PGA at MPs 0–30 are presented in Table5.Table 5
PGD, PGV, and PGA at MPs 0–30.
MP
PGD (m)
PGV (m/s)
PGA (m/s2)
EW
NS
UD
EW
NS
UD
EW
NS
UD
Input
Value
−0.0956
−0.1160
0.1074
0.1297
−0.1567
−0.1295
−1.368
1.011
0.7304
Time (s)
26.99
27.65
25.68
28.85
33.42
26.58
28.90
23.87
43.20
MP 0
Value
−0.1879
−0.2293
0.2162
−0.2260
−0.2487
−0.2581
−1.697
1.393
1.389
Time (s)
27.08
27.69
25.68
33.84
33.46
26.61
43.73
43.12
43.22
MP 1
Value
−0.1885
−0.2342
0.2116
−0.2199
−0.2601
−0.2475
−1.617
1.806
1.232
Time (s)
27.16
27.72
25.71
33.97
33.62
26.61
32.38
33.68
43.22
MP 2
Value
−0.1870
−0.2288
0.2133
0.2308
0.2436
−0.2577
−2.176
1.868
1.289
Time (s)
27.08
27.69
25.68
29.22
23.20
26.61
43.70
33.86
43.25
MP 3
Value
−0.1871
−0.2273
0.2129
0.2338
−0.2534
−0.2558
−1.967
1.470
1.313
Time (s)
27.08
27.69
25.66
29.22
33.46
26.58
43.73
43.09
43.22
MP 4
Value
−0.1888
−0.2327
0.2113
−0.2253
−0.2604
−0.2531
−1.718
1.789
1.216
Time (s)
27.14
27.72
25.68
33.94
33.62
26.61
45.65
33.65
43.30
MP 5
Value
−0.1882
−0.2306
0.2163
0.2400
−0.2874
−0.2596
−2.825
−1.699
1.517
Time (s)
27.06
27.69
25.68
29.25
33.46
26.61
33.70
43.46
43.22
MP 6
Value
−0.1907
−0.2375
0.2119
0.2504
−0.3139
−0.2546
−3.339
1.848
1.617
Time (s)
27.08
27.72
25.71
28.93
33.52
26.61
28.99
33.57
43.25
MP 7
Value
−0.1877
−0.2302
0.2133
0.2553
−0.2967
−0.2609
−2.427
−1.765
1.567
Time (s)
27.06
27.69
25.66
29.25
33.49
26.61
28.96
33.38
43.25
MP 8
Value
−0.1887
−0.2286
0.2130
0.2452
−0.2912
−0.2580
−2.577
−1.989
1.329
Time (s)
27.06
27.69
25.68
29.22
33.46
26.61
28.93
32.96
43.22
MP 9
Value
−0.1898
−0.2346
0.2116
0.2621
−0.3108
−0.2570
−2.999
2.026
1.499
Time (s)
27.11
27.72
25.68
29.28
33.52
26.63
28.99
23.97
43.25
MP 10
Value
−0.1913
−0.2310
0.2160
0.2786
−0.3143
−0.2595
−3.447
2.907
1.601
Time (s)
27.03
27.69
25.71
28.91
33.46
26.58
31.79
31.21
31.32
MP 11
Value
−0.1926
−0.2399
0.2095
0.2857
−0.3510
−0.2578
−4.231
3.623
1.693
Time (s)
27.08
27.72
25.71
28.93
33.52
26.63
31.85
35.00
43.28
MP 12
Value
−0.1881
−0.2304
0.2121
0.2840
−0.3397
−0.2608
−3.559
−3.691
−1.499
Time (s)
27.06
27.72
25.68
28.91
33.49
26.61
33.07
33.44
43.54
MP 13
Value
−0.1912
−0.2285
0.2122
0.2801
−0.3256
−0.2598
−4.153
3.132
1.358
Time (s)
27.03
27.69
25.71
28.91
33.46
26.61
33.70
34.98
43.22
MP 14
Value
−0.1885
−0.2336
0.2092
0.2906
−0.3538
−0.2583
−3.609
−3.392
1.720
Time (s)
27.08
27.72
25.74
28.96
33.52
−0.2583
33.09
23.76
43.30
MP 15
Value
−0.1937
0.2319
0.2081
0.3368
−0.3853
−0.2615
4.030
4.834
−2.002
Time (s)
27.08
23.70
25.71
28.96
33.54
26.66
34.29
31.29
29.14
MP 16
Value
−0.1928
0.2316
0.2106
0.3025
−0.3563
−0.2619
3.949
3.121
−2.268
Time (s)
27.08
23.68
25.71
29.28
33.54
26.63
31.77
35.03
29.14
MP 17
Value
0.1902
−0.2303
0.2136
0.3015
−0.3749
−0.2607
3.373
−3.682
1.705
Time (s)
33.12
27.72
25.74
28.93
33.52
26.63
33.97
32.11
29.46
MP 18
Value
0.1901
−0.2303
0.2133
0.2817
−0.3365
−0.2564
3.328
3.125
−1.727
Time (s)
33.12
27.72
25.68
28.93
33.49
26.61
34.26
36.07
33.49
MP 19
Value
−0.1858
−0.2299
0.2157
0.2627
−0.3108
−0.2565
−3.979
2.793
1.541
Time (s)
27.03
27.69
25.66
28.91
33.49
26.61
31.79
34.98
43.25
MP 20
Value
−0.1963
−0.2314
0.2155
0.2803
−0.3530
−0.2645
−4.061
−3.610
1.320
Time (s)
27.06
27.72
25.68
28.93
33.49
26.63
31.82
36.15
43.25
MP 21
Value
−0.1999
−0.2320
0.2128
0.2874
−0.3617
−0.2623
−3.993
3.608
1.637
Time (s)
27.06
27.72
25.68
28.93
33.52
26.61
33.73
33.60
43.25
MP 22
Value
−0.2003
−0.2309
0.2117
0.3059
−0.3758
−0.2625
−4.020
−4.031
1.855
Time (s)
27.08
27.72
25.68
28.96
33.52
26.58
48.12
34.07
43.22
MP 23
Value
−0.2002
0.2312
0.2106
0.3152
−0.3705
−0.2663
−4.345
4.241
−2.111
Time (s)
27.08
23.70
25.71
28.96
33.54
26.61
31.85
33.65
33.54
MP 24
Value
−0.1977
−0.2316
0.2068
0.3198
−0.3819
−0.2638
−4.269
−4.001
−2.247
Time (s)
27.08
27.74
25.68
28.96
33.54
26.63
33.76
50.36
29.14
MP 25
Value
−0.1893
−0.2302
0.2140
0.2748
−0.3117
−0.2597
−3.767
−3.104
1.631
Time (s)
27.03
27.69
25.66
28.91
33.46
26.61
31.79
36.12
43.25
MP 26
Value
−0.1896
−0.2306
0.2141
0.2770
−0.3129
−0.2610
4.401
3.233
1.637
Time (s)
27.03
27.69
25.68
28.91
33.46
26.61
34.21
31.24
43.22
MP 27
Value
−0.1888
0.2315
0.2154
0.2707
−0.3054
−0.2569
−3.752
−2.953
1.713
Time (s)
27.03
27.69
25.68
28.91
33.49
26.61
33.04
33.36
43.20
MP 28
Value
−0.1977
−0.2306
0.2156
0.2723
−0.3573
−0.2635
−3.680
−3.178
1.625
Time (s)
27.06
27.72
25.68
28.93
33.49
26.63
33.73
36.17
43.25
MP 29
Value
−0.1929
−0.2295
0.2138
0.2782
−0.3241
−0.2589
−3.652
2.763
−1.735
Time (s)
27.03
27.69
25.68
28.91
33.49
26.61
28.96
34.98
33.46
MP 30
Value
−0.1918
−0.2290
0.2146
0.2817
−0.3271
−0.2603
−3.743
−3.944
2.308
Time (s)
27.03
27.69
25.68
28.91
33.46
26.61
33.70
48.92
43.22Figure12 shows the PGA of the MPs on MSs 1, 2, and 3. The PGA amplifies with increases in the altitude of the MPs.Figure 12
Calculated PGA of the MPs on MSs 1, 2, and 3.Figure13 shows the calculated PGA of the MPs on MLs 1, 2, and 3. The PGA increases gradually as the seismic wave propagates upward, as shown in Figure 13(a). As shown in Table 5 and Figure 13(b), the PGA increases gradually from the centre to both ends of ML 2. Table 5 and Figure 13(c) illustrate that the change of PGA on ML 3 is small.Figure 13
Calculated PGA of the MPs on MLs 1, 2, and 3.
(a)
(b)
(c)The PGA varied with the elevation, and earthquake ground motions at the top of a valley are strengthened. This illustrates the notable topographic effect. As shown in Table5, the PGA, PGV, and PGD on the slope are larger than those on ML 3, and these values decrease to a normal level at even a small distance away from the slope body (Li [23]). The seismic damage is more severe on the slope.As shown in Figure13, the PGA magnification at the top of the slope is approximately a factor of 2, thus leading to serious collapse at the slope. The real 3D complex numerical model can simulate the phenomenon of landslide and collapse at specific locations (Li [23]). The uneven distributions of ground motions on a marine site caused by topographic effects are useful for the design and construction of underwater facilities.


## 4.1. Numerical Model

A 3D numerical model of a typical submarine site is shown in Figure2. The horizontal dimension of the model is 300×300 m and is divided equally into 100×100 elements in the X and Y directions. The sea level elevation is ±0.000 m, the elevation of the maximum height of the soil at the monitoring point (MP) 15 is −5.000 m, and the elevation of the model bottom is −59.226 m. The vertical dimension of the soil and bedrock (54.226 m) is divided into 16 elements. Therefore, the total number of elements is 160,000. The numbers in circles shown in Figure 2 are the monitoring points for acceleration, velocity, and displacement.Figure 2
3D numerical model and monitoring points for a complex submarine site.As shown in Figure2, three monitoring surfaces (MS) and three monitoring lines (ML) are selected to present the results. MS1, including MPs 0, 1, 2, 3, and 4, is located at the bottom of the model; MS2, including MPs 5, 6, 7, 8, and 9, is located at the top of the bedrock; MS3, including MPs 10, 11, 12, 13, and 14, is located at the top of the model. ML1 is composed of MPs 0, 5, and 10 in a vertical line; ML2 is composed of MPs 15–19, 10, and 20–24; ML3, which is vertical to ML2, is composed of MPs 25–27, 10, and 28–30. MPs 5 and 0 are the projection points of point 10 on the top of the bedrock and the input surface. MPs 1, 6, and 11 are on the same vertical line, and the same situation applies to points (2, 7, and 12), (3, 8, and 13), and (4, 9, and 14).


## 4.2. Material Parameters

Experimental data obtained from static and cyclic triaxial compression tests were used to determine the parameters of the marine soil [1–5, 17]. The material parameters of marine soil involved in our model are listed in Tables 1, 2, and 3. The comparison between the simulated stress-strain-dilatancy relations of soil sample under monotonic/cyclic loading in conventional triaxial tests is shown in Figures 3 and 4, respectively. The bedrock is assumed to be elastic, and its parameters are shown in Table 4.Table 1
Material parameters of the soil.
Parameter
Soil
Compression indexλ
0.026
Swelling indexκ
0.005
Critical state constantRf
2.890
Void ratioN (p′=98 kPa on the N.C.L.)
0.980
Poisson’s ratioν
0.320
Degradation parameter of overconsolidationm
2.200
Degradation parameter of structurem∗
0.100
Evolution parameter of anisotropybr
0.500
Parameter of suctionb
0.500
Void ratioNr (p′=98 kPa on the N.C.L.S.)
1.120Table 2
Initial conditions of the soil.
Parameter
Soil
Initial overconsolidation ratioOCR
2.200
Initial degree of structureR0∗
0.600
Initial anisotropyζ0
0.000
Initial void ratioe0 (p′=98 kPa)
0.900
Initial mean effective stressp′ (kPa)
98.000
Gravity densityγ (kN/m3)
1.940
Volume elasticity coefficient of waterKw (kPa)
2.23
0
e
6
Permeability coefficientk (m/s)
1.000e-7Table 3
Parameters of the moisture characteristics curve of the soil.
Parameter
Soil
Saturated degrees of saturationSrs
0.960
Residual degrees of saturationSrr
0.120
Parameter corresponding to drying AEV (kPa)Sd
320.200
Parameter corresponding to wetting AEV (kPa)Sw
13.600
Initial stiffness of scanning curve (kPa)ks0
77600.000
Parameter of shape functionc1
0.013
Parameter of shape functionc2
0.009
Parameter of shape functionc3
22.000Table 4
Material parameters of the elastic model for bedrock.
Parameter
Bedrock
Linear elastic modulusEc (kPa)
1.375
e
6
Poisson’s ratioν
0.250
Gravity densityγ (kN/m3)
2.200
Permeability coefficientk (m/s)
1.000e-10Figure 3
Stress-strain-dilatancy relations of soil sample in strain-controlled monotonic loading CD tests (confining stress: 98 kPa).Figure 4
Comparison between the calculated and experimental stress-strain relations of soil sample under undrained cyclic loading.


## 4.3. Input Ground Motion

The amplitude, spectral characteristics, and duration of the input ground motion have considerable influences on the seismic response of a marine site. The real three-component ground motions, shown in in Figure5(a), are selected as the input ground motions in the analysis. The original records were filtered and baseline-corrected; the velocity time history before and that after adjustment are shown in Figure 5. As shown in Figure 5, the portion in the blue dashed rectangular box was the selected input time histories. Before inputting the model, the velocity time history should be preliminarily reduced by half and transformed into stress time history (Jing et al. [34]). EW is in the east-west direction (positive: east), NS is in the north-south direction (positive: north), and UD is in the up-down direction (positive: up).Figure 5
Input ground motions before and after adjustment.
(a)
Original records of acceleration
(b)
Acceleration after adjustment
(c)
Velocity
(d)
Displacement


## 4.4. Boundary Condition

For static analysis, the four vertical boundaries were simply supported, and the bottom boundary was fixed. However, in the dynamic analysis, the four vertical boundary conditions were free-field boundaries (shown in Figure6), and the bottom boundary condition was a viscous damping boundary (Lysmer and Kuhlemeyer [35] and Jing et al. [34]).Figure 6
Artificial boundaries for dynamic analysis.


## 4.5. Results of the Static Analysis

The effects of water pressure and the pore water pressure cannot be neglected in the seismic response analysis of marine site. As shown in Figure7(a), the pore water pressure varies with the depth. The initial stress, which was generated by the self-weight stress, is shown in Figure 7(b).Figure 7
Initial pore water pressure and stress states (WLE: water depth above monitoring point 15; IGMSR: input ground motion scaling ratio).
(a)
Pore water pressure (Pa), WLE: +5.000 m, IGMSR: 1
(b)
Z-direction initial stress (Pa), WLE: +5.000 m, IGMSR: 1


## 4.6. Results of the Seismic Response Analysis

### 4.6.1. Spatial and Temporal Distribution of Displacement

Figure8 shows the contours of displacement at the following dynamic times: 5, 10, 20, 30, 40, 50, 60, and 67 s. The deformation increases from the bottom to top of the model. The residual deformations of the MPs on the top are considerably larger than the corresponding residual deformations at the foot of the hill. The numerical simulation results correspond well with the field survey.Figure 8
Contour of displacement at 5, 10, 20, 30, 40, 50, 60, and 67 s.
(a)
Time: 5 s, WLE: +5.000 m, IGMSR: 1
(b)
Time: 10 s, WLE: +5.000 m, IGMSR: 1
(c)
Time: 20 s, WLE: +5.000 m, IGMSR: 1
(d)
Time: 30 s, WLE: +5.000 m, IGMSR: 1
(e)
Time: 40 s, WLE: +5.000 m, IGMSR: 1
(f)
Time: 50 s, WLE: +5.000 m, IGMSR: 1
(g)
Time: 60 s, WLE: +5.000 m, IGMSR: 1
(h)
Time: 67 s, WLE: +5.000 m, IGMSR: 1The spatial and temporal distribution of the displacements shown in Figure8 illustrate that the slope has a large displacement. In the 3D simulation, we can identify the range of sliding, which is critical to our seismic design and reinforcement of marine slope sites.


### 4.6.2. Dynamic Time History

Typical acceleration, velocity, and displacement time histories at MPs 0, 5, 10, and 23 are shown in Figures9–11.Figure 9
Acceleration time histories of MPs 0, 5, 10, and 23.Figure 10
Velocity time histories of MPs 0, 5, 10, and 23.Figure 11
Displacement time histories of MPs 0, 5, 10, and 23.The PGD, PGV, and PGA at MPs 0–30 are presented in Table5.Table 5
PGD, PGV, and PGA at MPs 0–30.
MP
PGD (m)
PGV (m/s)
PGA (m/s2)
EW
NS
UD
EW
NS
UD
EW
NS
UD
Input
Value
−0.0956
−0.1160
0.1074
0.1297
−0.1567
−0.1295
−1.368
1.011
0.7304
Time (s)
26.99
27.65
25.68
28.85
33.42
26.58
28.90
23.87
43.20
MP 0
Value
−0.1879
−0.2293
0.2162
−0.2260
−0.2487
−0.2581
−1.697
1.393
1.389
Time (s)
27.08
27.69
25.68
33.84
33.46
26.61
43.73
43.12
43.22
MP 1
Value
−0.1885
−0.2342
0.2116
−0.2199
−0.2601
−0.2475
−1.617
1.806
1.232
Time (s)
27.16
27.72
25.71
33.97
33.62
26.61
32.38
33.68
43.22
MP 2
Value
−0.1870
−0.2288
0.2133
0.2308
0.2436
−0.2577
−2.176
1.868
1.289
Time (s)
27.08
27.69
25.68
29.22
23.20
26.61
43.70
33.86
43.25
MP 3
Value
−0.1871
−0.2273
0.2129
0.2338
−0.2534
−0.2558
−1.967
1.470
1.313
Time (s)
27.08
27.69
25.66
29.22
33.46
26.58
43.73
43.09
43.22
MP 4
Value
−0.1888
−0.2327
0.2113
−0.2253
−0.2604
−0.2531
−1.718
1.789
1.216
Time (s)
27.14
27.72
25.68
33.94
33.62
26.61
45.65
33.65
43.30
MP 5
Value
−0.1882
−0.2306
0.2163
0.2400
−0.2874
−0.2596
−2.825
−1.699
1.517
Time (s)
27.06
27.69
25.68
29.25
33.46
26.61
33.70
43.46
43.22
MP 6
Value
−0.1907
−0.2375
0.2119
0.2504
−0.3139
−0.2546
−3.339
1.848
1.617
Time (s)
27.08
27.72
25.71
28.93
33.52
26.61
28.99
33.57
43.25
MP 7
Value
−0.1877
−0.2302
0.2133
0.2553
−0.2967
−0.2609
−2.427
−1.765
1.567
Time (s)
27.06
27.69
25.66
29.25
33.49
26.61
28.96
33.38
43.25
MP 8
Value
−0.1887
−0.2286
0.2130
0.2452
−0.2912
−0.2580
−2.577
−1.989
1.329
Time (s)
27.06
27.69
25.68
29.22
33.46
26.61
28.93
32.96
43.22
MP 9
Value
−0.1898
−0.2346
0.2116
0.2621
−0.3108
−0.2570
−2.999
2.026
1.499
Time (s)
27.11
27.72
25.68
29.28
33.52
26.63
28.99
23.97
43.25
MP 10
Value
−0.1913
−0.2310
0.2160
0.2786
−0.3143
−0.2595
−3.447
2.907
1.601
Time (s)
27.03
27.69
25.71
28.91
33.46
26.58
31.79
31.21
31.32
MP 11
Value
−0.1926
−0.2399
0.2095
0.2857
−0.3510
−0.2578
−4.231
3.623
1.693
Time (s)
27.08
27.72
25.71
28.93
33.52
26.63
31.85
35.00
43.28
MP 12
Value
−0.1881
−0.2304
0.2121
0.2840
−0.3397
−0.2608
−3.559
−3.691
−1.499
Time (s)
27.06
27.72
25.68
28.91
33.49
26.61
33.07
33.44
43.54
MP 13
Value
−0.1912
−0.2285
0.2122
0.2801
−0.3256
−0.2598
−4.153
3.132
1.358
Time (s)
27.03
27.69
25.71
28.91
33.46
26.61
33.70
34.98
43.22
MP 14
Value
−0.1885
−0.2336
0.2092
0.2906
−0.3538
−0.2583
−3.609
−3.392
1.720
Time (s)
27.08
27.72
25.74
28.96
33.52
−0.2583
33.09
23.76
43.30
MP 15
Value
−0.1937
0.2319
0.2081
0.3368
−0.3853
−0.2615
4.030
4.834
−2.002
Time (s)
27.08
23.70
25.71
28.96
33.54
26.66
34.29
31.29
29.14
MP 16
Value
−0.1928
0.2316
0.2106
0.3025
−0.3563
−0.2619
3.949
3.121
−2.268
Time (s)
27.08
23.68
25.71
29.28
33.54
26.63
31.77
35.03
29.14
MP 17
Value
0.1902
−0.2303
0.2136
0.3015
−0.3749
−0.2607
3.373
−3.682
1.705
Time (s)
33.12
27.72
25.74
28.93
33.52
26.63
33.97
32.11
29.46
MP 18
Value
0.1901
−0.2303
0.2133
0.2817
−0.3365
−0.2564
3.328
3.125
−1.727
Time (s)
33.12
27.72
25.68
28.93
33.49
26.61
34.26
36.07
33.49
MP 19
Value
−0.1858
−0.2299
0.2157
0.2627
−0.3108
−0.2565
−3.979
2.793
1.541
Time (s)
27.03
27.69
25.66
28.91
33.49
26.61
31.79
34.98
43.25
MP 20
Value
−0.1963
−0.2314
0.2155
0.2803
−0.3530
−0.2645
−4.061
−3.610
1.320
Time (s)
27.06
27.72
25.68
28.93
33.49
26.63
31.82
36.15
43.25
MP 21
Value
−0.1999
−0.2320
0.2128
0.2874
−0.3617
−0.2623
−3.993
3.608
1.637
Time (s)
27.06
27.72
25.68
28.93
33.52
26.61
33.73
33.60
43.25
MP 22
Value
−0.2003
−0.2309
0.2117
0.3059
−0.3758
−0.2625
−4.020
−4.031
1.855
Time (s)
27.08
27.72
25.68
28.96
33.52
26.58
48.12
34.07
43.22
MP 23
Value
−0.2002
0.2312
0.2106
0.3152
−0.3705
−0.2663
−4.345
4.241
−2.111
Time (s)
27.08
23.70
25.71
28.96
33.54
26.61
31.85
33.65
33.54
MP 24
Value
−0.1977
−0.2316
0.2068
0.3198
−0.3819
−0.2638
−4.269
−4.001
−2.247
Time (s)
27.08
27.74
25.68
28.96
33.54
26.63
33.76
50.36
29.14
MP 25
Value
−0.1893
−0.2302
0.2140
0.2748
−0.3117
−0.2597
−3.767
−3.104
1.631
Time (s)
27.03
27.69
25.66
28.91
33.46
26.61
31.79
36.12
43.25
MP 26
Value
−0.1896
−0.2306
0.2141
0.2770
−0.3129
−0.2610
4.401
3.233
1.637
Time (s)
27.03
27.69
25.68
28.91
33.46
26.61
34.21
31.24
43.22
MP 27
Value
−0.1888
0.2315
0.2154
0.2707
−0.3054
−0.2569
−3.752
−2.953
1.713
Time (s)
27.03
27.69
25.68
28.91
33.49
26.61
33.04
33.36
43.20
MP 28
Value
−0.1977
−0.2306
0.2156
0.2723
−0.3573
−0.2635
−3.680
−3.178
1.625
Time (s)
27.06
27.72
25.68
28.93
33.49
26.63
33.73
36.17
43.25
MP 29
Value
−0.1929
−0.2295
0.2138
0.2782
−0.3241
−0.2589
−3.652
2.763
−1.735
Time (s)
27.03
27.69
25.68
28.91
33.49
26.61
28.96
34.98
33.46
MP 30
Value
−0.1918
−0.2290
0.2146
0.2817
−0.3271
−0.2603
−3.743
−3.944
2.308
Time (s)
27.03
27.69
25.68
28.91
33.46
26.61
33.70
48.92
43.22Figure12 shows the PGA of the MPs on MSs 1, 2, and 3. The PGA amplifies with increases in the altitude of the MPs.Figure 12
Calculated PGA of the MPs on MSs 1, 2, and 3.Figure13 shows the calculated PGA of the MPs on MLs 1, 2, and 3. The PGA increases gradually as the seismic wave propagates upward, as shown in Figure 13(a). As shown in Table 5 and Figure 13(b), the PGA increases gradually from the centre to both ends of ML 2. Table 5 and Figure 13(c) illustrate that the change of PGA on ML 3 is small.Figure 13
Calculated PGA of the MPs on MLs 1, 2, and 3.
(a)
(b)
(c)The PGA varied with the elevation, and earthquake ground motions at the top of a valley are strengthened. This illustrates the notable topographic effect. As shown in Table5, the PGA, PGV, and PGD on the slope are larger than those on ML 3, and these values decrease to a normal level at even a small distance away from the slope body (Li [23]). The seismic damage is more severe on the slope.As shown in Figure13, the PGA magnification at the top of the slope is approximately a factor of 2, thus leading to serious collapse at the slope. The real 3D complex numerical model can simulate the phenomenon of landslide and collapse at specific locations (Li [23]). The uneven distributions of ground motions on a marine site caused by topographic effects are useful for the design and construction of underwater facilities.


## 4.6.1. Spatial and Temporal Distribution of Displacement

Figure8 shows the contours of displacement at the following dynamic times: 5, 10, 20, 30, 40, 50, 60, and 67 s. The deformation increases from the bottom to top of the model. The residual deformations of the MPs on the top are considerably larger than the corresponding residual deformations at the foot of the hill. The numerical simulation results correspond well with the field survey.Figure 8
Contour of displacement at 5, 10, 20, 30, 40, 50, 60, and 67 s.
(a)
Time: 5 s, WLE: +5.000 m, IGMSR: 1
(b)
Time: 10 s, WLE: +5.000 m, IGMSR: 1
(c)
Time: 20 s, WLE: +5.000 m, IGMSR: 1
(d)
Time: 30 s, WLE: +5.000 m, IGMSR: 1
(e)
Time: 40 s, WLE: +5.000 m, IGMSR: 1
(f)
Time: 50 s, WLE: +5.000 m, IGMSR: 1
(g)
Time: 60 s, WLE: +5.000 m, IGMSR: 1
(h)
Time: 67 s, WLE: +5.000 m, IGMSR: 1The spatial and temporal distribution of the displacements shown in Figure8 illustrate that the slope has a large displacement. In the 3D simulation, we can identify the range of sliding, which is critical to our seismic design and reinforcement of marine slope sites.


## 4.6.2. Dynamic Time History

Typical acceleration, velocity, and displacement time histories at MPs 0, 5, 10, and 23 are shown in Figures9–11.Figure 9
Acceleration time histories of MPs 0, 5, 10, and 23.Figure 10
Velocity time histories of MPs 0, 5, 10, and 23.Figure 11
Displacement time histories of MPs 0, 5, 10, and 23.The PGD, PGV, and PGA at MPs 0–30 are presented in Table5.Table 5
PGD, PGV, and PGA at MPs 0–30.
MP
PGD (m)
PGV (m/s)
PGA (m/s2)
EW
NS
UD
EW
NS
UD
EW
NS
UD
Input
Value
−0.0956
−0.1160
0.1074
0.1297
−0.1567
−0.1295
−1.368
1.011
0.7304
Time (s)
26.99
27.65
25.68
28.85
33.42
26.58
28.90
23.87
43.20
MP 0
Value
−0.1879
−0.2293
0.2162
−0.2260
−0.2487
−0.2581
−1.697
1.393
1.389
Time (s)
27.08
27.69
25.68
33.84
33.46
26.61
43.73
43.12
43.22
MP 1
Value
−0.1885
−0.2342
0.2116
−0.2199
−0.2601
−0.2475
−1.617
1.806
1.232
Time (s)
27.16
27.72
25.71
33.97
33.62
26.61
32.38
33.68
43.22
MP 2
Value
−0.1870
−0.2288
0.2133
0.2308
0.2436
−0.2577
−2.176
1.868
1.289
Time (s)
27.08
27.69
25.68
29.22
23.20
26.61
43.70
33.86
43.25
MP 3
Value
−0.1871
−0.2273
0.2129
0.2338
−0.2534
−0.2558
−1.967
1.470
1.313
Time (s)
27.08
27.69
25.66
29.22
33.46
26.58
43.73
43.09
43.22
MP 4
Value
−0.1888
−0.2327
0.2113
−0.2253
−0.2604
−0.2531
−1.718
1.789
1.216
Time (s)
27.14
27.72
25.68
33.94
33.62
26.61
45.65
33.65
43.30
MP 5
Value
−0.1882
−0.2306
0.2163
0.2400
−0.2874
−0.2596
−2.825
−1.699
1.517
Time (s)
27.06
27.69
25.68
29.25
33.46
26.61
33.70
43.46
43.22
MP 6
Value
−0.1907
−0.2375
0.2119
0.2504
−0.3139
−0.2546
−3.339
1.848
1.617
Time (s)
27.08
27.72
25.71
28.93
33.52
26.61
28.99
33.57
43.25
MP 7
Value
−0.1877
−0.2302
0.2133
0.2553
−0.2967
−0.2609
−2.427
−1.765
1.567
Time (s)
27.06
27.69
25.66
29.25
33.49
26.61
28.96
33.38
43.25
MP 8
Value
−0.1887
−0.2286
0.2130
0.2452
−0.2912
−0.2580
−2.577
−1.989
1.329
Time (s)
27.06
27.69
25.68
29.22
33.46
26.61
28.93
32.96
43.22
MP 9
Value
−0.1898
−0.2346
0.2116
0.2621
−0.3108
−0.2570
−2.999
2.026
1.499
Time (s)
27.11
27.72
25.68
29.28
33.52
26.63
28.99
23.97
43.25
MP 10
Value
−0.1913
−0.2310
0.2160
0.2786
−0.3143
−0.2595
−3.447
2.907
1.601
Time (s)
27.03
27.69
25.71
28.91
33.46
26.58
31.79
31.21
31.32
MP 11
Value
−0.1926
−0.2399
0.2095
0.2857
−0.3510
−0.2578
−4.231
3.623
1.693
Time (s)
27.08
27.72
25.71
28.93
33.52
26.63
31.85
35.00
43.28
MP 12
Value
−0.1881
−0.2304
0.2121
0.2840
−0.3397
−0.2608
−3.559
−3.691
−1.499
Time (s)
27.06
27.72
25.68
28.91
33.49
26.61
33.07
33.44
43.54
MP 13
Value
−0.1912
−0.2285
0.2122
0.2801
−0.3256
−0.2598
−4.153
3.132
1.358
Time (s)
27.03
27.69
25.71
28.91
33.46
26.61
33.70
34.98
43.22
MP 14
Value
−0.1885
−0.2336
0.2092
0.2906
−0.3538
−0.2583
−3.609
−3.392
1.720
Time (s)
27.08
27.72
25.74
28.96
33.52
−0.2583
33.09
23.76
43.30
MP 15
Value
−0.1937
0.2319
0.2081
0.3368
−0.3853
−0.2615
4.030
4.834
−2.002
Time (s)
27.08
23.70
25.71
28.96
33.54
26.66
34.29
31.29
29.14
MP 16
Value
−0.1928
0.2316
0.2106
0.3025
−0.3563
−0.2619
3.949
3.121
−2.268
Time (s)
27.08
23.68
25.71
29.28
33.54
26.63
31.77
35.03
29.14
MP 17
Value
0.1902
−0.2303
0.2136
0.3015
−0.3749
−0.2607
3.373
−3.682
1.705
Time (s)
33.12
27.72
25.74
28.93
33.52
26.63
33.97
32.11
29.46
MP 18
Value
0.1901
−0.2303
0.2133
0.2817
−0.3365
−0.2564
3.328
3.125
−1.727
Time (s)
33.12
27.72
25.68
28.93
33.49
26.61
34.26
36.07
33.49
MP 19
Value
−0.1858
−0.2299
0.2157
0.2627
−0.3108
−0.2565
−3.979
2.793
1.541
Time (s)
27.03
27.69
25.66
28.91
33.49
26.61
31.79
34.98
43.25
MP 20
Value
−0.1963
−0.2314
0.2155
0.2803
−0.3530
−0.2645
−4.061
−3.610
1.320
Time (s)
27.06
27.72
25.68
28.93
33.49
26.63
31.82
36.15
43.25
MP 21
Value
−0.1999
−0.2320
0.2128
0.2874
−0.3617
−0.2623
−3.993
3.608
1.637
Time (s)
27.06
27.72
25.68
28.93
33.52
26.61
33.73
33.60
43.25
MP 22
Value
−0.2003
−0.2309
0.2117
0.3059
−0.3758
−0.2625
−4.020
−4.031
1.855
Time (s)
27.08
27.72
25.68
28.96
33.52
26.58
48.12
34.07
43.22
MP 23
Value
−0.2002
0.2312
0.2106
0.3152
−0.3705
−0.2663
−4.345
4.241
−2.111
Time (s)
27.08
23.70
25.71
28.96
33.54
26.61
31.85
33.65
33.54
MP 24
Value
−0.1977
−0.2316
0.2068
0.3198
−0.3819
−0.2638
−4.269
−4.001
−2.247
Time (s)
27.08
27.74
25.68
28.96
33.54
26.63
33.76
50.36
29.14
MP 25
Value
−0.1893
−0.2302
0.2140
0.2748
−0.3117
−0.2597
−3.767
−3.104
1.631
Time (s)
27.03
27.69
25.66
28.91
33.46
26.61
31.79
36.12
43.25
MP 26
Value
−0.1896
−0.2306
0.2141
0.2770
−0.3129
−0.2610
4.401
3.233
1.637
Time (s)
27.03
27.69
25.68
28.91
33.46
26.61
34.21
31.24
43.22
MP 27
Value
−0.1888
0.2315
0.2154
0.2707
−0.3054
−0.2569
−3.752
−2.953
1.713
Time (s)
27.03
27.69
25.68
28.91
33.49
26.61
33.04
33.36
43.20
MP 28
Value
−0.1977
−0.2306
0.2156
0.2723
−0.3573
−0.2635
−3.680
−3.178
1.625
Time (s)
27.06
27.72
25.68
28.93
33.49
26.63
33.73
36.17
43.25
MP 29
Value
−0.1929
−0.2295
0.2138
0.2782
−0.3241
−0.2589
−3.652
2.763
−1.735
Time (s)
27.03
27.69
25.68
28.91
33.49
26.61
28.96
34.98
33.46
MP 30
Value
−0.1918
−0.2290
0.2146
0.2817
−0.3271
−0.2603
−3.743
−3.944
2.308
Time (s)
27.03
27.69
25.68
28.91
33.46
26.61
33.70
48.92
43.22Figure12 shows the PGA of the MPs on MSs 1, 2, and 3. The PGA amplifies with increases in the altitude of the MPs.Figure 12
Calculated PGA of the MPs on MSs 1, 2, and 3.Figure13 shows the calculated PGA of the MPs on MLs 1, 2, and 3. The PGA increases gradually as the seismic wave propagates upward, as shown in Figure 13(a). As shown in Table 5 and Figure 13(b), the PGA increases gradually from the centre to both ends of ML 2. Table 5 and Figure 13(c) illustrate that the change of PGA on ML 3 is small.Figure 13
Calculated PGA of the MPs on MLs 1, 2, and 3.
(a)
(b)
(c)The PGA varied with the elevation, and earthquake ground motions at the top of a valley are strengthened. This illustrates the notable topographic effect. As shown in Table5, the PGA, PGV, and PGD on the slope are larger than those on ML 3, and these values decrease to a normal level at even a small distance away from the slope body (Li [23]). The seismic damage is more severe on the slope.As shown in Figure13, the PGA magnification at the top of the slope is approximately a factor of 2, thus leading to serious collapse at the slope. The real 3D complex numerical model can simulate the phenomenon of landslide and collapse at specific locations (Li [23]). The uneven distributions of ground motions on a marine site caused by topographic effects are useful for the design and construction of underwater facilities.


## 5. Conclusion

A new elastoplastic constitutive model with rotation hardening, which can systematically describe the monotonic and cyclic mechanical behaviours of typical marine soils combining the subloading, normal, and superloading yield surfaces, is employed in the seismic response analysis of a 3D marine site. New evolution equations for the stress-induced anisotropy and degree of overconsolidation are proposed. This model can describe the unified behaviour of unsaturated soil and saturated soil by using degree of saturation as independent state variable. By using the transform stress method, our model uniquely describes the overall mechanical properties of marine soils under general stress states, without changing the parameter values. Based on this model and three-phase field theory, an effective stress-based, fully coupled, explicit finite element–finite difference method is developed. The FEM and explicit integration method are applied in spatial and temporal discretization. An improved artificial boundary is adopted in the analysis. A finite deformation analysis is performed by introducing the Green-Naghdi rate tensor. This model is suitable for the static-dynamic analysis of complex 3D sites.A 3D numerical model for a typical marine sites is developed. The results indicate that the proposed method can accurately simulate a seismic disaster at a marine site. The operating status before the earthquake, the damage state after the earthquake, and the damage generation process are investigated. The displacement, velocity, and acceleration time histories suggest that the ground motion intensity will increase due to local uneven complex topography and site effects. The results also demonstrate the temporal and spatial distributions of landslides and collapse at the specific marine site. Local site and topographic effects control the damage level and pattern on the marine site. The ground motion is most easily focused at the top of the slope, where the wave finally reached. At that location, the ground motion is magnified approximately two times compared to normal values and will lead to serious collapse at the slope.



---

*Source: 3892193-2017-12-19.xml*