Geostatistical analysis of spatial variability of the liquefaction potential – Case study of a site located in Algiers (Algeria)

The cyclic stress ratio (CSR) is calculated according to the following formula^[27]: (1) CSR=τcycσV0′=0.65(σV0σV0′)(amaxg)rd(1MSF) {\rm{CSR}} = {{{\tau _{{\rm{cyc}}}}} \over {{{\sigma_{{\rm{v0}}}^\prime}}}} = 0.65\left( {{{{\sigma _{{\rm{v}0}}}} \over {{{\sigma_{{\rm{v}0}}^\prime}}}}} \right)\left( {{{{{\rm{a}}_{\max }}} \over {\rm{g}}}} \right){{\rm{r}}_{\rm{d}}}\left( {{1 \over {{MSF}}}} \right) where a_max is the maximum horizontal acceleration at the ground surface, g is the acceleration of gravity, σ_v0 and σ^′_v0 are the total and effective vertical stresses, respectively, r_d the stress reduction factor, and MSF is the Magnitude Scaling Factor.

To calculate the factor r_d, [31] proposed the following relationship: (2) rd=(1−0.4113z0.5+0.04052z+0.001753z1.5)(1−0.4177z0.5+0.05729z−0.0062z1.5+0.00121z2) {r_d} = {{\left( {1 - 0.4113{z^{0.5}} + 0.04052z + 0.001753{z^{1.5}}} \right)} \over {\left( {1 - 0.4177{z^{0.5}} + 0.05729z - 0.0062{z^{1.5}} + 0.00121{z^2}} \right)}} where z is depth below ground surface (m).

With regard to the Magnitude Scaling Factor (MSF), the lower bound equation, which was suggested by [31], was used as follows: (3) MSF=102.24Mw2.56 {MSF} = {{{{10}^{2.24}}} \over {M_w^{2.56}}} where Mw is the moment magnitude of the earthquake.

The cyclic resistance ratio (CRR) is estimated from CPT data following the procedure in [23], [24] as: (4) CRR={0.833(qc1N,cs1000)+0.05 if qc1N,cs<5093(qc1N,cs1000)3+0.08 if 50≤qc1N,cs<160 {\rm{CRR}} = \left\{ {\matrix{ {0.833\left( {{{{{\rm{q}}_{{\rm{c}}1{\rm{N}},{\rm{cs}}}}} \over {1000}}} \right) + 0.05\,{\rm{if}}\,{{\rm{q}}_{{\rm{c}}1{\rm{N}},{\rm{cs}}}} < 50} \cr {93({{{{{{\rm{q}}_{{\rm{c}}1{\rm{N}},{\rm{cs}}}}}}} \over {1000}})^3 + 0.08\,\,{\rm{if}}\,50 \le {{\rm{q}}_{{\rm{c}}1{\rm{N}},{\rm{cs}}}} < 160} \cr } } \right.

The resistance to penetration of a standard cone may be calculated using the following equation: (5) (qc1N)cs=Kc×qc1N {\left( {{q_{c1N}}} \right)_{cs}} = {K_c} \times {q_{c1N}} K_c is a correction factor that depends on I_c and defined by the equation of [23] as: (6) { Kc=1.0 if Ic≤1.64Kc=−0.403 Ic4+5.58 Ic3−21.63 Ic2+33.75Ic−17.88 if Ic>1.64 \left\{ {\matrix{ {\,{{\rm{K}}_{\rm{c}}} = 1.0\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,\,{\rm{if}}\,{{\rm{I}}_{\rm{c}}} \le 1.64} \hfill \cr {{{\rm{K}}_{\rm{c}}} = - 0.403\,{\rm{I}}_{\rm{c}}^4 + 5.58\,{\rm{I}}_{\rm{c}}^3 - 21.63\,{\rm{I}}_{\rm{c}}^2 + 33.75{{\rm{I}}_{\rm{c}}} - 17.88\,{\rm{if}}\,{{\rm{I}}_{\rm{c}}} > 1.64} \hfill \cr } } \right. On the other hand, the corrected peak resistance q_c1N (CPT) was determined using the following formula: (7) qc1N=(qc−σvoPat)(Patσv0′)n {q_{c1N}} = \left( {{{{q_c} - {\sigma _{v0}}} \over {{P_{at}}}}} \right){\left( {{{{P_{at}}} \over {{{\sigma_{v0}^\prime}}}}} \right)^n} The exponent n is defined using the formula of [24] as: (8) n=0.381(Ic)+0.05(σv0′Pat)−0.15 where n≤1 n = 0.381\left( {{I_c}} \right) + 0.05\left( {{{{{\sigma_{v0}^\prime}}} \over {{P_{at}}}}} \right) - 0.15\,{where}\,n \le 1 Note that P_at is the reference pressure (100 kpa), q_c is the peak resistance, and σ_v0 and σ^′_v0 are the total and effective vertical stresses, respectively.

The soil behavior type index I_c is defined by the expression proposed by [23]: (9) Ic=(3.47−logQ)2+(1.22+logF)2 {I_c} = \sqrt {{{\left( {3.47 - {log}Q} \right)}^2} + {{\left( {1.22 + {log} F} \right)}^2}} With (10) Q=[(qc−σV0)σV0′] Q = \left[ {{{\left( {{q_c} - {\sigma _{V0}}} \right)} \over {{{\sigma_{V0}^\prime}}}}} \right] (11) F=⌈fs(qc−qV0)⌉×100% F = \left\lceil {{{{f_s}} \over {\left( {{q_c} - {q_{V0}}} \right)}}} \right\rceil \times 100\% The quantity Q is known as the normalized tip resistance and F is the normalized friction ratio of the penetrometer sleeve.

The local factor of safety (FS) coefficient for the soil liquefaction is such that^[12]: (12) Fs=CRRCSR {F_s} = {{{CRR}} \over {{CSR}}} Note that liquefaction is predicted to occur if FS < 1; however, no liquefaction is expected to take place if FS > 1.^[12] It is worth mentioning that the soil becomes more resistant to liquefaction as the safety factor increases.

The first step consists of determining the safety factor (FS), and the second one is to find the liquefaction potential index (LPI), which is a measure of the total risk of liquefaction of soil to a depth of 20 m. The liquefaction potential index is defined as follows^[11,12]: (13) LPI=∫020F(z)w(z)dz {LPI} = \int_0^{20} {F\left( z \right)} w\left( z \right)dz According to [17] and [16], the liquefaction potential index can be calculated using the following expression: (14) LPI=∑i=1NLwiFLiHi {LPI} = \sum\limits_{i = 1}^{NL} {{w_i}{F_{Li}}{H_i}} where w_i(z) = 10 - 0.5z, and z denotes the depth in meters; Hi is the thickness of the discretized soil layers, NL is number of soil layers. Then, F_Li may be defined as [17]: (15) FLi={ 0for Fs(z)≥1,01−Fsfor Fs(z)≤1 {F_{Li}} = \left\{ {\matrix{ {\,\,\,\,\,\,\,\,\,\,\,0} \hfill & {{for}} \hfill & {\,{F_s}\left( z \right) \ge 1,0} \hfill \cr {1 - {F_s}} \hfill & {} \hfill & {{for}\,\,\,{F_s}\left( z \right) \le 1} \hfill \cr } } \right. Furthermore, [28] classifies a site, according to the importance of the soil liquefaction phenomenon, into five categories: Very low for LPI = 0, Low for 0 < LPI ≤ 2, Moderate for 2 < LPI ≤ 5, High for 5 < LPI ≤ 15, and Very high for LPI > 15.

Table 1 presents an example about the evaluation of the liquefaction potential index for a typical borehole (BH34). Figure 2 shows some illustrative graphs to show the major steps involved in calculating the liquefaction potential at one CPT sounding location. The liquefaction potential index was calculated for each CPT profile for an earthquake scenario with Mw = 6.8 and having a peak horizontal ground acceleration a_max = 0.3 g (Table 2). The CPT positions were defined by the Universal Transverse Mercator (UTM 31 Zone 31N) coordinate system.

Figure 2

The variogram is a tool that allows measuring the dissimilarity between values based on their separations. It describes the spatial continuity of the regionalized variable.^[8] The theoretical variogram is defined as: (16) γ(h)=12Var[Z(x)−Z(x+h)]=12E[(Z(x)−Z(x+h))2] \gamma\left( h \right) = {1 \over 2}{\rm{Var}}\left[ {Z\left( x \right) - Z\left( {x + h} \right)} \right] = {1 \over 2}{\rm{E}}[ {{{\left( {Z\left( x \right) - Z\left( {x + h} \right)} \right)}^2}}] where Z(x) is an observed value at a particular location, h is the distance between ordered data, Var is the variance of a random variable, E is the expectation operator of a random variable.

The experimental variogram, proposed by [18], is defined as follows^[8]: (17) γe(h,θ)=12N(h,θ)∑i=1N(h,θ)[Z(xi)−Z(xi+h)]2 {\gamma_e}\left( {h,\theta } \right) = {1 \over {2N\left( {h,\theta } \right)}}\sum\limits_{i = 1}^{N\left( {h,\theta } \right)} {{{\left[ {Z\left( {{x_i}} \right) - Z\left( {{x_i} + h} \right)} \right]}^2}} Note that x_i is an experimental point; (x_i + h) is an experimental point located at distance h from xi, in direction θ; Z (x_i) is the value measured at point x_i; Z (x_i + h) is the value measured at point (x_i + h); N (h, θ) is the number of pairs of observations in the direction considered θ and separated by distance h.

Ordinary kriging is an interpolation procedure that provides, from the available data, an optimal, linear and unbiased estimate of the property under study and whose estimation error is minimized.^[6] The interpolated value at the point x₀, denoted Z* (x₀), is given by: (18) Z*(x0)=∑i=1NλiZi(x) {Z^*}\left( {{x_0}} \right) = \sum\limits_{i = 1}^N {{\lambda _i}{Z_i}\left( x \right)} The weights λ_i are chosen so as to minimize the estimation variance. This should lead to the following ordinary Kriging equations ^{[19, 7, 29]}: (19) {∑i=1,nλiγ(xi,xj)+μ=γ(xj,x0)(i, j=1,2,...n)∑i=1Nλi=1 \left\{ {\matrix{ {\sum\limits_{i = 1,n} {{\lambda _i}\gamma \left( {{x_i},{x_j}} \right) + \mu = \gamma \left( {{x_j},{x_0}} \right)\left( {i,\,j = 1,2,...n} \right)} } \cr {\sum\nolimits_{i = 1}^N {{\lambda _i} = 1} } \cr } } \right. The kriging system of equation (19) will take the following matrix form: (20) (γ11γ(x1,x2)⋅γ(x1,xn)1γ(x2,x1)γ22⋅γ(x2,xn)1⋅⋅⋅⋅⋅γ(xn,x1)γ(xn,x2)⋅γnn111⋅10)(λ1λ2⋅λnμ)=(γ(x1,x0)γ(x2,x0)⋅γ(xn,x0)1) \left( {\matrix{ {{\gamma _{11}}} & {\gamma \left( {{{\rm{x}}_1},{{\rm{x}}_2}} \right)} & \cdot & {\gamma \left( {{{\rm{x}}_1},{{\rm{x}}_{\rm{n}}}} \right)} & 1 \cr {\gamma \left( {{{\rm{x}}_2},{{\rm{x}}_1}} \right)} & {{\gamma _{22}}} & \cdot & {\gamma \left( {{{\rm{x}}_2},{{\rm{x}}_{\rm{n}}}} \right)} & 1 \cr \cdot & \cdot & \cdot & \cdot & \cdot \cr {\gamma \left( {{{\rm{x}}_{\rm{n}}},{{\rm{x}}_1}} \right)} & {\gamma \left( {{{\rm{x}}_{\rm{n}}},{{\rm{x}}_2}} \right)} & \cdot & {{\gamma _{nn}}} & 1 \cr 1 & 1 & \cdot & 1 & 0 \cr } } \right)\left( {\matrix{ {{\lambda _1}} \cr {{\lambda _2}} \cr \cdot \cr {{\lambda _{\rm{n}}}} \cr \mu \cr } } \right) = \left( {\matrix{ {\gamma \left( {{{\rm{x}}_1},{x_0}} \right)} \cr {\gamma \left( {{{\rm{x}}_2},{x_0}} \right)} \cr \cdot \cr {\gamma \left( {{{\rm{x}}_{\rm{n}}},{x_0}} \right)} \cr 1 \cr } } \right) where n is the number of information available; x_i and x_j are measurement points; Z (x_i) and Z (x_j) are the values measured at points x_i and x_j; x₀ is the point (or volume) to be estimated; µ is the Lagrange parameter; γ (x_i, x_j) is the value of the variogram γ (h) for the distance h between x_i and xj; γ (x_j, x₀) is the value of the variogram γ (h) for the distance h separating x_i and x₀.

Kriging also makes it possible to evaluate the uncertainty associated with the estimator in the form of the kriging variance, which is given as: (21) σe2=∑i=1,nλjγ(xi,x)−γ¯(x,x)+μ \sigma _e^2 = \sum\limits_{i = 1,n} {{\lambda _j}\gamma \left( {{x_i},x} \right) - \bar \gamma \left( {x,x} \right) + \mu } where γ¯(x,x)γ¯(x,x) \bar \gamma \left( {x,x} \right)\bar \gamma \left( {x,x} \right) is the mean value of the variogram between two points belonging to the volume x. Note that if the volume is reduced to a point, then the distance between these points is zero.