Liquefaction Analysis using Shear Wave Velocity

The Andrus and Stokoe curves developed based on shear wave velocity case history databases, are the most widely used in the context of the Seed and Idriss simplified procedure as a deterministic model. Theses curves were developed from the database according to the calculate Cyclic Stress Ratio (CSR) proposed by Seed and Idriss in 1971 with the assumption that the dynamic cyclic shear stress (τd) is always less than the simplified cyclic shear stress (τr) deduced by Seed and Idriss based on their simplifying hypotheses (rd= τd / τr <1). Filali and Sbartai in 2017, showed that rd can in many cases be greater than 1, and they have proposed a correction for the CSR in the range where rd >1. In this paper, we will present a probabilistic study based on the Bayesian method for the evaluation of the liquefaction potential of a soil deposit using a case history database based on shear wave velocity measurement. The result of this analysis shows that by using the corrected version of the simplified method, the boundary curve is moved to a new position. Then, the objective of this study is to present an adjusted mathematical model which characterizes the new position of the boundary curve (CRR) and a new formulation for computing the probability of liquefaction based on the probabilistic shape of the CRR curves using the corrected and the original version of the simplified method.


Introduction
After the earthquakes of Alaska (1964) and Niigata in Japan (1964), Seed and Idriss [1] developed a simplified procedure based on in-situ tests to evaluate the liquefaction potential which is defined by a safety factor calculated by the ratio between the Cyclic Resistance Ratio and the Cyclic Stress Ratio (CRR/CSR). Thereafter, this procedure was modified and improved, in particular by Seed [2], Seed and Idriss [3], Seed et al. [4]. Youd et al. [5,6] in their contribution have modified the expression of the stress reduction factor (rd) to extend it whatever the depth of the soil deposit, Akhila et al. [7] have used an artificial intelligence techniques to predict the cyclic resistance ratio for clean sands. The contribution of Kuo et al. [8] in the improvement of this method were summarized in a proposed empirical simplified method to evaluate the liquefaction potential, Guoxing et al. [9] have developed from a liquefaction case history database a new mathematical model to predict the CRR curves. This procedure is based on simplifying hypothesis by considering the soil column as a rigid body with the assumption that the actual peak shear stress (τd) induced at depth, h, is always less than that predicted by the simplified procedure (τr) of Seed and Idriss (rd= τd / τr <1). Thus, Filali and Sbartai [10], in their study, showed that the dynamic cyclic shear stress (CSRD) can in many cases be greater than the Simplified Shear Stress (CSR) according to the used earthquake. This result (rd>1) was found in the study conducted by Farrokhzad [11] for many sites at significant depth and in the work presented by Sun et al. [12] at shallow depth for a few sites. Therefore, rd, can be greater than 1 (rd>1), in this case, this procedure cannot be considered as conservative, thus, the simplified procedure of Seed and Idriss [1] can't be applied because it's based on the assumption that rd<1, and all the modifications and improvements made in the literature are based on this assumption. For this reason, Filali and Sbartai [10], in order to generalize the use of the simplified procedure, have proposed a corrector factor in order to adjust the simplified CSR in the range where rd>1 which correspond to a maximum acceleration of the earthquake less than 0.30g (amax<=0.30g).
In this paper, we will present a probabilistic analysis of liquefaction potential based on the proposed correction [10] in order to define the cyclic resistance ratio (CRR) curves used to characterize the boundary between liquefied and non-liquefied regions. The residual strength of liquefied soil, or liquefied shear strength, is defined as the shear strength mobilized at large displacement after liquefaction is triggered in a saturated soil. As the liquefaction potential is governed by this parameter, several works have been published in the literature, such as the study of Fadhil and Ali [13] in which they have performed a series of laboratory tests to improve the parameters of the shear strength. The structure of this article starts with an introduction in which a literature review related to liquefaction potential is presented. Followed by the presentation of the deterministic model of the CRR and CSR used in this study, the next section treats the used probabilistic analysis and the obtained result, followed by two case studies and a comparison with previous studies to valid the obtained results, and closed by a conclusion that summarizes the results and objective of this study.

Deterministic Model
The approach of Seed and Idriss [1] is the most widely used procedure in practice for estimating the liquefaction resistance of sandy soils. To represent the ground motions caused by earthquakes with one single parameter, a simplified procedure has been developed by Seed and Idriss [1] and updated in Youd et al. [6]. The resistance to liquefaction is evaluated by comparing a property index of the soil to the Cyclic Stress Ratio (CSR) given by the following equation for a magnitude earthquake adjusted to 7.5: Where σv = the vertical total stress of the soil at the depth studied, σ′v = the vertical effective stress of the soil at the depth studied, amax = the peak horizontal ground surface acceleration, g = the acceleration of gravity and rd = the shear stress reduction factor. After Filali and Sbartai [10], as the assumption rd<1 is verified only when amax>0.30g, in other word, when amax<0.30g which correspond to rd>1, the deformable and rigid body are not adjusted in accordance with the assumption on which is based the simplified procedure, and in order to generalize the use of the simplified method by adjusting the deformable and rigid body whatever the used earthquake, the authors have proposed a new earthquake corrector factor, RC, in the range where amax<=0.30g in order to adjust the dynamic and simplified results when rd>1 and ensure the reliability of the simplified method by giving the most conservative case for all earthquakes. The proposed correction [10] is defined by an earthquake corrector factor, RC, which is the ratio between the dynamic and the simplified shear stress expressed as follows:

Cyclic Resistance Ratio (CRR)
Andrus and Stokoe [14] have collected a case history data based on Vs measurement, at over 50 sites (124 test arrays) and field performance data from 20 earthquakes, including a total of 193 cases for liquified and non-liquified sites. From this database, they have developed a bounding curve for fine content (FC) ≤5%, 20% and ≥35% and proposed a correlation between Vs1 and CRR expressed as: Where Vs1* is the critical value of Vs1 equal to 220 m/s for FC≤5%, 210 m/s for FC≈20% and 200 m/s for FC≥35% , a and b are curve fitting parameters equal to 0.03 and 0.9 respectively to adjust the curves in the limit which separate liquefied and non-liquefied cases for earthquake magnitude of 7.5. Andrus et al. [15] For clean sand, the proposed Andrus and Stokoe [16] relationship based on Vs-CRR curves for non-cemented soil, Holocene age with different percentages of fines shown in the following equation:.
Where (Vs1) Cs is the overburden stress-corrected shear wave velocity defined as follows: Where Vs1 is the overburden stress-corrected shear wave velocity of sandy soils, pa is the reference stress of 100 kPa, and Kcs is a fines content (FC) correction factor. Kcs can be estimated by the relationships proposed by Juang et al. [17]: The CRR7.5 should be corrected for the earthquake magnitude, overburden pressure, and static shear [3, 18,19].
Where MSF is the magnitude scaling factor, and Kσ and Kα are factors for overburden and initial static stress ratio corrections, respectively. These factors are calculated by the formulae recommended by Boulanger and Idriss [19].

Magnitude Scaling Factor (MSF)
Several equations have been proposed for the assessment of MSF according to the earthquake moment magnitude [3,20]. Idriss [20] proposed the magnitude scaling factor as:

Overburden Correction Factor, Kσ
The overburden correction factor Kσ can be estimated by the relationship proposed by Boulanger and Idriss [19]: Where the coefficient Cσ can be expressed in terms of corrected shear wave velocity.

Static shear stress correction factor kα
To take into account the influence of static shear stresses on CRR, Seed et al. [18] have proposed a correction factor Kα to correct the CRR. Several researches were conducted by Idriss and Boulanger [21,22]. The author believes that these results can be used. As the soil layers are considered horizontal, the value of Kα in this study is kept equal to 1

Bayesian Mapping Function (BMF)
Since the deterministic safety factor (Fs) is the most widely used in Geotechnical practice, it is interesting to relate it to the probability of liquefaction in order to facilitate the use of the probabilistic approach for engineers for making an adequate decision. Juang et al. [23] have proposed a mapping function approach which linked the deterministic Fs to the probability of liquefaction; this approach has been refined by Juang et al. [24,25]. In this approach, the conditional probability of liquefaction for a given site is deduced from the information contained in the case history database [17,25] in according to the following equation: Where fL(FS) and fNL(FS) are the probability density functions of the calculated FS for the sets of liquefied cases and non-liquefied cases, respectively. Based on the obtained Equation 14, the probability of liquefaction is calculated for each of the 225 cases in the database using the original and the corrected version of the simplified procedure.

Original procedure of Seed and Idriss (1971)
The variation of the probability of liquefaction against the deterministic safety factor (Fs) calculated using the original version of the simplified procedure [1] is plotted in Figure 1 with that obtained by Juang et al. [26]. The set of the 225 points can be fitted in terms of mapping function which linked PL to Fs defined by the following equation: The deterministic curve model is defined by Fs=1. Thus, the Andrus et al. [15] curves can be characterized with a probability of liquefaction of 30% based on Byesian mapping model, this result is similar to that obtained by Juang et al. [26] from the same database. From Equation 15, we can plot for a given value of PL the CRR boundary curves presented in Figure 2. This figure shows that the value of VS1CS converges to 215m/s for high values of CSR.

Corrected version of the simplified procedure (Filali and Sbartai, 2017)
The safety factor calculated using Equation 4, for the CSR and Equation 7, for the CRR is used to recalculate the probability of liquefaction. By fitting the set of points presented on Figure 3, the Mapping function can be expressed by the relationship below: In this equation, a value of Fs=1 corresponds to the deterministic curve model. Therefore, for this case, the Andrus et al [15] curves can be characterized with a probability of liquefaction of 11.3% based on Bayesian Mapping model. From Equation 16, we can plot for a given value of PL, the CRR boundary curves presented in Figure 4. This figure shows that the value of VS1CS converge to 215m/s for high values of CSR, this result is also the same as that obtained by Juang et al. [15] according to their Bayesian mapping model.
The figure also shows that the boundary curve proposed by Andrus and Stokoe [16] is characterized by a PL=0.113, and in according to the corrected version of the simplified procedure, it is not conservative because it cannot be considered as a boundary curve, which separate the liquefied and non-liquefied cases and must be adjusted to the curve corresponding to PL=0.30 which represents the true boundary between the two zones.
By fitting this curve using the Andrus and Stokoe [16] model, the cyclic resistance ratio, CRR7.5, can be expressed by the following equation:  By comparing this equation with that proposed by Andrus and Stokoe [16], we can say that only the curve fitting parameters have changed (a=0.03433, b=4.369). This result is reasonable, because, according to the corrected version of the simplified procedure, the values of the cyclic stress ratio, CSR, have changed, therefore, the boundary between the liquefied and nonliquefied cases may also change and the curve fitting parameters must be adjusted to the new position of the boundary. Since, the mathematical model of the true boundary is defined, we must recalculate the safety factor and the probability of liquefaction for all cases in the database using Equation 4 and 17. By the same manner, the mapping function deduced by fitting the set of points presented in Figure 5 can be expressed as follows: This equation shows that the deterministic boundary curve, which correspond to FS=1 is characterized by a probability of liquefaction of 24% instead 30%. The set of probabilistic boundary curves deduced from Equations 17 and 18 are plotted in Figure 6. Then, the deterministic design decision is always made based on the safety factor which indicates that the liquefaction occur or no according to a reference value by choosing the most conservative case. The liquefaction boundaries plotted in Figure 6 show that the Andrus and Stokoe [16] CRR curve is characterized by a probability of 6% using the Bayesian mapping function with the deterministic model given by Equation 17 based on the corrected simplified method which correspond to a deterministic safety factor (FS) of 1.52, while the adjusted model proposed in this study shown in Equation 17 is related to a probability of 24%, which correspond to FS=1. Then, according to these results, the more conservative case is always given by the corrected simplified method. Thus, the set of curves shown in Figure 6 indicates a non-liquefaction for the zone above the boundary curve of PL=90% and a liquefaction for the zone below the boundary curve of PL=6%. The zone between PL=6% and PL=90% is an intermediate zone in which 6% and 24% represent the lower and the marginal probabilities, above 24% the risk of liquefaction increase with the probability of liquefaction. To define the severity of the liquefaction potential using a probabilistic analysis, Juang et al. [27] have proposed a liquefaction likelihood classification which can be used for probabilistic design decision using the corrected version of the simplified method. Comparison with previous studies The correlation between shear wave velocity and liquefaction resistance has been studied by several authors. Based on field performance data from sites in Imperial Valley, Robertson et al. [28] developed a liquefaction resistance curve where the shape was based on an analytical result. Kayen et al. [29] have conducted a probabilistic analysis to provide an unbiased assessment of Vs-based in situ soil liquefaction triggering potential, and assess the probability of liquefaction triggering for use in performance-based engineering applications. Guoxing et al. [9], based on an updated calibration using an expanded database for evaluation of the soil liquefaction potential at nuclear power plant project sites with extremely high-risk potential have proposed a deterministic empirical liquefaction triggering correlation with its probabilistic version. A comparison between these liquefaction resistance correlations and the adjusted model proposed in this study are plotted on Figure 7. This figure shows clearly that the best fit is given by the corrected version of the simplified method which materialize the true liquefaction boundary expressed by Equation 17, the four other correlations in Figure 7 cannot be considered as boundary curves because according to the corrected version of the simplified method, the values of CSR for all cases in the database for which amax≤0.30g are adjusted through a corrector factor, RC, defined by Equation 3. Therefore, the plotted set of points of the case history database is translated upwards, which leads to move the boundary curve to a new position different to that defined by the original simplified method. Then, in this figure, we have kept the original position of the plotted curves in order to show the effect of the proposed correction on these curves.

Treasure Island Site
Treasure Island is a man-made island located between San Francisco and Oakland. More than 29 million cubic yards of finely to medium-grained sand have been dredged from borrow areas in the San Francisco Bay and used as fill  In this example, we will evaluate the liquefaction potential with the original and the corrected version of the simplified procedure in order to define which of the two methods gives the more conservative case. Then, the cyclic stress ratio CSR is calculated by using both Equations 1 and 4; for the estimation of the cyclic resistance ratio, we will use Equation 5 adjusted to fines content ≤5% with fitting parameters of Andrus and Stokoe [16] and Equation 17.The peak ground acceleration amax value used for the calculation of CSR is taken equal to 0.1129g. The depth to the ground water table is kept 1.5m relative to the ground surface. The average value of the unit weight is taken equal to 17.6 KN/m 3 above the water table and 19.2 KN/m 3 below the water table below ,the third and the fifth layers (Clayey sand) may be non-liquefiable by Chinese criteria [16]. In Figure 8, are shown the profiles of safety factor according to the depth computed by the original and the corrected version of the simplified procedure, the profile of the dynamic FS is deduced from a dynamic analysis with lumped mass discussed in Filali and Sbartai [10] by using equation 17 to estimate the cyclic resistance ratio.  Figure 8 shows that the more conservative case is given by the corrected version of the simplified procedure, and the profile of the corrected safety factor is very close to the dynamic profile. These results indicate that the maximum shear stress given by the corrected version which is almost equal to that computed from a dynamic analysis is always for this case greater than the shear stress estimated by the original simplified method, which implies that the stress corrector factor, rd, is greater than 1, and to confirm this, we have conducted a dynamic analysis using Shake91_input software [30] in which the Loma Prieta earthquake is simulated by the DIAM accelerogram applied at the bottom of the soil profile. In this analysis, we have calculated the maximum shear stress for soil profile using Shake91_input and the simplified method with the original and the corrected version using the maximum acceleration of DIAM accelerogram which is 0.1129g. The results are presented in Figure 9. This figure shows clearly that the maximum shear stress computed by the original simplified method is less than that given by the dynamic analysis conducted using Shake91_input (rd>1, amax<0.30g) while the corrected version of the simplified method gives values greater than or equal to those of the dynamic method (rd≤1). Then, for this site, the liquefaction potential evaluation must be conducted using the corrected version of the simplified method, because the original version cannot be applied since rd is greater than 1.

Petrochemical Zone of Skikda Site (Algeria)
Based on the request of the National Petroleum Refining Company of Skikda department (NAFTEC), the laboratory has performed a geophysical investigation with three down-hole tests. The study site is located within the industrial zone of Skikda, it has a flat topography. The down-hole test SC02 detected the presence of a sandy horizon, reddish to brownish which extends up to depth 20 m and saturated with a mean diameter D50 varying between 0.11 and 1 mm. The average value of the unit saturated weight is taken between 19.6 and 20.5 KN/m 3 . The water table is assumed on the ground surface. The magnitude of the earthquake is in the range of 6.8 and the maximum acceleration at the surface is equal to 0.122g, the site is classified in the zone II according to the Algerian earthquake code RPA 2003. The profile of soil and shear wave velocity chosen in this study are shown on the Figure 10. In Figure 11, are shown the profiles of safety factor according to the depth computed by the original and the corrected version of the simplified procedure, the profile of the dynamic FS is deduced from a dynamic analysis performed with Shake_input software [30] in which the dynamic cyclic stress ratio (CSRD) was expressed as the ratio of the maximum shear stress and the vertical effective stress. For this site, the conclusion is the same as the Treasure Island site. To confirm this, we have conducted a dynamic analysis using Shake91_input software [30] in which the Boumerdes earthquake of 21/05/2003 is simulated by the Azazga station accelerogram EW component applied at the bottom of the soil profile. In this analysis, we have calculated the maximum shear stress for soil profile using Shake91_input and the simplified method with the original and the corrected version using the maximum acceleration of the used accelerogram which is 0.122g. The results are presented in Figure 12. This figure shows clearly that the maximum shear stress computed by the original simplified method is less than that given by the dynamic analysis conducted using Shake91_input (rd>1, amax<0.30g) while the corrected version of the simplified method gives values greater than or equal to those of the dynamic method (rd≤1).

Conclusion
In this study, we have critically compared the models of the probability of liquefaction obtained by the original simplified method [1] and the corrected version of this method [10] by using a Bayesian mapping function based on shear wave velocity test. The results show that the boundary curve is characterized in one hand, by PL=0.30 which correspond to Fs=1 by using the original simplified method, and in other hand, by PL=0.24 which correspond to Fs=1 by using the corrected version of this method with the Andrus and Stokoe [16] shape of the CRR expressed by Equation 7. Then, the proposed model for the CRR curve of Andrus and Stokoe [16] must be adjusted to the new boundary in accordance with the corrected version of the simplified method because the boundary curve is obtained by plotting CSR against Vs from the case history data, and, as the CSR have changed for all sites in the database where amax<=0.30g, the boundary curve must also change and may be readjusted. This readjustment is materialized by the proposed Equation 17 from which the CRR curve is positioned on the true boundary which separate the liquefied and non-liquefied zones according to the corrected version of the simplified method. This correction is only valid for clean sand, then, other sands where FC>5% must be adjusted to clean sand according to VS1CS in order to be able to use the proposed correction. The case studies of Treasure Island and the petrochemical zone of Skikda sites show clearly that by using the corrected version of the simplified method the assumption rd<1 is always verified whatever the used earthquake and is not verified when the original simplified method is used for the two sites.

Conflicts of Interest
The authors declare no conflict of interest.