Seismic Fragility Analysis of Retaining Walls Dependent on Initial Conditions

Abstract

Fragility curves of retaining walls constitute an efficient tool for the estimation of seismic risk and can be utilized for prevention from potential damage or for immediate decision-making. In this work, fragility curves for cantilever retaining walls of three different heights are proposed, considering cohesionless soil materials. The seismic response of the soil-wall system, in terms of permanent vertical ground displacement of the backfill soil and permanent horizontal displacement of the wall’s base, is estimated by conducting non-linear time history analyses, through the 2D finite element simulation method. Five initial conditions are investigated regarding the value of the global factor of safety (F_S) under static conditions. An initial value of F_S equal to 1.5 is considered for dry conditions. If the presence of the water table is taken into account, the corresponding F_S drops to values ranging from 1.4 to 1.1. Parameters that characterize seismic intensity are evaluated based on criteria, in order to identify the intensity measures that best correlate with the system’s response. Three damage states are adopted, corresponding to minor, moderate, and extensive damage. The approach of combined damage criteria is also investigated. Finally, fragility curves are derived demonstrating the degree of dependency on initial conditions.

1. Introduction

The impact of seismic events on retaining walls along roadway axes can be of major importance, as induced damage can affect the orderly road operation directly and in the long term, not to count fatalities, human injuries, and loss of properties. According to the Geotechnical Extreme Events Reconnaissance (GEER) database, deformations and failures of retaining structures have been reported in historical earthquake events of the past, like the Kobe, Chi-Chi, Iquique, and Cephalonia earthquakes [1,2,3,4], but also in more recent events, such as the earthquake sequence in central Italy in 2016 [5], and Kahramanmaras (Turkey) in 2023 [6]. Therefore, a seismic fragility assessment of such elements is a critical factor that contributes to the prevention and mitigation of the direct loss corresponding to the physical damages of the road network, and the indirect loss related to the long-term reduced functionality [7].

The term of seismic fragility refers to the tendency of an element or system to perform inadequately when subjected to a seismic event. Fragility curves are utilized as a probabilistic measure to assess the seismic performance of an element exposed to risk. They describe the probability of reaching or exceeding predefined damage states, as a function of the intensity of ground shaking. Different methods can be applied for the derivation of fragility curves, namely empirical, based on expert judgment, numerical, and hybrid methods. Empirical fragility curves are derived by collecting data from post-earthquake observations. The difficulties involved in this approach include the comprehensive amount of data required from many seismic events and the estimation of spatial distribution of seismic intensity at various locations. Fragility curves based on expert judgement (e.g., HAZUS, 2004 [8]) are based on experts’ opinion regarding the possible damage of a structure, induced by certain levels of seismic intensity. This approach is characterized by a high level of subjectivity difficult to quantify [9]. Recently, numerical analyses have been widely used, as they provide the ability to model and estimate the seismic response of any type of structure. Finally, hybrid methods result from the combination of the previously mentioned approaches, as they use empirical data from past earthquakes in order to optimize relationships obtained by analytical or judgment-based methods.

In seismic fragility analyses, the possible damage states of the structure are bounded by some predefined limit states. Damage levels, as well as damage indices, are generally linked to the seismic-induced displacements of the structure. In the case of retaining walls, a first reference to failure criteria by defining the horizontal displacement of the wall as a percentage of its height has been made by Prakash et al. (1995) [10]. They proposed the percentage 2% of wall height as a criterion for permissible horizontal displacements, whereas 10% of wall height as a failure criterion. Following this proposal, Huang et al. (2009) [11], investigated the permissible displacements of conventional retaining walls in the scope of internal friction angle mobilization along the potential failure line of the cohesionless backfill. They conducted shaking table tests, as well as seismic displacement analyses utilizing Newmark’s sliding block theory. The results demonstrated that for horizontal displacements up to 2% of wall height, the backfill soil retains its maximum strength, while for displacements greater than 5% of wall height, the friction angle reaches a reduced value, and the soil undergoes residual deformations. The above-mentioned criteria have been introduced to seismic fragility analyses, by means of a damage index. Zamiran & Osouli (2018) [12] adopted the horizontal displacement of the wall as a damage index, assuming the percentages 2%, 5% and 8% of wall height as limit values corresponding to the exceedance of minor, moderate, and extensive damage levels. Cosentini & Bozzoni (2022) [13] and Seo et al. (2022) [14] adopted the limit values of 2%, 5% and 10% for the three damage states, respectively.

Other parameters adopted in the literature as damage indices for fragility analyses of retaining walls are the peak vertical ground displacement of the backfill and the angle of rotation of the wall. Argyroudis et al. (2013) [15] and Seo et al. (2022) [14], used the threshold values of 0.05 m, 0.15 m, and 0.40 m as vertical displacement limit states for minor, moderate, and extensive damage levels, according to the criteria proposed by the research program SYNER-G (2013) [16]. In those cases, the damage states were qualitatively correlated with the serviceability level of the road. Furthermore, the rotation angle was adopted for gravity walls by Cosentini & Bozzoni (2022) [13], calculated as the arc tangent of the ratio of the horizontal displacement to the total wall height. The predefined limit states were equal to the damage criteria proposed by PIANC (2001) [17] for quay walls, i.e., 3° for minor, 5° for moderate, and 8° for extensive damage.

So far, much research has been carried out regarding the derivation of fragility curves for retaining walls along roadways, by using the finite element method (FEM) or finite difference method (FDM) software. Argyroudis et al. (2013) [15] proposed fragility curves for bridge abutments of 6 m and 7.5 m height for soil classes C and D, according to Eurocode 8 [18] by conducting dynamic analyses in FEM software PLAXIS 2D [19]. The proposed fragility curves were compared to observed damage during the Niigata-Chuetsu Oki 2007 earthquake. Zamiran & Osouli (2018) [12] evaluated the seismic response of a 6 m height cantilever retaining wall, considering three different configurations of the backfill cohesion, i.e., 0 kPa, 10 kPa, and 30 kPa. FDM software FLAC 2D [20] was used to determine the seismic response of the wall. The developed fragility curves denoted to what level a difference in cohesion can affect the probability of damage. Cosentini & Bozzoni (2022) [13] proposed fragility curves for two different shapes of gravity retaining walls of 3.6 m, 8 m, and 12 m height, and two configurations of the backfill inclination (0° and 30°). Dynamic analyses were carried out with FDM software FLAC 2D [20]. Following an investigation, peak ground acceleration (PGA) and peak ground velocity (PGV) were selected as the seismic intensity measures that best correlated with the wall’s response, in terms of permanent horizontal displacement and rotation of the wall. Seo et al. (2022) [14], evaluated seismic fragility of 4 m height cantilever retaining walls, for 0°, 10°, and 20° slope angle of the backfill. Considering the four soil classes of the South Korean classification system, site response analyses were carried out to estimate the corresponding ground motions characteristics. The seismic response of the soil-wall system was computed by conducting 2D FDM analyses (FLAC 2D [20]). Fragility curves were proposed as a function of the peak ground acceleration (PGA) and cumulative absolute velocity (CAV).

In the existing research on seismic fragility assessment of retaining walls, only dry backfill conditions have been considered, therefore assuming that the material is well drained. In practice, malfunctioning of drains due to a blockage of the drainage pipe and inadequate maintenance is usually observed. In these cases, the presence of water leads to the increase of lateral pressure against the wall, resulting in larger displacements and a decrease of the safety factor (F_s) as the water table rises to higher levels [21]. This study focuses on the development of novel fragility curves for cantilever retaining walls with cohesionless foundation soil and backfill materials, considering five different initial conditions regarding the value of global factor of safety (F_s) under static conditions. The variation in the values of F_s is due to the change of the water table level behind the wall. In particular, the dimensions of three cantilever retaining walls with different heights are determined, provided that they have a value of F_s = 1.5 under dry conditions. Subsequently, the presence of water table in different levels behind the wall is considered, resulting in the drop of F_s to values equal to 1.4, 1.3, 1.2, and 1.1. Hence, the effect of different initial conditions to the fragility of the soil-wall system is examined, adopting as damage criteria the permanent vertical ground displacement of the backfill soil, the permanent horizontal displacement of the wall’s base, and the approach of combined criteria. An investigation to identify the seismic intensity parameters that satisfactorily correlate with the system’s response is also carried out.

2. Materials and Methods

The objective of this article is the generation of fragility curves for cantilever retaining walls, while investigating to what extent the presence of water table in various levels affects the seismic fragility. Towards this aim, the working stages included (a) the generation of the FE soil-wall system models by selecting the appropriate geometry and parameters of the constitutive models; (b) the 2D dynamic analyses considering five different initial conditions; and (c) the derivation of fragility curves after estimating the requisite parameters of the fragility function. A thorough description of the working stages is presented below.

In order to estimate the seismic response of the cantilever retaining wall, the FEM software PLAXIS 2D [19] was utilized. The initial step involved the creation of the soil-wall model, as shown in Figure 1. The 2D model consists of the cantilever retaining wall, with a geometry within the range of common practice, the backfill soil behind the wall, at the toe of the wall and at its foundation, the soil material that is an extension of the backfill and extends to one meter below the footing level (layer 1), and the underlying soil layer with a thickness of 7 m (layer 2). Bedrock is simulated at the base of the model, so that the compliant base boundary condition required for the dynamic analyses can be applied. Soil layers were considered as non-cohesive materials, with the properties presented in

Table 1. Main parameters of soil and concrete.

	Layer 1	Layer 2	Backfill	Wall
φ (°)	30	36	32	-
ψ (°)	0	6	2	-
γunsat (kN/m3)	19.50	21	20	25
γsat (kN/m3)	21	22.5	21.5	-
Vs (m/s)	200	290	230	2215

.

A total model width of 60 m was selected, to ensure that the influence of vertical boundaries is limited. Regarding the dimensions of the wall, a parameterization was performed for three different heights of 3 m, 6 m, and 9 m. Opting for wall dimensions that would ensure stability under dry conditions, the “safety calculation” analysis provided by PLAXIS 2D [19] was utilized. This type of analysis implements the phi/c reduction method, that reduces the shear strength parameters tanφ and c until failure occurs. Hence, trial investigations were carried out to identify the exact wall dimensions of the three models R.W.3, R.W.6, and R.W.9, that would result in a value of the global factor of safety under dry conditions equal to 1.5. A schematic view of the retaining wall is depicted in Figure 2, while the selected dimensions of R.W.3, R.W.6, and R.W.9 are presented in

Table 2. Dimensions of the retaining wall models R.W.3, R.W.6, and R.W.9.

Dimensions (m)	R.W.3	R.W.6	R.W.9
H	3.00	6.00	9.00
B	2.30	4.70	7.00
HR	2.50	5.40	8.10
tb1	0.20	0.40	0.50
tb2	0.10	0.20	0.40
ts1	0.25	0.30	0.30
ts2	0.30	0.60	0.90
bh	1.90	3.92	5.83
bt	0.40	0.78	1.17

.

A retaining wall can experience a seismic event anytime during its life. It is not unusual that several years after the construction of a retaining wall (especially at provincial or national mountainous road networks) drainage systems become less efficient, and therefore, the initial conditions of safety can be seriously reduced. This phenomenon is simulated with an increase of the water level behind the retaining wall. To investigate the dependency of seismic fragility on initial conditions, four additional FE models of R.W.3, R.W.6, and R.W.9 were generated, considering the presence of water table at different levels, so that the value of F_s will drop to 1.4, 1.3, 1.2, and 1.1 accordingly. Regarding the pore water pressure, the assumption of hydrostatic conditions has been made. The precise water levels behind the wall (measured from the wall base level), as shown in

Table 3. Water table levels behind the walls and the corresponding values of global F_s, regarding the three models R.W.3, R.W.6, and R.W.9.

Fs	Level of Water Table behind the Wall (m) 1
	R.W.3	R.W.6	R.W.9
1.5	-	-	-
1.4	0.1	0.4	0.5
1.3	0.3	1.0	2.0
1.2	1.0	2.2	3.5
1.1	1.5	3.0	4.5

¹ Assuming the wall base foundation level as a reference level.

, were determined for the three retaining walls individually by conducting trial safety analyses implementing the phi/c reduction method. It is observed that for the unfavorable case of F_s = 1.1, the water table rises to a level equal to 50% times the wall height, for all three cases of retaining walls. Figure 3 indicatively presents the water levels corresponding to the four initial conditions, for the case of R.W.6.

The defined boundary conditions for the static analyses include a fully fixed base of the model, and normally fixed boundaries at the vertical edges constraining horizontal movement. For the dynamic analyses, free field boundary conditions were defined for lateral boundaries, enabling free field motion, whilst absorbing the reflected waves as they are connected to the main domain through viscous dampers. Concerning the base of the model, the boundary condition of the compliant base was adopted for the dynamic analyses. The compliant base consists of the combination of a viscous boundary in order to have a minimum wave reflection at the base of the model and allows the input of the earthquake motion.

To simulate soil’s behavior, the Hardening Soil small (HSsmall) constitutive model was implemented. HSsmall is an advanced elastic-plastic model that accounts for the soil’s high stiffness in the small strain range, the strain dependency of stiffness, and the stress dependency of stiffness and strength. In HSsmall, the Mohr–Coulomb failure criterion is introduced. The nonlinear behavior of soil is determined by a shear hardening (deviatoric) yield surface, and a cap yield surface that explains the plastic volumetric strain mostly observed in softer types of soil. The soil stiffness is simulated by means of the secant stiffness modulus E₅₀, corresponding to 50% of the ultimate deviatoric stress at the primary loading stage of the triaxial compression test, the oedometric tangent modulus, E_oed, and the elastic modulus for unloading and reloading, namely ${\mathrm{E}}_{\mathrm{ur}}$. During unloading and reloading, the behavior of soil is considered elastic, while the soil stiffness decreases nonlinearly as the number of load cycles increases. In HSsmall, the shear modulus reduction curve is characterized by the small-strain shear modulus G₀, and the shear strain γ_0.7, at which the secant shear modulus has reduced to 72.2% of G₀. The maximum value of G₀, was calculated utilizing the relationship:

$${\mathrm{G}}_0={\mathsf{\rho }\ast \mathrm{V}}_S{}^2$$

(1)

where $\mathsf{\rho }$ is the soil density in kg/m³. The decay of the shear modulus with increasing shear strain, is approximated in PLAXIS 2D by applying relationship (2) that is based on the Hardin-Drnevich [22] equation. According to Bringkreve et al. (2007) [23], the comparison of curves of different types of soil indicated that the value of γ_0.7 is between 1–2 × 10⁻⁴, whereas the value of G₀ can range from 10 × G_ur for soft soils to 2.5 × G_ur for harder types of soil, where G_ur is the shear modulus for unloading and reloading.

$${\mathrm{G}}_{\mathrm{S}}=\frac{{\mathrm{G}}_0}{1+0.385\frac{\mathsf{\gamma }}{{\mathsf{\gamma }}_{0.7}}}$$

(2)

Within the very small strain range, the response of the soil is assumed to be linear elastic. Thus, a Rayleigh damping of 5% was introduced to the model at the target frequencies f₁ = 6.7 Hz and f₂ = 19 Hz, to account for the realistic viscous characteristics of soil in the dynamic analyses. The target frequencies correspond to the first and the second natural frequencies of the soil deposit, as demonstrated by the Fourier transfer function between the surface and the bedrock level.

The values of the parameters required to determine the constitutive model HSsmall were selected from the experimental triaxial and oedometer test results found in [24]. These parameters are the reference stiffness moduli ${\mathrm{E}}_{50}^{\mathrm{ref}}$, ${\mathrm{E}}_{\mathrm{oed}}^{\mathrm{ref}},$ and ${\mathrm{E}}_{\mathrm{ur}}^{\mathrm{ref}}$ calculated for reference stress p^ref = 100 kPa, the secant shear modulus at very small strains ${\mathrm{G}}_0^{\mathrm{ref}}$, the shear strain γ_0.7 at which the secant shear modulus has reduced to 72.2% of its initial value, the coefficient m that defines the amount of stiffness dependency on stress level, and the Poisson ratio at unloading-reloading v_ur. The soil and wall properties adopted for the determination of the constitutive models are demonstrated in

Table 4. Soil and concrete properties of the dynamic analyses model.

	Layer 1	Layer 2	Backfill	Retaining Wall
Constitutive model	HSsmall	HSsmall	HSsmall	Linear Elastic
${\mathrm{E}}_{50}^{\mathrm{ref}}$ (kN/m2)	24,000	46,800	31,600	$\mathrm{E}=\mathrm{30,000}\ast {10}^3$
${\mathrm{E}}_{\mathrm{oed}}^{\mathrm{ref}}$ (kN/m2)	24,000	46,800	31,600	-
${\mathrm{E}}_{\mathrm{ur}}^{\mathrm{ref}}$ (kN/m2)	66,000	187,200	94,800	-
${\mathrm{G}}_0^{\mathrm{ref}}$ (kN/m2)	78,000	176,600	110,800	$\mathrm{G}=\mathrm{12,500}\ast {10}^3$
m	0.50	0.50	0.50	-
vur	0.20	0.20	0.20	v = 0.20
γ0.7	0.00015	0.0002	0.00015	-

.

The selection of the appropriate element size for mesh discretization is a prerequisite for the proper simulation of wave propagation through soil strata. According to international practice, the maximum dimension l_max of each element should not exceed 10% of the wavelength propagating in the medium. Considering λ_min as the wavelength, Vs as the velocity of the shear wave, and f_max as the maximum frequency of the seismic excitation, the maximum permissible element dimension of each soil material was calculated by the relation:

$${\mathrm{l}}_{\max }=\frac{{\mathsf{\lambda }}_{\min }}{10}=\frac{Vs}{10{ \ast \mathrm{f} }_{\max }}$$

(3)

Considering f_max ≈ 12 Hz, for layer 1, the maximum permissible element dimension is equal to 1.67 m, for layer 2, l_max = 2.42 m, for the backfill, l_max = 1.92 m, and for the bedrock, l_max = 6.67 m. The model was discretized, using 15-node triangular elements, that provide a fourth order interpolation for displacements and the numerical integration involves twelve stress points. To properly model the soil-wall interaction, interface elements were used, defined by five pairs of nodes. The interface characteristics were determined based on the adjacent soil by setting the strength reduction factor (R_inter) equal to 0.8. At the first computational stage, initial vertical horizontal effective stresses were calculated. On a second stage, elastic-plastic deformation analysis was carried out. At this point, the backfill was simulated to be constructed per stages of one meter, starting from the wall foundation level. This is a realistic approach as it allows to account for the generated deformations at each individual stage of construction.

In order to estimate the seismic response of the three wall models, eight real acceleration time histories from five different seismic events were selected, as the strongest seismic motions of the Hellenic database recorded on rock conditions (Margaris et al., 2021 [25]) (

Table 5. Characteristics of input motions used for the dynamic analyses (Margaris et al., 2021 [25]).

Seismic Event	Μw	Repi (km)	VS30 (m/s)	Record	PGA (g)	Tp (s)	DS5–95 (s)
Athens 1999	6	18.8	1020	ATH4_L	0.1186	0.10	6.51
Athens 1999	4.6	21.23	1020	ATH4_L2	0.0275	0.16	2.60
Cephalonia 2014	6.1	37.27	1183	VSK1_E	0.0783	0.14	4.22
				VSK1_N	0.0968	0.12	4.78
Evia 2014	5.3	22.4	872	VLSA_E	0.196	0.06	5.70
				VLSA_N	0.257	0.06	6.22
Crete 2015	6.1	61.12	877	ZKRA_E	0.0494	0.12	15.08
				ZKRA_N	0.056	0.34	15.73

). The normalized acceleration response spectra of the selected records, as well as the mean normalized spectrum ± one standard deviation, are presented in Figure 4. Aiming to correlate the induced displacements with the increasing seismic intensity under free field conditions, the input time histories were scaled to PGA values equal to 0.1 g, 0.2 g, 0.3 g, 0.4 g and 0.5 g. The upper scaling limit of PGA = 0.5 g has been chosen as a probable seismic intensity level in Greece, and precisely in the Region of East Macedonia and Thrace, according to seismic hazard maps proposed in the work of Sotiriadis et al. (2023) [26]. Therefore, 600 dynamic analyses (3 retaining wall models × 40 input motions × 5 initial conditions) and two additional sets of 40 2D soil response analyses were carried out to estimate the values of seismic intensity measures under free field conditions, considering and ignoring the presence of water table, respectively. It should be noted that the aim of this work is the proposal of fragility curves that have a broad implementation area and can be subsequently utilized in conjunction with seismic hazard analysis results concerning specific sites, in order to estimate seismic risk of the corresponding areas. Therefore, record selection was not performed based on spectrum compatibility, as this would have led to limitations regarding the site applicability of the proposed fragility curves. Instead, real acceleration time histories recorded on rock conditions were selected and scaled, following the incremental dynamic analysis method. A similar approach regarding ground motion selection has been followed in the literature by Che et al. (2020) [27], Zamiran and Osouli (2018) [12], and Lesgidis et al. (2017) [28].

To calibrate the soil model and verify the FE seismic analyses results, the response of a one-dimensional (1D) soil column created in PLAXIS was compared with the results obtained from the 1D site response analysis software STRATA [29]. By adopting the soil properties of the 2D model, a 1D soil column was simulated in PLAXIS after applying the boundary condition of tied degrees of freedom at the lateral boundaries. Thus, the nodes at the lateral boundaries of the model would undergo the same displacements and the simulation of 1D wave propagation could be achieved. Moreover, the fundamental frequencies of the three soil-retaining wall models of 3 m, 6 m, and 9 m height were estimated by generating the ratios of the power spectra of the wall stem over the power spectra at the bedrock.

In this study, the permanent horizontal displacement of the wall base (U_x), and the permanent vertical ground displacement of the backfill material (U_y), were selected as damage indices (DI) to be compared with predefined thresholds that express different damage levels. For the horizontal displacement, three threshold values of U_x, expressed as a function of the wall height (0.02 × H m, 0.05 × H m, 0.10 × H m) were assumed [13,14]. Regarding the vertical ground displacement of the backfill, the damage criteria proposed by the research program SYNER-G [16] were adopted, considering three limit values of U_y (0.05 m, 0.15 m, 0.40 m). Based on the different limit states three damage levels were determined, namely minor, moderate, and extensive damage state, that were qualitatively correlated with the road serviceability level according to Argyroudis et al. (2013) [15], as presented in

Table 6. Definition of damage states for retaining walls.

Damage State	Ux (m)	Uy (m)	Serviceability of the Road
DS1: Minor	0.02 × H	0.05	Open road with reduced speeds or partially closed during repair works
DS2: Moderate	0.05 × H	0.15	Closed or partially closed road during repair works
DS3: Extensive	0.10 × H	0.40	Closed road during repair works

.

The seismic intensity was expressed by means of intensity measures (IM), which characterize the ground motion while showing a good correlation with the induced damage. Utilizing the outcomes of the dynamic analyses, the ln(IM)−ln(DI) dispersion plots were generated. Subsequently, the linear regression expressions were calculated in order to simulate the relationship between the different seismic intensity measures (independent variables) and the damage indices in terms of horizontal displacement of the base of the wall and vertical displacement of the backfill material (dependent variables). In addition, the approach of combined criteria was investigated, to estimate fragility due to the most critical damage mechanism. To derive the envelope fragility curves, new dispersion plots were generated utilizing both the outcomes of U_x and U_y as dependent variables. Specifically, the computed values of permanent horizontal displacements of the wall base were plotted together with the computed values of permanent vertical displacements of the backfill, versus the corresponding intensity measures. This resulted in the generation of linear regression expressions, accounting for both damage indices. The intensity measures selected for evaluation were the free field peak ground acceleration (PGA_FF), as the most commonly used measure for describing ground motion, the free field peak ground velocity (PGV_FF), and two parameters related to the energy content of the ground motion, namely the free field Arias intensity (IA_FF), and the free field cumulative absolute velocity (CAV_FF). The evaluation of the intensity measures was performed based on the criteria of practicality, efficiency, and proficiency [13]. Efficiency expresses the degree of dispersion of the structure’s response with respect to seismic intensity, and it is evaluated by the lognormal standard deviation of seismic demand (β_D) or by the coefficient of determination (R²). Practicality expresses whether the intensity measure is directly correlated with the damage, evaluated by the slope of the linear relationship that describes the evolution of damage with the increasing seismic intensity. Proficiency is defined as the ratio of efficiency to practicality, expressing the composite result of the two criteria. A lower value of the ratio of lognormal standard deviation to the slope of the linear relationship indicates a more proficient intensity measure.

After defining the damage states and selecting the optimal seismic intensity measures, we proceeded to the calculation of probabilities of exceeding the three predefined damage states. This was achieved by applying the lognormal cumulative distribution function:

$${\mathrm{P}}_{\mathrm{f}}\left(\mathrm{ds}\ge {\mathrm{ds}}_{\mathrm{i}}|\mathrm{IM}\right)=\mathsf{\Phi }\left[\frac{1}{{\mathsf{\beta }}_{\mathrm{tot}}}\ast \ln \left(\frac{\mathrm{IM}}{{\mathrm{IM}}_{\mathrm{mi}}}\right)\right]$$

(4)

where ${\mathrm{P}}_{\mathrm{f}}$ is the probability of exceeding a damage state $\mathrm{ds}$ for a given level of seismic intensity, which is determined by the seismic intensity measure $\mathrm{IM}$, $\mathsf{\Phi }$ is the standard cumulative probability function, ${\mathrm{IM}}_{\mathrm{mi}}$ is the median threshold value of the intensity measure required to cause the ith damage state, and ${\mathsf{\beta }}_{\mathrm{tot}}$ is the total lognormal standard deviation. For the estimation of the total lognormal standard deviation, three different sources of uncertainty were taken into account, associated with the definition of damage states $({\mathsf{\beta }}_{\mathrm{ds}})$, the capacity of the element $({\mathsf{\beta }}_{\mathrm{C}}$), and the seismic input motion $({\mathsf{\beta }}_{\mathrm{D}}$), as described by Equation (5):

$${\mathsf{\beta }}_{\mathrm{tot}}=\sqrt{{\mathsf{\beta }}_{\mathrm{ds}}^2+{\mathsf{\beta }}_{\mathrm{C}}^2+{\mathsf{\beta }}_{\mathrm{D}}^2}$$

(5)

Finally, fragility curves were proposed. Concerning the approach of combined criteria, the envelope fragility curves were generated as the composite result of two sets of fragility curves, i.e., a set adopting the limit state values of the U_x criterion, and a set adopting the limit state values of the U_y criterion. The comprehensive assessment of the results indicated the impact of initial conditions and wall dimensions on the seismic fragility of typical cantilever retaining walls, as well as the most critical between the failure mechanisms under investigation.

3. Results

In this section, the results derived from the aforementioned methodology are presented with reference to the validation of the soil FE model, the estimation of the models’ fundamental frequencies, the evaluation of the intensity measures, and the estimation of the fragility parameters required to generate the fragility curves.

3.1. Validation of the Finite Element Soil Model

In order to calibrate the soil model and verify the validity of the FE seismic analyses results, the response of a 1D soil column created in PLAXIS was compared with the results, obtained from the 1D site response analysis software STRATA [29]. Initially, a linear elastic analysis was performed in STRATA by applying the weaker seismic motion ATH4_L2 (PGA = 0.0275 g), since the soil response is almost linear in the range of very small strains. The results demonstrated the response at the soil surface, and the Fourier transfer function between the response at the soil surface and the input motion indicated the fundamental natural frequency of the soil column f₁ = 7.6 Hz and the second one f₂ = 19 Hz. Thereafter, adopting the linear elastic constitutive model in PLAXIS, a 1D soil response analysis was performed, using 5% Rayleigh damping at frequencies 7.6 and 19 Hz. Figure 5 and Figure 6 depict that both analyses yield the same results of acceleration time histories and response spectra at the ground surface, as well as the transfer functions which also show a good convergence.

Subsequently, a comparative evaluation of soil response was carried out adopting the equivalent linear approach in STRATA and the non-linear HSsmall constitutive model in PLAXIS. To account for the hysteretic behavior of soil during the seismic loading in the equivalent linear analyses, shear modulus reduction and damping curves of Darendeli and Stokoe (2001) [30] were utilized. Figure 7 shows the response spectra and acceleration time histories at the ground surface, as obtained from the equivalent linear and non-linear analyses, by applying the input motions ZKRA_N (PGA = 0.0571 g), and VLSA_E (PGA = 0.1998 g). Since the deviations of the equivalent linear and non-linear approach depend on the degree of non-linearity of the soil, it is concluded that the PLAXIS 1D-wave propagation results show satisfactory convergence with those of STRATA, and therefore the validity of the soil simulation using the HSsmall constitutive model is verified.

3.2. Estimation of Fundamental Frequencies

The root mean square of the power spectral density ratios, demonstrated the fundamental frequencies of the three retaining wall models, at 6 Hz for R.W.3, and around 4 Hz for the cases of R.W.6 and R.W.9. These results are depicted in Figure 8, for the input motion VSK_N scaled at a value of PGA equal to 0.1 g.

3.3. Evaluation of Intensity Measures

The parameters that characterize seismic intensity PGA_FF, PGV_FF, IA_FF, and CAV_FF were evaluated for their appropriateness in correlating with the induced damage, in terms of permanent horizontal displacement of the wall, permanent vertical displacement of the backfill, and the combination of the two damage criteria. The results of the investigation are expressed by rankings of practicality, efficiency, and proficiency, as presented in

Table 7. Evaluation of intensity measures by practicality, efficiency, and proficiency criteria, in terms of permanent horizontal displacement of the wall base U_x, permanent vertical displacement of backfill U_y, and combining U_x + U_y, for the case of F_s = 1.5.

Fs = 1.5		Ux			Uy			Ux + Uy
		Practicality	Efficiency	Proficiency	Practicality	Efficiency	Proficiency	Practicality	Efficiency	Proficiency
R.W.3	PGAFF	2.7299	0.4216	0.1544	2.3251	0.3184	0.1369	2.5225	0.3801	0.1507
	PGVFF	1.7062	0.4497	0.2636	1.4226	0.3582	0.2518	1.5644	0.3955	0.2528
	IAFF	1.0495	0.2458	0.2342	0.8962	0.1573	0.1755	0.9728	0.2241	0.2204
	CAVFF	1.6296	0.3439	0.2111	1.3567	0.3052	0.2249	1.4931	0.3356	0.2248
R.W.6	PGAFF	2.7404	0.4596	0.1677	2.2796	0.3857	0.1692	2.5100	0.4344	0.1731
	PGVFF	1.8473	0.3206	0.1735	1.5142	0.3061	0.2021	1.6808	0.3298	0.1962
	IAFF	1.0763	0.2364	0.2196	0.9036	0.1794	0.1985	0.9899	0.2259	0.2282
	CAVFF	1.6679	0.3474	0.2083	1.3854	0.2969	0.2143	1.5267	0.3374	0.2210
R.W.9	PGAFF	2.5576	0.4612	0.1803	2.0728	0.4049	0.1954	2.3157	0.4537	0.1959
	PGVFF	1.7610	0.2841	0.1613	1.4445	0.2506	0.1735	1.6027	0.2994	0.1868
	IAFF	1.0140	0.2453	0.2419	0.8356	0.2078	0.2486	0.9248	0.2635	0.2849
	CAVFF	1.5770	0.3338	0.2117	1.2956	0.2857	0.2205	1.4363	0.3377	0.2351

and

Table 8. Evaluation of intensity measures by practicality, efficiency, and proficiency criteria, in terms of permanent horizontal displacement of the wall base U_x, permanent vertical displacement of backfill U_y, and combining U_x + U_y, for the case of F_s = 1.1.

Fs = 1.1		Ux			Uy			Ux + Uy
		Practicality	Efficiency	Proficiency	Practicality	Efficiency	Proficiency	Practicality	Efficiency	Proficiency
R.W.3	PGAFF	2.5073	0.3826	0.1526	2.267	0.3614	0.1594	2.3872	0.4027	0.1687
	PGVFF	1.6840	0.3504	0.2081	1.5186	0.3084	0.2028	1.6013	0.3294	0.2057
	IAFF	0.9773	0.17259	0.1766	0.8881	0.1492	0.1681	0.9327	0.1888	0.2024
	CAVFF	1.5189	0.2072	0.1364	1.3747	0.2116	0.1539	1.4468	0.2603	0.1799
R.W.6	PGAFF	2.6330	0.4273	0.1623	2.3927	0.3909	0.1634	2.5129	0.4321	0.1720
	PGVFF	1.8279	0.3177	0.1738	1.6521	0.3053	0.1848	1.7400	0.3239	0.1862
	IAFF	1.0284	0.1708	0.1661	0.9329	0.1700	0.1822	0.9807	0.2202	0.2245
	CAVFF	1.5909	0.2704	0.1700	1.4388	0.2637	0.1833	1.5149	0.3014	0.1989
R.W.9	PGAFF	2.6371	0.4453	0.1689	2.5499	0.4140	0.1624	2.5350	0.4434	0.1749
	PGVFF	1.8492	0.3144	0.1700	1.7617	0.3113	0.1767	1.7719	0.3316	0.1872
	IAFF	1.0378	0.1733	0.1670	0.9917	0.1885	0.1901	0.9966	0.2075	0.2082
	CAVFF	1.6092	0.2658	0.1652	1.5310	0.2669	0.1743	1.5429	0.2882	0.1868

, for the two extreme cases of initial conditions, i.e., F_s = 1.5, and F_s = 1.1. Figure 9 also depicts the results plotted per criterion for each case of initial condition. It is observed that PGA_FF is the most practical and proficient intensity measure, while in terms of efficiency it generally shows values of lognormal standard deviation < 0.46. In the case of IA_FF, despite the maximum values of efficiency, the parameter shows the lowest degree of practicality and proficiency. CAV_FF performs well for all three criteria, while for the PGV_FF parameter, a generally satisfactory degree of efficiency, practicality, and proficiency is obtained.

3.4. Fragility Function Parameters

The implementation of the lognormal cumulative distribution function requires the estimation of the median threshold values ${\mathrm{IM}}_{\mathrm{mi}}$ for the three damage states, as well as the total lognormal standard deviation β_tot. The median threshold values ${\mathrm{IM}}_{\mathrm{mi}}$, were estimated for each damage state by the plots ln(IM) − ln(DI) that describe the evolution of damage with the increasing ground motion intensity. In

Table 9. Fragility function parameters for the cases of F_s = 1.5, 1.3, and 1.1, for the IM PGA_FF.

PGAFF		Fs = 1.5			Fs = 1.3			Fs = 1.1
		R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9
Ux	DS1 = 2% H	0.6157	0.7606	0.8814	0.4693	0.5701	0.6610	0.3541	0.4428	0.5181
	DS2 = 5% H	0.8623	1.0883	1.2297	0.6496	0.7982	0.9526	0.5104	0.6272	0.7334
	DS3 = 10% H	1.1126	1.4271	1.5821	0.8307	1.0296	1.2561	0.6729	0.8160	0.9539
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.654	0.6791	0.6802	0.609	0.6439	0.6497	0.6296	0.6577	0.6695
Uy	DS1 = 0.05 m	0.5921	0.5184	0.4895	0.4395	0.4047	0.3767	0.3594	0.3390	0.3342
	DS2 = 0.15 m	0.9496	0.8393	0.8316	0.7115	0.6432	0.6243	0.5835	0.5366	0.5250
	DS3 = 0.40 m	1.4480	1.2906	1.3348	1.0939	0.9727	0.9803	0.8994	0.8085	0.7857
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.5928	0.6315	0.6434	0.6127	0.6198	0.6265	0.6169	0.6346	0.6491
Ux + Uy	DS1 = 0.05 m	0.5832	0.5305	0.5168	0.4397	0.4093	0.3872	0.3433	0.3276	0.3261
	DS2 = 0.15 m	0.9015	0.8218	0.8307	0.6766	0.6300	0.6191	0.5439	0.5073	0.503
	DS3 = 0.40 m (0.3 m) *	1.1866	1.2147	1.2689	0.888	0.9260	0.9414	0.7271	0.7495	0.7406
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6281	0.6623	0.6752	0.6184	0.6348	0.6409	0.642	0.6608	0.6683

* Limit value of DS3 for the case of R.W.3.

, the values of ${\mathrm{IM}}_{\mathrm{mi}}$ and β_tot are indicatively presented for the three cases of F_s = 1.5, 1.3, and 1.1, concerning the intensity measure PGA_FF and the three damage indices, while

Table A1. Fragility function parameters for the case of F_s = 1.5.

Fs = 1.5		PGAFF (g)			PGVFF (cm/s)			CAVFF (cm/s)
		R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9
Ux	DS1 = 2% H	0.6157	0.7606	0.8814	32.3536	43.1532	53.416	1518.491	2125.69	2707.35
	DS2 = 5% H	0.8623	1.0883	1.2297	55.3547	72.6084	86.599	2664.446	3800.48	4654.50
	DS3 = 10% H	1.1126	1.4271	1.5821	83.097	107.6291	124.811	4076.93	5898.29	7012.81
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.654	0.6791	0.6802	0.6725	0.5939	0.5751	0.6069	0.6088	0.6012
Uy	DS1 = 0.05 m	0.5921	0.5184	0.4895	30.5131	24.6776	22.7852	1428.17	1140.607	1042.22
	DS2 = 0.15 m	0.9496	0.8393	0.8316	66.0500	50.9801	48.7478	3209.64	2520.741	2433.45
	DS3 = 0.40 m	1.4480	1.2906	1.3348	131.6141	97.4354	96.1269	6613.49	5116.797	5188.05
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.5928	0.6315	0.6434	0.6321	0.5862	0.5593	0.5858	0.5815	0.5759
Ux + Uy	DS1 = 0.05 m	0.5832	0.5305	0.5168	29.7201	25.5531	24.6296	1389.64	1184.26	1136.40
	DS2 = 0.15 m	0.9015	0.8218	0.8307	59.9842	49.1258	48.8827	2900.37	2432.02	2441.86
	DS3 = 0.40 m (0.3 m) *	1.1866	1.2147	1.2689	93.4247	88.0525	90.1441	4613.89	4623.59	4833.86
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6281	0.6623	0.6752	0.6573	0.5989	0.5828	0.6022	0.6032	0.6034

* Limit value of DS3 for the case of R.W.3.

,

Table A2. Fragility function parameters for the case of F_s = 1.4.

Fs = 1.4		PGAFF (g)			PGVFF (cm/s)			CAVFF (cm/s)
		R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9
Ux	DS1 = 2% H	0.5026	0.6129	0.7211	25.1296	33.6937	42.4850	1185.82	1665.457	2184.628
	DS2 = 5% H	0.6952	0.8568	1.0234	41.9035	55.4081	70.8132	2094.93	2958.293	3954.065
	DS3 = 10% H	0.8886	1.1041	1.3338	61.6929	80.7217	104.222	3222.06	4568.636	6193.854
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6044	0.6293	0.6343	0.66	0.6022	0.586	0.5936	0.5864	0.5786
Uy	DS1 = 0.05 m	0.4795	0.4278	0.4049	23.3239	19.7387	18.3320	1092.284	899.90	825.12
	DS2 = 0.15 m	0.7568	0.6830	0.6631	48.1952	39.6561	37.5181	2439.71	1999.71	1886.05
	DS3 = 0.40 m	1.1372	1.0370	1.0300	92.1325	73.9295	71.1092	4999.468	4079.05	3945.43
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.5753	0.6038	0.61	0.6274	0.5881	0.5686	0.5721	0.5657	0.5564
Ux + Uy	DS1 = 0.05 m	0.4750	0.4369	0.4245	22.9807	20.3812	19.6257	1074.08	932.69	891.44
	DS2 = 0.15 m	0.7229	0.6732	0.6682	44.6704	38.7713	37.9846	2246.74	1954.08	1915.08
	DS3 = 0.40 m (0.3 m) *	0.9422	0.9903	1.0018	67.9425	68.8411	68.4928	3579.13	3781.93	3790.34
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.595	0.6206	0.6302	0.6486	0.60	0.59	0.5883	0.5804	0.576

* Limit value of DS3 for the case of R.W.3.

,

Table A3. Fragility function parameters for the case of F_s = 1.3.

Fs = 1.3		PGAFF (g)			PGVFF (cm/s)			CAVFF (cm/s)
		R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9
Ux	DS1 = 2% H	0.4693	0.5701	0.6610	22.5841	30.1715	37.4722	1053.48	1467.58	1864.34
	DS2 = 5% H	0.6496	0.7982	0.9526	37.4521	49.3309	63.9797	1850.35	2593.36	3413.58
	DS3 = 10% H	0.8307	1.0296	1.2561	54.9103	71.5552	95.8949	2833.42	3989.40	5394.14
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.609	0.6439	0.6497	0.6479	0.5945	0.6107	0.5794	0.5823	0.5681
Uy	DS1 = 0.05 m	0.4395	0.4047	0.3767	20.5171	18.2550	16.4840	942.97	822.66	738.69
	DS2 = 0.15 m	0.7115	0.6432	0.6243	42.4034	36.0423	34.4060	2143.86	1794.78	1689.12
	DS3 = 0.40 m	1.0939	0.9727	0.9803	81.0744	66.1562	66.3659	4463.28	3601.48	3534.67
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6127	0.6198	0.6265	0.6087	0.5869	0.5912	0.5763	0.5657	0.5515
Ux + Uy	DS1 = 0.05 m	0.4397	0.4093	0.3872	20.4644	18.5798	17.1428	942.10	838.28	770.95
	DS2 = 0.15 m	0.6766	0.6300	0.6191	39.6285	34.9430	34.0207	1976.96	1736.26	1669.75
	DS3 = 0.40 m (0.3 m) *	0.888	0.9260	0.9414	60.13	61.4141	62.7325	3155.68	3326.13	3328.83
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6184	0.6348	0.6409	0.634	0.5941	0.6039	0.5851	0.5772	0.5629

* Limit value of DS3 for the case of R.W.3.

,

Table A4. Fragility function parameters for the case of F_s = 1.2.

Fs = 1.2		PGAFF (g)			PGVFF (cm/s)			CAVFF (cm/s)
		R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9
Ux	DS1 = 2% H	0.4193	0.5087	0.5869	19.0396	25.5276	31.3448	870.61	1209.09	1532.96
	DS2 = 5% H	0.5872	0.7128	0.8273	31.8545	41.5207	51.3266	1544.27	2119.65	2712.24
	DS3 = 10% H	0.7576	0.9200	1.0727	47.0165	59.9891	74.5351	2382.36	3241.12	4176.20
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6171	0.6607	0.6591	0.6297	0.5931	0.5877	0.5615	0.575	0.574
Uy	DS1 = 0.05 m	0.4138	0.3786	0.3637	18.6337	16.6714	15.7620	852.11	741.20	693.78
	DS2 = 0.15 m	0.6535	0.5962	0.5776	37.5950	32.0919	30.6399	1853.20	1571.14	1491.41
	DS3 = 0.40 m	0.9826	0.8943	0.8730	70.3526	57.5871	55.4641	3708.28	3072.67	2953.49
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6008	0.6441	0.6382	0.6169	0.5896	0.581	0.5561	0.5673	0.5638
Ux + Uy	DS1 = 0.05 m	0.4021	0.3732	0.3634	17.8484	16.3341	15.7456	811.00	723.19	692.24
	DS2 = 0.15 m	0.6174	0.5725	0.5619	34.4210	30.2726	29.4416	1681.70	1471.07	1425.63
	DS3 = 0.40 m (0.3 m) *	0.8092	0.8389	0.8290	52.0929	52.5143	51.4782	2664.3	2772.99	2717.18
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6198	0.6588	0.8290	0.634	0.5986	0.5886	0.5707	0.5785	0.5732

* Limit value of DS3 for the case of R.W.3.

and

Table A5. Fragility function parameters for the case of F_s = 1.1.

Fs = 1.1		PGAFF (g)			PGVFF (cm/s)			CAVFF (cm/s)
		R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9	R.W.3	R.W.6	R.W.9
Ux	DS1 = 2% H	0.3541	0.4428	0.5181	17.5285	20.9055	26.1938	664.5898	961.29	1244.70
	DS2 = 5% H	0.5104	0.6272	0.7334	30.2031	34.5115	42.9926	1214.915	1709.97	2199.65
	DS3 = 10% H	0.6729	0.8160	0.9539	45.5838	50.4254	62.5436	1917.5	2643.68	3383.94
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6296	0.6577	0.6695	0.6105	0.5924	0.5906	0.5412	0.5684	0.5662
Uy	DS1 = 0.05 m	0.3594	0.3390	0.3342	15.2181	14.1858	13.9660	681.658	615.79	604.98
	DS2 = 0.15 m	0.5835	0.5366	0.5250	31.3723	27.5836	26.7082	1515.791	1321.42	1273.16
	DS3 = 0.40 m	0.8994	0.8085	0.7857	59.8476	49.9442	47.6461	3093.871	2612.75	2473.90
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.6169	0.6346	0.6491	0.6037	0.5858	0.589	0.5429	0.5653	0.5668
Ux + Uy	DS1 = 0.05 m	0.3433	0.3276	0.3261	14.2198	13.5224	13.5081	631.57	582.65	581.79
	DS2 = 0.15 m	0.5439	0.5073	0.503	28.2392	25.4251	25.1111	1349.58	1203.27	1185.77
	DS3 = 0.40 m (0.3m) *	0.7271	0.7495	0.7406	43.5355	44.6757	43.6785	2179.04	2299.05	2239.16
	${\mathsf{\beta }}_{\mathrm{tot}}$	0.642	0.6608	0.6683	0.6265	0.6054	0.60	0.5637	0.5838	0.5771

* Limit value of DS3 for the case of R.W.3.

in Appendix A include the fragility function parameters regarding all five initial conditions, in terms of PGA_FF, PGV_FF, and CAV_FF, and U_x, U_y, and U_x + U_y. Figure 10 shows the case of R.W.6 for F_s = 1.2, considering the permanent vertical displacement of the backfill as damage index and the free field PGV as intensity measure. It should be noted that in cases where the calculated displacements were not sufficient to cause extensive damage, the linear relationship was extrapolated [13]. The uncertainties associated with the damage states can be addressed utilizing the maximum likelihood method, as presented in the work of Karakas et al. (2022) [31]. In this work, the values of the uncertainties ${\mathsf{\beta }}_{\mathrm{ds}}$ and ${\mathsf{\beta }}_{\mathrm{C}}$ were set equal to 0.4 and 0.3 according to HAZUS (2004) [8]. The coefficient of uncertainty that expresses the deviation of the structure’s response at each intensity level due to the different input motions (${\mathsf{\beta }}_{\mathrm{D}})$, was calculated as the root mean square error of the differences between the damage index values predicted by the linear regression line, and the observed values that resulted from the dynamic analyses.

Regarding the investigation of combined criteria, the comparative assessment of fragility curves generated by adopting (a) the limit values of U_x criterion, and (b) the limit values of U_y criterion indicated that higher probabilities of damage arise from the second case, with the exception of R.W.3 for DS3. This is a reasonable conclusion due to the fact that the vertical displacement criterion requires generally lower limit state values than the horizontal displacement criterion. However, this statement does not apply to the case of R.W.3 for DS3, where the U_x threshold required to cause extensive damage (0.30 m) is lower than the corresponding U_y threshold (0.40 m). Hence, accounting for this exclusion, it is concluded that the envelope fragility curves are derived by adopting the limit values of the U_y criterion, indicating that the permanent vertical displacement of the backfill is a more critical mechanism for seismic induced damage on cantilever retaining walls.

3.5. Fragility Curves

Based on the implementation of the computed fragility curves parameters, IM_mi, concerning the three damage states, and β_tot (Table 9, and Table A1, Table A2, Table A3, Table A4 and Table A5 in Appendix A) to the lognormal cumulative distribution function (Equation (4)), the probabilities of exceeding each damage state for a given level of seismic intensity were calculated. Hence, the corresponding fragility curves for the three models R.W.3, R.W.6, and R.W.9 were proposed, referring to the five cases of initial conditions defined by the different values of global F_s. Figure 11 illustrates the fragility curves for F_s = 1.5, 1.3, and 1.1, in terms of PGA_FF − DIs, while the fragility curves regarding all five initial conditions, in terms of PGA_FF, PGV_FF, and CAV_FF, and U_x, U_y, and U_x + U_y, are illustrated on Figure A1, Figure A2, Figure A3, Figure A4 and Figure A5, in Appendix B.

The proposed fragility curves in terms of permanent vertical displacement of the backfill (U_y) denote that retaining walls R.W.6 and R.W.9 show approximately the same degree of fragility at each damage state, with divergence in damage probability lower than 5%. On the contrary, reduced fragility levels are presented for R.W.3, showing about 8–20% lower probabilities of damage than R.W.6. This difference between R.W.3 and the higher walls seems to be decaying with the decrease of F_s; while shorter deviation is generally observed in the cases of minor damage (DS1). The fragility curves proposed by adopting the combined criteria (U_x + U_y) denote overall similar remarks with the ones proposed by adopting U_y criterion, with the exception of extensive damage state (DS3) where all three walls show almost the same degree of fragility, regardless of their height. On the other hand, fragility curves calculated considering the permanent horizontal displacement of the base (U_x) as a damage index illustrate a different pattern on damage probability regarding the three models R.W.3, R.W.6, and R.W.9. Particularly, higher probabilities of damage are presented for walls of lower height. This pattern is retrieved from the definition of limit states as a percentage of the wall height. Therefore, R.W.3 requires lower values of horizontal displacement to exceed each damage state (0.06 m, 0.15 m, 0.30 m), as opposed to the corresponding values for R.W.6 (0.12 m, 0.30 m, 0.60 m), and R.W.9 (0.18 m, 0.45 m, 0.90 m). Considering that the horizontal displacement of the base represents the failure mechanism of sliding, whereas the vertical displacement of the backfill might result from the composite effect of sliding and rotation, it can be concluded that U_y is a more critical damage index than U_x.

For the comparative assessment of fragility curves referring to different initial conditions, Figure A6, Figure A7, Figure A8, Figure A9, Figure A10, Figure A11, Figure A12, Figure A13 and Figure A14, in Appendix C are presented. Herein, Figure 12 indicatively illustrates the cases of PGA_FF − U_x + U_y. Each figure corresponds to a different case of IM−DI, while each plot refers to an individual soil-wall model at a certain damage state. As expected, lower values of F_s result to higher values of damage probability. In cases where the presence of water table is considered, i.e., F_s = 1.1, 1.2, 1.3, 1.4, there is a relatively small change in the degree of fragility between two sequencing initial conditions. On the contrary, fragility curves corresponding to dry initial conditions, i.e., F_s = 1.5, show considerably lower probabilities of damage than the following curves of F_s = 1.4. In the case of R.W.3, DS2 for PGA_FF = 0.5 g, the estimated probabilities of damage for F_s = 1.1, 1.2, 1.3, and 1.4 are equal to 45%, 37%, 32%, and 27%, respectively; while for F_s = 1.5 the probability of damage drops to 17%, demonstrating that there is no proportionality in the degradation of fragility curves.

4. Discussion

In this study, novel fragility curves for typical cantilever retaining walls have been developed, with the aim of investigating the impact of initial conditions on seismic fragility. Three retaining wall models of the same typology but with different heights were simulated in order to compute their seismic response, considering dry conditions and also the presence of water table at four different levels, leading to five discrete values of global safety factor ranging from 1.1 to 1.5. Damage probability was determined in terms of permanent horizontal displacement of wall base, permanent vertical displacement of backfill, and utilizing the outcomes of both horizontal and vertical displacements, with respect to increasing seismic intensity.

At an early stage of the investigation, the estimation of fundamental frequencies of the three soil-wall models resulted in the same value of f₁ = 4 Hz, for R.W.6 and R.W.9, indicating a potential similar dynamic response for the two models. This was later verified by the fragility curves proposed in terms of U_y and U_x + U_y, with the two walls showing approximately the same degree of fragility. Additionally, R.W.3 showed lower probabilities of damage, demonstrating that with a decrease in height, the induced seismic deformation reduces. On the other hand, a reversed pattern was observed in the case of fragility curves proposed in terms of U_x, resulting from the definition of limit state displacements as a percentage of wall height. Consequently, it is evident that the selection of the damage index, as well as the definition of damage states, plays a significant role on the formulation of fragility curves.

To calculate an envelope curve representative of the most unfavorable damage scenario at each intensity level, the approach of combining two damage criteria was implemented. The results denoted the vertical displacement of backfill as the critical damage index, leading to higher probabilities of damage compared to horizontal displacement, with the exception of only one case (R.W.3 − DS3). Taking into account that the horizontal displacement of the wall base is correlated with the failure mechanism of sliding, while the vertical displacement of backfill might result from the contribution of both sliding and overturning, the aforementioned remark can be regarded reasonable.

In order to validate the proposed fragility curves, a comparative assessment with existing curves from the literature was made. The research that has been carried out so far concerns different typologies of retaining walls and soil properties, leading to variable results. Hence, the proposed fragility curves by Seo et al. [14] regarding a 4m height cantilever retaining wall model with cohesionless soil corresponding to site category S2 of the South Korean classification system, were chosen for this comparison. Figure 13, depicts the fragility curves by Seo et al. (dashed), along with the newly proposed curves (continuous) that were retrieved after interpolation between the R.W.3 and R.W.6 curves for F_s = 1.5, in terms of PGA_FF − U_y, and PGA_FF − U_x. In the case of U_y, the two sets of fragility curves corresponding to minor damage show a maximum divergence in damage probability around 10%, for PGA = 0.3–0.4 g. In the case of U_x, the sets of curves corresponding to DS1 show a maximum divergence around 12%, for the larger PGA values. For both damage indices, the proposed curves of this work illustrate 10–13% higher damage probabilities for moderate PGA values in DS2, whereas the same divergence is observed in DS3 for larger PGA values. The lower damage probabilities concerning the curves from literature, might be attributed to the presence of a shear key at the base of the retaining wall model, resulting in a higher sliding resistance.

The comparative representation of fragility curves derived for different values of safety factor demonstrated the effect of initial conditions to the degree of seismic fragility. In cases of heavy rainfall or retaining walls with malfunctioning/insufficient drainage, the implementation of fragility curves corresponding to F_s = 1.1 would be recommended in order to cover all intermediate initial conditions. In cases of dry conditions, fragility curves corresponding to F_s = 1.5 shall be implemented. Future research could further examine in what way factors other than the level of water table will affect the degree of fragility, assuming the same initial conditions of F_s = 1.1–1.5. The proposed fragility curves constitute a reliable tool that can be utilized in conjunction with seismic hazard analysis results in seismic risk assessment. The results of this investigation were derived concerning cohesionless soil materials, and typical geometries of cantilever retaining walls. For retaining structures of different typology, relevant fragility curves should be considered or developed. Regarding pore water pressure, the assumption of hydrostatic conditions has been made. For a more realistic simulation of water conditions, groundwater flow could be taken into account in future research work.