In the past decades, seaports suffered large structural disorders, economic and social losses, mostly due to earthquakes. Among the damaged port structures are pile-supported wharves founded in sloping rockfill or through the rock dikes. These structures are commonly used worldwide, and have been the subject of several experimental and numerical studies. Yang^{[1]} has studied the seismic response of the Seventh Street Terminal Port of Oakland using FLAC software,^{[2]} and used some modifications to take into account the effect of sloping ground. Several dynamic centrifuge tests were conducted at the University of California, Davis (UCD).^{[3–10]} All centrifuge models represented case histories, consisting of rock dike configurations and wharf deck supported by piles. These tests were carried out in order to better understand the seismic performance of wharf piles. Hwang et al.^{[11]} performed shaking table model tests, to reproduce the dynamic behavior of a gravity quay wall and a pile-supported wharf damaged during the 1995 Kobe earthquake. The results of shaking table model tests were compared with the field measurements and with the results of the previous model testing. McCullough^{[3]} and Dickenson and McCullough^{[12]} presented centrifuge modeling results of pile supported wharves and pile in sloping rockfill. These results were used to validate two geotechnical seismic analysis methods: Newmark sliding block method,^{[13]} and numerical modeling using the finite difference code FLAC.^{[2]} The authors also performed a numerical analysis of the Seventh Street Terminal (Berths 35 through 38), at the Port of Oakland subjected to the 1989 Loma Prieta earthquake and compared the displacement results to in-situ observations.^{[14,15]} Takahashi and Takemura^{[16]} carried out a centrifuge model test of the Takahama Wharf damaged during the 1995 Hyogo-ken Nambu earthquake, and performed a parametric study to establish the effect of the thickness of the liquefiable sand layer on the results, and to study the piled deck-caisson dynamic interaction. Dickenson et al.^{[17]} provided a synthesis of modeling considerations and summary of the 2D FLAC numerical modeling results of the Berth 404 at Pier 400 of Port of Los Angeles, and considered pile-deck connections and pile-rockfill interaction. The results were validated with the case history data and the results from large-scale tests. Huertas and Romanel^{[18]} investigated the coupled hydro-mechanical behavior of soils and studied the influence of the liquefaction on the stability of pile supported wharf in rockfill dike, considering the nonlinear elasto-plastic constitutive model UBCSand. Dickenson et al.^{[19]} investigated the seismic performance of pile-supported wharves subjected to long-duration ground motions and considered the evaluation of kinematic loads on pile and soil-pile interaction in rockfill dikes. Su et al.^{[20]} conducted 3D finite element analysis of pile-supported wharf and considered the effect of the ground slope and pile connectivity along the wharf deck. Souri et al.^{[21]} conducted physical modeling of wharf piles in sloping rockfill subjected to liquefaction-induced lateral spreading, and studied kinematic demands and full superstructure inertia in order to estimate bending moments of piles at different depths. Vytiniotis et al.^{[22]} studied the effect of prefabricated vertical drains installed behind the crest of the slope to control the lateral spreading in the fill and to reduce seismic damage of the wharf structure.

Recently, a large volume of research works used seismic fragility curves as a useful tool to evaluate the seismic performance of pile-supported wharves.^{[23–30]}

Several port structures use batter piles to provide lateral resistance and to limit deflections from lateral loads, such as ship berthing and seismic loads. In the past, engineers contest the use of batter piles due to the poor performance under earthquake, caused by the lack of proper design of the batter piles head connections. Actually, by resolving this problem, batter piles have smaller deformations than vertical ones. Roth et al.^{[31]} studied the role of batter piles in the seismic performance of wharves. Schlechter et al.^{[32]} carried out a series of large-scale centrifuge models, with and without batter piles of typical pile-supported wharf configurations. Razavi et al.^{[33]} presented some batter pile damages during severe earthquakes, and performed seismic analysis of a typical pier with and without batter piles. Heidary-Torkamani et al.^{[28, 29]} developed the seismic fragility curves of wharves with batter piles, considering different engineering demand parameters (EDPs) for each damage state, and also performed a sensitivity analysis in order to evaluate the influence of geotechnical parameter uncertainties in the seismic performance of the wharves. Li et al.^{[34,35]}) studied dynamic response characteristics and failure process of vertical and batter pile-supported wharf structure, established on the same conditions, and subjected to two different earthquake levels.

The soil is a complex material and exhibits nonlinearly behavior when it is subjected to very low strain levels Atkinson and Sallfors,^{[36]} as for example, far-field deformations around deep excavation and tunnel, and soils subjected to seismic ground motion. Researchers developed different soil constitutive models with various degrees of complexity, and it is important to choose the most representative soil behavior model, as it has a significant impact on the results.

A great amount of investigation was conducted to determine the influence of soil constitutive models at small-strain levels for soils around deep excavation and tunnel; however, less researches were done in the case of soils under seismic loading. Truty^{[37]} studied the deformations near the excavations in Berlin sand, using Mohr-Coulomb (MC), Hardening Soil (HS) and Hardening Soil with Small-strain stiffness (HSS) models, and presented the soil deformation, the bending moment and the deflection of the wall. Hamrouni et al.^{[38]} highlighted the impact of soil constitutive models on the prediction of the displacements of the surrounding soil in the execution of a shallow tunnel. Benz et al.^{[39]} described soil stiffness at small-strains and presented the available small-strain stiffness models. The authors also performed four numerical applications on: triaxial test, spread footing problem and excavation problem 2D and 3D, using the finite element method with and without small-strain stiffness formulation. Obrzud^{[40]} presented Hardening Soil standard and Hardening Soil Small-strain models in finite element program, and exposed the comparison result with Mohr-Coulomb modeling of excavation in Berlin sand, tunneling and shallow footing in Texas sand. Besseling^{[41]} evaluated the sensitivity of the site response to variation of the soil constitutive models. Four soil constitutive models implemented in PLAXIS software were considered: the linear elastic (LE), the Mohr-Coulomb (MC), the Hardening Soil (HS) and the Hardening Soil with Small-strain stiffness (HSS) models. Based on the comparison of site responses, the HSS model was adopted for modeling of the soil deposit in the study of the dynamic jetty analysis. Vakili et al.^{[42]} presented the numerical approach for adaptive constitutive modeling of soils using the finite element method, and gave results using this approach in numerical soil element tests (Oedometer and Triaxial tests) and staged tunneling construction. Vakili et al.^{[43]} evaluated the results for the vertical displacements of four selected nodes in shallow and deep tunnel projects using three soil constitutive models. Op de Kelder^{[44]} analyzed the horizontal deformations on the sheet pile wall as a result of excavation of a building pit, using the Hardening Soil and the Hardening Soil Small strain stiffness models implemented in the FE code PLAXIS 2D. For the calculation of shear modulus at very small strains (_{0}) and the threshold shear strain (_{0,7}) of HSS model, the author employed the correlation of Alpan^{[45]} and the correlation of Benz and Vermeer.^{[46]} Thereafter, the obtained numerical results were compared to the inclinometer measurements. Alpan^{[45]} gave curves that describe the relationship between static _{stat}_{dyn}_{stat}_{dyn}_{0}). _{0} is the Young’s modulus at very small strains. According to Wichtmann and Triantafyllidis,^{[47]} it can be assumed that the correlation of Alpan estimated _{stat}_{50} (_{stat}_{50}). For the correlation of Benz and Vermeer, it is recommended to interpret _{stat}_{ur}_{stat}_{ur}_{50}). Szerző and Batali^{[48]} used different constitutive laws for the soil and interface elements to analyze the numerical settlement results of piled raft foundations supporting a circular gas tank, and compared these to field measurements and two simplified methods: the Poulos-Davis-Randolph method^{[49]} and the Piglet software.^{[50]} Akbari Hamed[51] compared the experimental results of Triaxial and Oedometer tests to the numerical results in two finite element software, using various soil constitutive models, and also performed a comparison between the measured data and the numerical analysis results of the asphaltic core rockfill dam (Dam-X, Quebec). Yeganeh and Fatahi^{[52]} examined the effects of soil plasticity type and soil-structure interaction in the seismic analyses of the typical high rise buildings, constituted by 20-story reinforced concrete moment-resisting building with a mat foundation, using the explicit finite difference program FLAC3D.

Considering the few researches to evaluate the impact of choosing soil constitutive model on the seismic results and, in addition to the complexity of seismic behavior of wharf structures with batter piles in sloping rockfill, this study aimed to capture the more realistic seismic behavior of pile-supported wharf with batter piles in cut-slope rock dike. For this, three soil constitutive models with increasing levels of complexity were considered using finite element method. The techniques of the finite element modeling will be given in the following sections.

A two-dimensional finite element modeling was generated using the software PLAXIS 2D.^{[53]} The details of the modeling adopted in the dynamic analysis are listed below:

The geometry model boundaries are taken relatively far away to avoid direct influence of the boundary conditions Hamrouni et al.^{[54]} Horizontal fixity has been introduced at the vertical boundaries, while the bottom boundary of the model has been subjected to total fixity (both horizontal and vertical fixities).

A prescribed displacement is selected at the bottom boundary of the model in order to simulate the dynamic loading. Standard earthquake boundaries are specified at the right and left side of the model geometry, to avoid the reflection of the seismic waves on the model boundaries and perturbations.^{[55]} The relaxation coefficients _{1} and _{2}were kept in default values 1 and 0.25, respectively.

The mesh coarseness was fine and refined around the dike with 15 node triangular elements. The influence of the mesh must be tested for each analysis. In Plaxis,^{[56]} the Average Element Size (AES) is calculated from the outer geometry dimensions (_{min}_{max}_{min}_{max}_{c}^{[57]} provided a guideline to use AES, which must be smaller than approximately one-tenth to one-eighth of the wavelength (

The damping matrix [_{R}_{R}^{[58,59]} as presented in equation (4):
^{[60]} The displacement and the velocity at the point in time _{N}_{N}^{[61]} that uses the Newmark scheme to solve the modified equations of motion shown in equation (9):
_{N}_{N}_{R}^{[62]}

The standard setting of the PLAXIS code for time integration is the damped Newmark scheme with _{N}_{N}

The time step used in a calculation is constant and it is presented in equation (13):
^{[63,64]}

To avoid unreliable numerical response in a finite element model, a careful evaluation of the critical time step _{crit}^{[58,62]}

The configuration of the wharf chosen for this study is the JCB01 model,^{[10]} selected from a series of centrifuge models carried out at the Center for Geotechnical Modeling of California University, Davis Campus (UC Davis).

This model consists of a wharf deck supported by seven rows of vertical piles included two pairs of batter piles, embedded in cut slope rock dike configuration. The dike is constructed by a thin layer of rockfill overlaid by a slope of loose sand (40% relative density). The side slope was 2 horizontal to 1 vertical (

In numerical analysis, material behavior of the soils and structural elements are represented by material constitutive models.

The sandy soils used in the centrifuge model were Nevada sand.^{[65]} The rock was imported from a quarry on Catalina Island, off the coast of southern California, and crushed to an appropriate size in order to represent the rock-pile interaction.

The geotechnical parameter values used in numerical modeling are presented below according to the three soil constitutive models applied in this study:

The most common material model used is the linear-elastic perfectly plastic model (Mohr-Coulomb model). It requires only five input parameters: Young’s modulus ^{[66]} The stress-strain relationship is considered as elastic-perfectly plastic (Fig. 2(1)).

The soil is composed of granular particles, and is generally modeled in numerical analysis as a continuum. The individual particle size of the soil is in the order of 30 to 40 times less than the pile diameter or the dimensions of the foundation. In the case of rockfill used in the marine structures, the median rockfill size was 2.5 to 3 times smaller than the pile diameter. Then, the rock particles act as individual particles and less as a continuum. Diaz et al.^{[67]} introduced a concept of pseudo cohesion in rockfill so as to modify the modeled continuum behavior to match the measured behavior. McCullough^{[3]} added 15 kPa pseudo-cohesion for the rockfill in order to consider the individual rock particle interaction with the pile and took into account the difference between upslope and downslope spring stiffness. This concept was also adopted by Martin,^{[68]} and McCullough and Dickenson.[69] Kawamata^{[70]} considered the lateral reaction of the rockfill as a combination of reactions due to both friction and interlocking between large rock particles.

The FE code PLAXIS can handle cohesionless soils with cohesion equal to zero, but some options do not perform well. The material models manual of PLAXIS^{[71]} advises to enter at least a small value of cohesion (use ^{[3]}

The input parameters of soil used for the Mohr-Coulomb model are summarized in Table 1, and for more details refer to Boland et al.,^{[10]} McCullough,[3] Heidary-Torkamani et al.^{[28,29]} and Amirabadi et al.^{[26]}

Soil parameters used for the MC model.

° | 37 | 33.2 | 45 | |

° | 21 | 7 | 0 | |

kN/m^{2} | 0.5 | 0.5 | 15 | |

- | 0.3 | 0.3 | 0.2 | |

kN/m^{2} | 6.032 × 10^{4} | 1.638 × 10^{4} | 2.232 × 10^{4} |

The dry mass density was 1633 kg/m^{3}, 1522 kg/m^{3} and 1611 kg/m^{3} respectively for dense sand, loose sand and rockfill. The material damping of soils was simulated with 5% of damping ratio (Rayleigh damping).

The permeability was 0.043 × 10^{–3} m/s, 0.066 × 10^{–3} m/s and 0.01 m/s respectively for dense sand, loose sand and rockfill. For water conditions, the ground was saturated. A water depth of 12.4 m, and volumetric weight of 10.3 kN/m^{3} (salted water) were defined.

The Hardening Soil model, developed by Schanz et al.,^{[72]} describes a hyperbolic stress – strain curve (Fig. 2(2)). It has the same failure criterion as the Mohr-Coulomb model and uses three types of soil stiffness for three loading conditions, namely triaxial loading, oedometer loading and unloading/reloading, as presented in equations (14), (15) and (16):
^{ref}^{ref}^{2}).
^{[73]} The amount of stress dependency is given by the power

The HS model was implemented into PLAXIS finite element code. In many practical cases, it is appropriate to determine the stiffness moduli by using the default setting of PLAXIS:
^{[71]}

The soil parameters required in the Hardening Soil model for the materials used in this study, in addition to the soil parameters of MC model, are presented in Table 2.

Soil parameters used for HS model

kN/m^{2} | 6.032 × 10^{4} | 1.638 × 10^{4} | 2.232 × 10^{4} | |

kN/m^{2} | 4.826 × 10^{4} | 1.310 × 10^{4} | 1.786 × 10^{4} | |

kN/m^{2} | 1.810 × 10^{4} | 4.914 × 10^{4} | 6.696 × 10^{4} | |

^{ref} | kN/m^{2} | 100 | 100 | 100 |

- | 0.5 | 0.5 | 0.5 | |

- | 0.3982 | 0.4524 | 0.2929 | |

_{f} | - | 0.9 | 0.9 | 0.9 |

The Hardening Soil model with small-strain stiffness is a modification of the Hardening Soil model. It takes account the degradation of shear stiffness with shear strain during loading (developed by Benz, 2007).^{[74]}

The model’s stiffness parameters in a triaxial test were presented in Fig. 2(3). The stress-strain curve for small-strains was described by a simple hyperbolic law, analogous to the hyperbolic law for larger strains.^{[75]}

The HSS model uses almost entirely the same parameters as HS model, except two additional parameters: The initial or very small-strain shear modulus _{0}, and the shear strain level _{0,7} at which the secant shear modulus _{S} degrades to about 72% of initial shear modulus _{0}.

Several works have demonstrated that when increasing strain amplitude, soil stiffness decays nonlinearly.^{[36]} The stiffness reduction curves exhibit a characteristic S-shape in logarithmic scale (Fig. 3), only at very small-strain levels (less than 10^{–6}) where the soil describes linear elastic behavior.

The stress dependency is given by the following equation (17).
^{–6}).

The very small-strain shear modulus _{0} of non-cohesive soils can be determined using different methods, with different efforts and accuracy.^{[47]} In the absence of experimental data for the determination of _{0} parameter, approximations from empirical equations can be appropriate using the relation given by Hardin and Black^{[78]} (equation (18). Similarly, in the lack of experimental data for the determination of the threshold shear strain _{0,7}, correlations are also available, for example, a correlation between the threshold shear strain _{0,7} and the Plasticity Index PI,^{[79]} or using the relationship of Hardin-Drnevich[80] as shown in equation (19):
_{0} is the earth pressure coefficient at rest, and _{1}^{‘} is the effective vertical stress (negative pressure).

In dynamic applications, the HSS model introduces hysteric material damping and shows a typical hysteretic behavior (Ishihara, 1996),^{[81]} in contrast to the Hardening Soil model.

The values of the two additional soil parameters
_{0,7} used in this analysis are presented in Table 3:

Additional input parameters of HSS model.

_{0.7} | 4.918 × 10^{–4} | 3.629 × 10^{–4} | 4.246 × 10^{–4} | |

kN/m^{2} | 1.051 × 10^{5} | 9.108 × 10^{4} | 1.038 × 10^{5} |

The constitutive parameters of the MC model have to be calibrated to advanced models, using two soil element tests available in the two-dimensional finite element code PLAXIS. These tests are Oedometer test (one-dimensional compression) and Triaxial test. The following conditions are adopted:

In the Oedometer test, the soil sample is subjected to four phases: two loading phases from 0 to 100 kPa and

200 kPa, followed by two unloading phases to 100 kPa and 50 kPa.

In Triaxial test, a confining pressure _{3} of 100 kPa is applied.

The calibration is limited to the variation of elasticity modulus (^{ref}

According to the comparison between curves for both Oedometer and Triaxial tests, and for all the different soil materials, the best consistency between MC model and advanced models is observed when the value of the Young’s modulus in the MC model is assumed to be equal to

The structural elements of the JCB01 centrifuge model were made from aluminum section. The piles were modeled using aluminum (6061-T6) tubing that has 0.636 m diameter, 0.0508 m wall thickness and 70 GPa modulus of elasticity. The wharf deck has 28.1 m transverse length, 12 m longitudinal length, 0.255 m thickness and 27.23 kN/m^{3} density. More properties were given in Boland et al.^{[10]} and McCullough.^{[3]}

The main problem when modeling with PLAXIS 2D is the transition from three to two dimensions. Then, the wharf deck supported by three rows of piles was modeled as plane strain. The row of piles was transformed into continuous wall elements defined per meter length, called plane strain piles (^{[82–84]}

The normal stiffness for plane strain pile (_{psp}_{psp}_{psp}_{p}_{p}_{p}_{p}_{p}_{r}_{p-row i}

In this work, the material behavior of the structural elements (pile and wharf deck) is considered to be linear elastic and modeled by using plate elements. Interaction between ground (soils, rockfill) and plane strain piles is described using interface elements.

The material of the interface element could be modeled by a new material, which has the same strength properties of the surrounding soil (or rockfill) multiplied by a factor, called strength reducing factor for interface _{int}_{int}_{i}_{i}_{soil}_{soil}

The input parameters used in the modeling of the structural elements are presented in Table 4. For more details, refer to studies by Boland et al.,^{[10]} McCullough,[3] Heidary-Torkamani et al., ^{[28,29]} and Amirabadi et al.[26]

Pile and deck properties.

kN/m | 162.75 × 10^{4} | 17.85 × 10^{6} | |

kN.m^{2}/m | 70.49 × 10^{3} | 96.724 × 10^{3} | |

kN/m/m | 0.633 | 6.9436 | |

- | 0.2 | 0.2 |

The numerical model is subjected to the 1989 Loma Prieta earthquake time history recorded at the Oakland Outer Harbor Station (Fig. 5):

In order to verify the numerical model using the finite element method of the PLAXIS 2D, the displacement of soils around the rock dike and the pore pressure are determined and compared to the numerical results predicted using the two-dimensional ^{[28]}

The results of the total displacement of soils around the cut-slope rock dike for the three soil constitutive models are presented in Fig. 6, 7 and 8, respectively, for MC, HS and HSS models.

The maximum total displacement of soil observed for all three soil used constitutive models is located in the upper part of the loose sand embankment. The value of this maximum total displacement is obtained between 0.14 m and 0.15 m with MC and HS models, and between 0.08 m and 0.09 m with HSS model.

Heidary-Torkamani et al.^{[28]} have obtained 0.08973 m for the maximum vector displacement (Fig. 9). This displacement is included in the interval of the displacement determined with HSS model. It can be concluded that the results obtained from a 2D plane strain model (as seen in the paragraph 4.3 above), by using the HSS model give adequate results for soil and wharf displacements.

Concerning the pore pressure contour results obtained in the numerical model with PLAXIS 2D, the contours have the same shape for all three soil used constitutive models. So, it was plotted here with just one pore pressure contours (Fig. 10). The pore pressure around the zone of embedded pile is approximately 0.2 MPa. This result coincides with the total pore pressure contours presented by Heidary-Torkamani et al.^{[28]} in Fig. 11.

It can be justified that when the soils are subjected to dynamic loading, the MC, HS and HSS models of the finite element code PLAXIS 2D and the MC model of the finite difference code FLAC 2D^{[28]} give the same pore pressures, due to the fact that, all these models do not take into account the accumulation of pore pressures or liquefaction.^{[71]}

Fig. 12 presents the total pile displacements under the 1989 Loma Prieta earthquake for the three soil constitutive models.

According to the pile deflection shown in Fig 12, the total displacements at the base of the piles are small and are at about the same for all the models. Thereafter, these total displacements gradually increase until maximum values are reached in the upper part of the piles, and become different for the three models.

This is due to the fact that at depth, the piles are embedded in good soil (dense sand layer), able to take up the efforts of wharf and undergo small deformations during the earthquake. While in the dike embankment, liquefaction occurs and induces large deformations especially at the embankment top.

Analyzing the maximum total pile displacements, it can be noted that for the first three piles of the wharf and also for the batter piles, the maximum total displacements are reached at the piles head. For the other piles, the maximum total displacements are particularly localized at a height in the range of 20.2 m to 22.2 m.

The largest of the maximum total displacements for the different piles of the wharf is observed for all models at the last pile (pile 7), caused by the significant deformation of loose sand at the top of the embankment, as shown in Fig. 6, 7 and 8. The values obtained are 0.139 m, 0.131 m and 0.077 m for the MC, HS and HSS models, respectively.

Comparing the results of this study with those of the numerical analysis using the finite difference method (Fig. 9, where the maximum pile displacement equals 0.08606 m), the MC and HS models give a maximum total pile displacement overestimated by 61.51% and 52.2% respectively, while for the HSS model, the maximum total displacement is underestimated by 10.52%.

It can be concluded that the results obtained with the HSS model are close to the results of numerical analysis conducted by Heidary-Torkamani et al.,^{[28]} and the MC model results are larger than the results of advanced models (HS and HSS models). This is due to the fact that the MC model uses a constant Young modulus for soil and only Rayleigh damping, in contrast to the advanced models, which use three different types of stiffness modulus (
^{[73]}

Table 5 shows the displacement time histories at the top and at the toe of the cut-slope rock dike.

Lateral and vertical displacement time histories at the top and at the toe of the rock dike.

The top of the cut-slope rock dike undergoes a significant lateral displacement caused by the deformations of loose sand at the embankment top. The maximum lateral displacement calculated by MC and HS models is 89.77% and 46.59%, respectively, higher than the result of HSS model.

At the end of seismic loading, the MC and HS models recorded a permanent lateral displacement of 85.5% and 65.21%, respectively, higher than that obtained by the HSS model.

When looking at the vertical displacement of the dike top, the MC and HS models provide a maximum settlement of 135.29% and 52.94%, respectively, larger than that predicted by the HSS model. For the permanent vertical displacement, the results obtained with MC and HS models increase relative to HSS model by 84.61% and 61.53%, respectively.

Concerning the toe of the cut-slope rock dike during the earthquake, the MC and HS models have a maximum lateral displacement of 51.51% and 39.39%, respectively, greater than in the HSS model. For the permanent lateral displacement, the MC and HS models recorded a result of 76% and 62%, respectively, higher than that estimated by the HSS model.

It can be noted that the rock dike toe has a tendency to heave. The maximum vertical displacement time history calculated for MC, HS and HSS models are 0.075 m, 0.031 m and 0.015 m, respectively, and occurred from the twentieth second during the earthquake loading. Beyond that time, the maximum vertical displacements stabilize. The MC and HS models give a maximum vertical displacement of 400% and 106.66%, respectively, higher than that observed by the HSS model.

Comparing the lateral displacements of the dike top with those of its toe, it can be concluded that the top moves more than the toe of the dike by 67%, 40.21% and 33.33% for the maximum lateral displacement, and by 45.45%, 40.74% and 38% for the permanent lateral displacement calculated with the MC, HS and HSS models, respectively. This difference is due to the significant deformation of the loose sand at the embankment top, as presented in Fig. 6, 7 and 8.

Concerning the predicted vertical displacement, the rock dike top tends to settle contrary to its toe, which heaves.

It can also be noted, that all the results obtained at the top and at the toe of the rock dike using the MC and HS models are relatively higher than the HSS model results.

The two-dimensional finite element method was used to study the sensitivity of the seismic response of pile-supported wharf to the variation of the soil constitutive models. Indeed, three constitutive models with different accuracies were adopted for this analysis, namely the MC model, the HS and HSS advanced models.

The JCB01 model selected from a series of centrifuge models carried out at the Center for Geotechnical Modeling of California University (UC Davis) was chosen for the study.^{[3,10,26,28,29]}

This model has a configuration of a typical wharf of Western United States ports. It consists of a deck supported by seven rows of vertical piles including two pairs of batter piles, embedded in cut slope rock dike overlaid an embankment of loose sand. The model was subjected to the 1989 Loma Prieta earthquake.

The wharf was modeled using the plane strain model. The concept of pseudo cohesion in rockfill was introduced. The analysis was carried out in time-domain and modification of the initial Newmark scheme defined with ^{[61]} was used. The constitutive parameters of soil were calibrated through numerical soil element tests.

In this work, the performance of numerical modeling in the estimation of wharf structural displacement and pore pressure were compared to the numerical results of Mohr Coulomb model in ^{[28]} The verification shows that the model based on higher stiffness under small-strain (HSS model) comes closest to the results of Heidary-Torkamani et al.^{[28]} contrary to the MC and HS models.

This is due to the fact that the MC model considers a constant Young’s modulus for soil and uses only Rayleigh damping, contrary to the advanced models (HS and HSS), which uses three different types of stiffness modulus (
_{0,7}, and incorporates the hysteretic damping of materials in addition to Rayleigh damping.

The use of approximations from empirical equations developed by Hardin and Black^{[78]} in order to determinate _{0} parameter for HSS model, gives appropriate result.

The displacement of piles in sloping rockfill and the displacement time histories of the rock dike at the top and at the toe were also highlighted. It can be observed that the earthquake induced significant lateral movement at the top of the cut-slope rock dike, due to deformations at the upper part of the embankment during the liquefaction of loose sand. This lateral ground displacement affects the behavior of the wharf and the displacement of piles, specially the piles closest to the top of the rock dike.

At the base of the piles, the total displacements are small, and then gradually increase until reaching maximum values at the top of the first three piles (seaward) as well as for the batter piles. But for the other piles, the maximum values of the total displacement are reached at the upper part of the piles when deformation of loose sand is important. It can be also seen that the batter piles have performed adequately.

According to the analysis results, the MC and HS models overestimate the maximum total displacement in piles by 61.51% and 52.2%, respectively, whereas for the HSS model, the maximum total displacement is underestimated by 10.52%, relative to the result of the numerical analysis established with the finite difference method.^{[28]}

The displacement time history at the dike top is greater than its displacement at the toe due to lateral movement at the embankment top. A difference of about 67%, 40.21% and 33.33% is recorded for the maximum lateral displacement for MC, HS and HSS models, respectively, and a difference of about 45.45%, 40.74% and 38% is estimated regarding its permanent lateral displacement for MC, HS and HSS models, respectively. Concerning predicted vertical displacement, the top of the rock dike has a tendency to settle contrary to its toe, which heaves.

The influence of the soil constitutive models used in this study is demonstrated on the results of piles and cut-slope rock dike displacements. It is concluded that the MC model is not appropriate for modeling seismic behavior of pile-supported wharf structures, and it is recommended to use an advanced model with nonlinear behavior under small-strain and hysteretic damping.

The response of the piles embedded into the rock dike is dictated by the resulting soil deformation mechanisms. The choice of the appropriate soil constitutive model has significant effects on the numerical results, and depends on many factors related to the type of analysis to perform, expected precision of predictions and available soil parameters.

The use of 2D finite element analysis with plane-strain and HSS model is an efficient tool for predicting the seismic behavior of pile-supported wharf. However, some suggestions and perspectives have emerged from this work. First, it would be interesting to use a model better suited to liquefaction analysis and allowing to determine the generation of excessive pore pressure during an earthquake (such as the UBCsand model), and to compare its results with those obtained in this paper.

Secondly, it can be noticed that the use of complex soil constitutive model (as advanced models) is computationally expensive, especially in dynamic calculation. And in order to remedy this, it is advisable to use adaptive constitutive modeling. That can considerably reduce the computational effort. This approach consists of using different soil constitutive models with different levels of complexity, but from the same hierarchical soil model family. Then, more sophisticated models are assigned to the regions that are directly influenced by the geotechnical problem. On the other hand, simpler models are assigned to other regions with low concentration of stress and strain. Lastly, it will be interesting to evaluate and to compare the results of this work with advanced 3D finite element analysis, focusing on the differences between the two methods.

#### Soil parameters used for HS model

kN/m^{2} | 6.032 × 10^{4} | 1.638 × 10^{4} | 2.232 × 10^{4} | |

kN/m^{2} | 4.826 × 10^{4} | 1.310 × 10^{4} | 1.786 × 10^{4} | |

kN/m^{2} | 1.810 × 10^{4} | 4.914 × 10^{4} | 6.696 × 10^{4} | |

^{ref} | kN/m^{2} | 100 | 100 | 100 |

- | 0.5 | 0.5 | 0.5 | |

- | 0.3982 | 0.4524 | 0.2929 | |

_{f} | - | 0.9 | 0.9 | 0.9 |

#### Additional input parameters of HSS model.

_{0.7} | 4.918 × 10^{–4} | 3.629 × 10^{–4} | 4.246 × 10^{–4} | |

kN/m^{2} | 1.051 × 10^{5} | 9.108 × 10^{4} | 1.038 × 10^{5} |

#### Lateral and vertical displacement time histories at the top and at the toe of the rock dike.

#### Soil parameters used for the MC model.

° | 37 | 33.2 | 45 | |

° | 21 | 7 | 0 | |

kN/m^{2} | 0.5 | 0.5 | 15 | |

- | 0.3 | 0.3 | 0.2 | |

kN/m^{2} | 6.032 × 10^{4} | 1.638 × 10^{4} | 2.232 × 10^{4} |

#### Pile and deck properties.

kN/m | 162.75 × 10^{4} | 17.85 × 10^{6} | |

kN.m^{2}/m | 70.49 × 10^{3} | 96.724 × 10^{3} | |

kN/m/m | 0.633 | 6.9436 | |

- | 0.2 | 0.2 |

^{th} ECSMFE, Florence 3, 1991, pp. 915–956.

Reviews on Finite Element Modeling Practices of Stone Columns for Soft Soil Stabilization Beneath an Embankment Dam Proposal of concept for structural modelling of hybrid beams Resonance of a structure with soil elastic waves released in non-linear hysteretic soil upon unloading Settlement Analysis of a Sandy Clay Soil Reinforced with Stone Columns The evolution of the shape of composite dowels Parametric study of the earth dam's behaviour subjected to earthquake Impact of longwall mining on slope stability – A case study