Abstract
This study analyzes the impact of the number of ground motions on the seismic fragility of a high earth-rockfill dam and the estimation of reasonable fragility parameters based on a sufficient number of earthquake records. In this paper, the vertical deformation is obtained using the three-dimensional finite element program DYNE3WAC combined with the Pastor–Zienkiewicz–Chan model and Biot dynamic consolidation theory. The relative seismic settlement rate is considered the damage index for the seismic fragility analysis of the dam. The fragility curves of the high earth-rockfill dam are determined by the multiple stripe analysis (MSA) method. A set of seismic waves is chosen based on the spectrum in the Chinese hydraulic structure seismic code. With an increasing number of earthquake records, the coefficients of variation (COV) of the mean and standard deviation (STD) of the relative seismic settlement rate decrease and tend to stabilize when the number of earthquake records reaches 34. The estimated fragility parameters and are constant when the number of earthquake records exceeds 34. The requisite number of earthquake records for an accurate fragility estimation is determined by analyzing the lower and upper confidence intervals for the estimated and . The 95% and 90% confidence interval can accurately estimate the fragility of a high earth-rockfill dam when the number of ground motion records reaches 36 and 32, respectively. The results of the fragility analysis demonstrate that the DYNE3WAC program and MSA method can provide an effective basis for determining fragility curves. Furthermore, the sensitivity analysis of earthquake records is essential for assessing the seismic fragility of high earth-rockfill dams.
1. Introduction
High earth-rockfill dams are an important part of the China’s infrastructure and represent sources of massive amounts of sustainable, renewable and clean hydropower energy [1]. China is one of the countries most affected by natural disasters, especially earthquakes, due to its specific location and large territory. The 1999 Chi-Chi earthquake, the 2008 Wenchuan earthquake, and the 2010 Yushu earthquake have demonstrated that the dams are vulnerable to earthquakes. The great Wenchuan earthquake that occurred in China caused earthquake damage to more than 2500 earth-rockfill dams. High earth-rockfill dams with heights of 200-300 meters are located in areas prone to strong earthquakes [2, 3]. Therefore, seismic performance assessments must be performed to assess the safety of these dams.
Seismic fragility analysis is one of the most effective methods for evaluating seismic performance and can effectively estimate the risk to dams affected by earthquakes; thus, this technique has gradually becomes an important method in seismic safety evaluation [4]. Fragility analysis can predict the probability of damage when structures reach or exceed a certain limit state probability under different strengths of seismic action and can employ a fragility curve to describe the probability distributions of the limit states [5, 6]. The fragility curves can be used to describe the probability of selected structure responses exceeding a certain critical level with a ground motion intensity measure (IM). To construct a fragility curve, incremental dynamic analysis (IDA) [7] and multiple strip analysis (MSA) [8] have been extensively applied as seismic analysis tools. The IDA method is employed to estimate the limit state capacity and seismic performance by performing a series of nonlinear time history analyses for a suite of ground motion records. In contrast to IDA, the MSA approach is a nonlinear dynamic analysis method and involves a sufficient set of stripes that are used to analyze and evaluate structural seismic responses under different IM levels. Baker [9] found that MSA produces more efficient fragility estimates than IDA for a given number of structural analyses and proposed an assessment approach. During the past several decades, many studies have been conducted to research the seismic fragility of arch dams and concrete dams [10, 11]. However, few papers have analyzed the seismic fragility of high earth-rockfill dams. Pang and Xu et al. [12] recently evaluated the seismic fragility of concrete-faced rockfill dams (CFRDs) subjected to the acceleration response spectra of 11 earthquake waves. However, these fragility analyses did not consider a large number of earthquake records, which might affect the parameter estimation for the fragility function. Several studies have been performed to promote the accurate estimation of the probability of exceeding a certain limit state in terms of earthquake IMs. Cimellaro et al. [13] investigated the impact of the number of ground motions on elastic shear-type buildings using a multidegree of freedom (MDOF) system and concluded that the minimum number of earthquake records needed to accurately estimate the parameters of the fragility functions is 23. Baker and Cornell [14] investigated a method of selecting ground motions and concluded that selecting earthquake records according to spectral shape can reduce the variance and bias of structural responses.
Deformation and stability are the two major aspects of the earthquake destruction patterns of earth-rockfill dams and are typically regarded as seismic performance assessment indexes [15–17]. Generally, the deformation index, especially the relative seismic settlement ratio [18, 19], is the most common seismic assessment index for earth-rockfill dams. The relative seismic settlement ratio is an important physical parameter when evaluating the seismic performance of earth-rockfill dams [20–22]. Furthermore, the index can illustrate the safety of the earth-rockfill dam and can be easily found in seismic records. Therefore, it is reasonable and feasible to regard the relative seismic settlement rate as an evaluation index of the seismic performance of earth-rockfill dams.
This paper aims to investigate the effect of the number of earthquake records on the fragility of high earth-rockfill dams and to perform a sensitivity analysis of the fragility parameters. The seismic fragility assessment of the dam is studied based on the MSA method. First, 60 seismic waves are selected according to the site conditions and the spectrum in China’s hydraulic structure seismic code. The relative seismic settlement rate of the high earth-rockfill dam is selected as the damage index. The three-dimensional elastoplastic finite element program DYNE3WAC (DYNamic Earthquake Analysis Program 3D Window Version for Academic), which is the 3D version of SWANDYNE II [23], can effectively calculate geotechnical engineering under static and dynamic conditions. The DYNE3WAC program and MSA method can reasonably determine the fragility curves of high earth-rockfill dams. Considering the number of ground motion inputs, the appropriate number of seismic records can be determined according to the fragility equation parameters. The coefficient of variation (COV) of the mean and standard deviation (STD) of the relative seismic settlement rate and the estimated and are analyzed. The estimated fragility parameters and tend to be constant when the number of earthquake records increases. The confidence interval can accurately estimate the fragility of a high earth-rockfill dam when the number of earthquake records is sufficiently high.
2. Seismic Fragility Method Based on MSA
Performance-based earthquake engineering (PBEE) methods, such as IDA and MSA, have been proposed by many scholars. The difference between IDA and MSA lies in the record selection. In the MSA method, ground motion records are selected to fit all predefined IM levels. Furthermore, the MSA method is a sufficient set of single stripes used to analyze and evaluate structural seismic responses at different IM levels [24]. Compared with IDA, MSA is less time-consuming and relatively accurate. The seismic fragility curves provide the conditional probabilities of the structural response reaching or exceeding certain limit states that correspond to the seismic capacity under different earthquake intensities. The probability that the structure response exceeds specified limit states can be described as a function of the IM of ground motion. For fragility curves, peak ground acceleration (PGA) is commonly exploited as a convenient IM. The variation in the IM values of ground motions is commonly assumed to be related to the damage state of the selected structure following a lognormal distribution. The fragility function is generally defined as a lognormal cumulative distribution function and is shown in where is the probability that a ground motion with results in a structure that attained a damaged state, is the standard normal cumulative distribution function, and and are the median and STD of , respectively.
The MSA approach does not suggest that estimates of IM amplitudes for all the selected ground motions result in the exceedance of limit states. To estimate the fragility parameters, the likelihood function can be defined as follows:where is the number of levels, is a product over all levels, and and are the number of times the limit state is exceeded and the number of ground motions with , respectively.
The parameters and can be estimated by maximizing this likelihood function. It is equivalent and numerically easier to maximize the logarithm of the likelihood, as shown in
2.1. Seismic Records
An important step in seismic fragility analysis is the selection of a representative set of earthquake motions at different levels of ground motion intensity. For ground motion selection, the magnitude, distance, and soil condition should be considered. There are uncertainties in the magnitude and distance, so it is reasonable to select earthquake records within appropriate ranges. The general consensus is to select ground motions whose magnitude and distance are within a close range of the target magnitude and distance. The dynamic peak acceleration with a transcendental probability of 2% in 100 years is selected as the seismic response spectrum. Shome [25] notes that the use of 10~20 seismic waves can produce an accurate estimation of the structure's earthquake response. However, due to the complex structure of high earth-rockfill dams, the geological conditions are different from those of ordinary buildings.
Because randomness is an inherent characteristic of earthquake ground motions [26], seismic randomness analysis should be researched and applied. The seismic performance assessment of the earth-rockfill dam employs fragility analysis based on 10–20 ground motion records from the PEER database, which has tens of thousands of ground motions [27, 28]. To analyze the relationship between the seismic fragility of high earth-rockfill dams and the number of earthquake records, the number of ground motions is selected to be 3 times the upper limit of Shome’s recommended value. Therefore, the selection of the ground motion input should take full account of the randomness of earthquakes. The response spectra based on the site conditions of a high earth-rockfill dam are selected as the target response spectra. Furthermore, 60 actual seismic records that agree well with the target response spectra based on site conditions are selected in PEER. These records are obtained from 60 different earthquakes with magnitudes ranging from 5.0 to 7.0 and with epicentral distances ranging from 10 to 55 km. The ground records are scaled and matched to the 5% damped target spectrum using SeismoMatch software. The distributions of PGAs, magnitudes and epicentral distances of these seismic waves are shown in Figure 1(a), and the seismic response spectrum of acceleration is shown in Figure 1(b):
(a) Epicentral distance, magnitude, and the PGA
(b) The curves of earthquake acceleration response spectra
2.2. Damage Index
To obtain the fragility curves for high earth-rockfill dams subjected to seismic loads, it is important to define the damage state depending on the demand parameters of the earth-rockfill dam. Deformation is the primary aspect of the seismic performance assessment that should be considered when evaluating the seismic safety of earth-rockfill dams [29, 30]. Permanent seismic deformation has been selected as the seismic design standard in the seismic safety evaluation of earth-rockfill dams [16, 24]. The relative seismic settlement rate can comprehensively consider the earthquake-induced deformations of dams. Kong and Pang et al. [31] considered relative settlement ratios (crest vertical settlement values/height of the dam) of 0.4%, 0.7%, and 1% of the dam crest as the assessment limits and analyzed the fragility when the studied dam exhibited minor, moderate and severe failure. The three limit states corresponding to relative settlement ratios of 0.4%, 0.7%, and 1.0% of the dam crest were chosen as the minor, moderate, and severe, respectively. To evaluate the seismic fragility in terms of damage levels based on the relative settlement ratio of an earth-rockfill dam, damage grades should be proposed. In this study, three limit states with relative settlement ratios of 0.4%, 0.7%, and 1.0% were selected, and they correspond to three damage grades of minor (limit state 1, LS1), moderate (limit state 2, LS2) and severe (limit state 3, LS3), respectively.
3. Seismic Fragility of High Earth-Rockfill Dams
3.1. Pastor–Zienkiewicz–Chan (PZC) Model
Using the well-known generalized plasticity model, the dynamic response and the relative settlement rate can be used to effectively and reasonably analyze the material of a high earth-rockfill dam. The generalized plasticity framework was first introduced by Zienkiewicz and Mroz [32]. The main advantage of the generalized plasticity theory is that the yield surface and plastic potential surface are not defined explicitly, and the consistency law is not required to determine the plastic modulus. The elastoplastic tensor in the generalized plasticity theory is as follows:where represents the elastic stiffness tensor, is the loading direction vector, is the plastic flow direction vector, and is the plastic modulus under loading (L) or unloading (U) conditions.
A more detailed introduction of the PZC model appears in the Appendix.
3.2. Finite Element Model (FEM) and Method
In recent years, many high earth-rockfill dams have been built to satisfy the growing demand for water resources and hydropower in Southwest China. However, these dams are distributed in areas that experience strong ground motions. Hence, seismic fragility is a matter of great concern in dam safety. In this paper, the Nuozhadu earth-rockfill dam located on the Lancang River in Yunnan Province, Southwest China, was analyzed. The dam crest has a longitudinal length of 630 m. The reservoir has a storage capacity of 237.03×108 m3 with a normal storage water level of 812 m. The upstream and downstream sides of the dam slope have grades of 1:1.9 and 1:1.8, respectively. The 2D finite typical section of the earth-rockfill dam, with a maximum height of 261.5 m and a top width of 18 m, is shown in Figure 2.
The earth-rockfill dam was constructed in 2008 and completed at the end of 2012. The Nuozhadu dam is the highest earth-rockfill dam in China. The nonlinear finite element program DYNE3WAC is applied to perform static, seepage and dynamic analyses. The water loading is applied to the core wall of the earth-rockfill dam. The permeability coefficient of the clay is far lower than the permeability coefficient of the rockfill I and rockfill II units. Therefore, the core wall of the earth-rockfill dam has a good antiseepage effect. The material parameters of the Nuozhadu dam are shown in Table 1.
The initial stress field and phreatic surface can be used for the dynamic analysis of the dam. The two-dimensional finite element model (FEM) and the material parameters of the earth-rockfill dam are added in the manuscript according to the phreatic surface in Figure 3.
According to the phreatic line in Figure 3, the clay and rockfill of the dam above the phreatic surface which were located in the unsaturated zone were discretized by displacement elements. The clay and rockfill of the dam below the phreatic surface were modeled with displacement-pore water pressure coupling elements. The base and lateral boundaries were assumed to be impervious. The nodes at the normal water level were traction free with zero prescribed pore pressure. Figure 4 shows the finite element mesh of the Nuozhadu earth-rockfill dam with 21034 nodes and 20389 elements.
The PZC model parameters in the paper were referenced from Dong’s research which analyzed the permanent displacement of the Nuozhadu earth-rockfill dam with the same model [33]. The materials of the Nuozhadu dam can be analyzed by the PZC model, and the material parameters are shown in Table 2.
Following the Biot formulation governing the behavior of porous media and neglecting relative fluid acceleration, the equation of motion can be readily written with the solid displacement u and fluid pore pressure p as
A more detailed introduction of the Biot formulation appears in the Appendix.
3.3. Dynamic Analysis Results
The relative seismic settlement rate of the high earth-rockfill dam can be obtained according to the selected seismic wave and dynamic analysis. The earthquake deformation exceeds the relative settlement rate of the dam top by 1% when the PGA is 1.0 g, so the amplitude PGA range is 0.1 g to 1.0 g. Owing to a large number of dynamic finite element calculations and seismic fragility analysis, the ASUS workstation with 128 G of RAM and a main frequency of 3.3 GHz was used to calculate the dynamic response of the dam. The time history curves of the relative seismic settlement rate for a PGA of 0.3 g are shown in Figure 5.
PGA data ranging from 0.1 g to 1.0 g with increments of 0.05 g are input to perform nonlinear dynamic analyses for 60 earthquakes. According to the MSA, statistical information about the relative settlement rate can be obtained. The numerical analysis is performed to simulate the deformation of the high earth-rockfill dam. The seismic fragility curves shown in Figure 6 clearly illustrate the fragility.
As observed from Figure 6, the probability of exceeding the three damage levels varies according to the PGA. When the PGA is 0.3 g, the probability of exceeding the moderate level is only 3.3%, while the probability of exceeding the minor level is 47.2%.
When the PGA is 0.4 g, the probability of exceeding the minor level is 79.4%, while the probabilities of exceeding the moderate and severe damage levels are 18.6% and 1.6%, respectively. Additionally, when the PGA is 0.8 g, the probability that the earthquake load causes the relative seismic settlement rate of the high earth-rockfill dam to exceed the minor level is 100%, while the probabilities of exceeding the moderate and severe damage levels are 91.2% and 56.3%, respectively. It can be concluded that the quantitative results for the performance of the earth-rockfill dam are helpful and useful to dam owners for managing uncertainties.
4. Fragility Parameter Analysis
From the set of 60 earthquake records, a subset can be chosen randomly to obtain the mean and STD of the maximum relative seismic settlement rate. The values of the mean and STD are selected randomly from 1000 groups without repetition, and 1000 mean values are obtained. Then, the COV of the mean and STD can be calculated. To obtain the relationship between the relative seismic settlement rate and the number of ground motions, 60 selected ground motions were applied and analyzed. Figure 7(a) shows the COV of the means of 1000 randomly chosen groups when the number of ground motions is 4. The relationship between the COV of the mean of the maximum relative seismic settlement rate and different numbers of earthquake records are shown in Figure 7(b).
(a)
(b)
Figure 7(b) shows an obvious reduction in the COV of the mean in response to an increase in the number of earthquakes. The COV in the initial stages from 4 to 10 earthquakes decreases rapidly and then declines more gradually until it is nearly stable at 34 earthquakes. Figure 8(a) shows the COV of the STD for 1000 randomly chosen groups when the number of ground motions is 4. The relationship between the COV of the STD of the relative seismic settlement rate and different numbers of earthquake records is shown in Figure 8(b) by repeating the same procedure for up to 60 records.
(a)
(b)
In Figure 8(b), the trend of the COV of the STD with the number of records is the same as that of the COV of the mean, initially decreasing then becoming stable as the number of earthquakes increases. Furthermore, each dot in Figures 7(a) and 8(a) corresponds to one of the 1000 random simulations for the same number of ground motions.
4.1. Effect of the Number of Earthquake Records on the Fragility Parameters
The effect of the number of earthquake records on the fragility function parameter ( and ) estimation is investigated. Each number of ground motion records contains 1000 random groups without repetition with the same members in each group. Then, the fragility function parameters for each group can be estimated based on (1) and (3). Therefore, there are 1000 estimated values and 1000 estimated values for each number of earthquake records, and the related mean can be obtained. The MSA parameters corresponding to different ground motion records were obtained after a large number of nonlinear finite element calculations. Figure 9(a) shows the seismic fragility parameters of 1000 random groups when the number of ground motions is 4.
(a)
(b)
Figure 9(b) shows that the estimated value decreases as the number of earthquakes increases and tends to be constant when the number of earthquakes exceeds 34. Figure 10(a) shows the seismic fragility parameters of 1000 random groups when the number of the ground motions is 4.
(a)
(b)
Figure 10(b) demonstrates that the estimated value increases as the number of earthquake records increases and tends to be constant when the number of earthquake records exceeds 34. Comparison of the estimated parameters with different numbers of earthquake records reveals that when the number of earthquake records exceeds 34, stable parameters can be estimated. Notably, the stable parameters in the different numbers of records may have different confidence intervals. Therefore, to obtain a sufficient number of records for accurate parameter estimation, more in-depth research is needed. Each dot in Figures 9 and 10 corresponds to 1000 random simulations based on the number of ground motions. From the results of Figures 9 and 10, the relationship between the fragility function parameters ( and ) and the numbers of ground motions can be analyzed.
4.2. The 95% Confidence Interval for Fragility Parameters
The effect of the number of earthquake records on the 95% confidence interval for estimated parameters is investigated to determine the number of records needed to obtain accurate fragility estimation. For each number of earthquake records, the same 1000 possible groups are selected, and the confidence interval of the estimation value with the same 95% probability for each group is developed. According to the concept of the fragility function, the lower and upper bounds of the 95% confidence intervals for the estimated and in each group can be obtained using Eqs. (7) and (8) based on the t and distributions: where and are the estimated values for the th group when the number of earthquake records is .
Therefore, there are 1000 values for the lower and upper 95% confidence intervals for each parameter for a constant number of earthquake records, and the related mean values can be obtained. The relationship between the mean values of the lower and upper 95% confidence intervals and different numbers of earthquake records are shown in Figure 11.
(a)
(b)
Notably, the values in Figure 11(a) are the exponents of for the related values of the lower and upper confidence intervals. Figure 11(a) clearly illustrates that the range of the 95% confidence interval decreases and remains constant as the number of earthquake records increases to 36. In Figure 11(b), the values of exhibit the same trend with the number of records as those of , i.e., decreasing and tending to be stable as the number of earthquakes increases. However, the sufficient number of earthquake records cannot be determined from only the example in Figure 11. When the scope of the 95% confidence interval for both estimated and values are less than 0.2, the related number of earthquake records can accurately estimate the fragility of dams subjected to seismic loads. Therefore, it is concluded that the required number of earthquake records for an accurate fragility estimation of dams subjected to seismic loads is 36.
Compared to the 95% confidence intervals for the estimated and values, the 90% confidence levels for the estimated parameters can illustrate the relationship between the seismic fragility parameters and the number of ground motions. According to the concept of the fragility function, the lower and upper bounds of the 90% confidence intervals for the estimated and values in each group can be obtained using Eqs. (9) and (10) based on the t and distributions:where and are the estimated values for the th group when the number of earthquake records is .
The relationship between the mean values of the lower and upper bounds of the 90% confidence intervals and different numbers of earthquake records are shown in Figure 12.
(a)
(b)
From Figure 12, the 90% confidence levels for the estimated parameters have the almost same trend as those in Figure 11, and the seismic fragility parameters are stable when the number of earthquake records reaches 32.
5. Conclusions
This study uses the three-dimensional dynamic finite program DYNE3WAC combined with the PZC elastoplastic constitutive model and the Biot consolidation theory to determine the vertical deformations in a high earth-rockfill dam. The relative crest settlement rate is used to estimate the damage state for the high earth-rockfill dam. The seismic fragility analysis based on the MSA method can reasonably assess the seismic performance of the high earth-rockfill dam. Considering the randomness of ground motion records, the parameters of seismic fragility can comprehensively reflect the sensitivity of the analysis to the number of earthquake waves.(1)In this paper, the elastoplastic finite element program is used to analyze the dynamic responses of the high earth-rockfill dam and reasonably calculate the vertical deformations of the crest. The relative crest settlement ratio of the high earth-rockfill dam is established as the damage index, reflecting the seismic fragility. The dam failure levels based on the relative crest settlement ratio are established. The seismic fragility curves of the minor, moderate, and severe failure levels are described by the MSA approach.(2)The MSA approach can reasonably describe the fragility curves and estimate fragility function parameters, which are influenced by the number of ground motion records. The influence of the number of earthquake records on the seismic fragility parameters of the high earth-rockfill dam is assessed based on the seismic fragility parameters. The effect of the number of earthquake records based on the fragility parameters and the 95% and 90% confidence intervals for the fragility parameters are accurately determined.(3)The COV of the mean and STD of the relative seismic settlement clearly indicate that the value of these indexes decrease and tend to be stable as the number of earthquake records increases. The estimated fragility parameters and tend to be constant when the number of earthquake records exceeds a certain number. The number of earthquake records required for accurate fragility estimation is determined by analyzing the 95% and 90% confidence intervals for the fragility parameters, which tend to be stable when the number of earthquakes reaches 36 and 32, respectively.
Appendix
A. PZC Elastoplastic Model
The PZC model can simulate the mechanical behaviors of sand but cannot simulate the behavior of rockfill material, which is more complex and contains particles larger than sand. The model can be formulated in terms of the three invariants of the effective stress tensor p, q, and :where p, q, and represent the mean effective stress, deviatoric stress and Lode’s angle, respectively; is the first stress invariant of the stress tensor; and and represent the second and third deviatoric stress invariants of the deviatoric stress tensor, respectively.
The original PZC model assumes plane-strain conditions and cannot simulate the three-dimensional response of earth-rockfill dams. The Hoek elasticity defines the shear and bulk moduli (G and K), which are dependent on the stress level (p). The mean effective stress (p) is normalized by atmospheric pressure (). The shear and bulk moduli can be expressed aswhere is the current void ratio, and are the shear and bulk modulus numbers, respectively, and is equal to 101.325 kPa.
The stress ratio and dilatancy are related as follows:where and represent the incremental plastic volumetric and plastic deviatoric strains, respectively, is a model parameter, is the slope of the critical state line on the p-q plane, and is the stress ratio. Accordingly, dilatancy will equal zero when the stress ratio equals . According to the angle of internal friction at the critical state , Lode’s angle and can be expressed as
The loading direction vector is the gradient vector of the yield surface, and the plastic flow direction vector is the gradient vector of the potential surface. In the PZC model, the flow rule is assumed to be nonassociated. During the loading process, the plastic flow direction vector for loading can be expressed as
In the case of unloading, the unloading plastic flow direction vector is defined as
Generally, the plastic flow of soil has a similar expression as follows:where and and are two constant model parameters.
The loading plastic modulus is proposed as follows:where is the plastic modulus in loading, is the plastic modulus number, can limit the possible state, accounts for phase transformation, considers soil degradation, is the stress ratio parameter, and are material model constants, and is the accumulated plastic deviatoric strain .
The reloading plastic modulus is expressed as
In this paper, to consider the history of past loading conditions, the plastic modulus is modified by introducing aswhere is a model constant that has to be calibrated to provide the best prediction of the loading-reloading experiments. is a mobilized stress function expressed by
B. FE Implementation
Dynamic control equations of the interaction between fluid and solid can be transformed into spatially discretized FEM formulation based on the Galerkin weighted residual method and (5) and (6) can be discretized into where is the mass matrix, is the solid stiffness matrix, is the coupling matrix, is the seepage matrix, and is the fluid-compressibility matrix.
is the fluid-phase load vector.
is the solid-phase load vector.
Equation (A.12) is solved using a finite difference technique called the generalized Newmark (GNpj) single-step time integration scheme. At time , the acceleration , velocity , and displacement are given by
where and are coefficients of integration chosen in the range of 0 to 1.
The rate of pore pressure change and pore pressure are expressed as
The integration constants for the generalized Newmark (GNpj) [34] have been changed. Equation (B.1) is solved using a finite difference technique called the generalized Newmark (GNpj) single-step time integration scheme. The acceleration, velocity, displacement, rate of pore pressure change, and pore pressure, respectively, are given by formulas (B.2)-(B.3), where , , and are coefficients of integration chosen in the range of 0 to 1. For the schemes to be unconditionally stable, the conditions and must be satisfied. Usually some algorithmic (numerical) damping is introduced using such values as [35]: .
In terms of incremental displacements and pore pressure as primary unknowns, the matrix form of the discretized FEM formulation can be written as follows:
Data Availability
The data used to support the findings of this study are included within the article.
Conflicts of Interest
The authors declare that they have no conflicts of interest regarding the publication of this paper.
Acknowledgments
The authors would like to gratefully acknowledge Professor A. H. C. Chan for making the use of SWANDYNE program available to the authors. The research was supported by the National Major Scientific Research Program in 13th Five-Year Plan (Grant no. 2016YFB0201001).