This thesis explores the application of a Bayesian approach to hyperparameter optimization in surrogate modeling for geotechnical engineering problems. Surrogate modeling, particularly employing Gaussian Processes and Kriging, has become an essential tool for accelerating complex numerical simulations in geotechnical engineering. Traditional Maximum Likelihood Estimation (MLE) approaches to hyperparameter optimization, although straightforward, often overlook the inherent uncertainties in model hyperparameters, potentially leading to sub-optimal prediction accuracy and poor generalization.

The research presented in this thesis investigates the feasibility and benefits of applying a Bayesian inference approach to the tuning of hyperparameters in surrogate models. This approach allows for a probabilistic treatment of hyperparameters, providing a comprehensive quantification of uncertainty. The Markov Chain Monte Carlo sampling method, specifically the No-U-Turn Sampler, was employed to sample from the posterior distributions of hyperparameters, addressing the challenges posed by their non-Gaussian nature and non-linear relationship with model outputs.

Three case studies of varying complexity from geotechnical engineering practice are examined to compare the Bayesian approach against traditional MLE in terms of hyperparameter determination, prediction accuracy, uncertainty quantification, and computational efficiency. The findings suggest that the Bayesian approach, while computationally more intensive, could potentially offer more accurate predictions in terms of the mean-squared error and provide a deeper understanding of uncertainty, which is crucial for risk-informed decision-making in geotechnical engineering.

The study concludes that Bayesian hyperparameter optimization in surrogate modelling holds significant potential for improving the robustness and reliability of predictions in geotechnical engineering, particularly in applications involving complex dependencies and where a thorough understanding of uncertainty is crucial. Further research is recommended to enhance the computational efficiency of the Bayesian method and to explore its integration with multi-point enrichment strategies for practical engineering applications.

The Random Finite Element Method (RFEM) is a robust stochastic method for slope reliability analysis that incorporates the spatial variability of soil properties. However, the extensive computational time associated with the direct Monte Carlo simulation limits its practical application. To overcome this problem, this study investigates the use of machine learning (ML) models as surrogate models for the RFEM in both 2D and 3D contexts. It investigates the performance of three ML models in predicting slope stability by means of the factor of safety (FoS) based on a generated random field of the undrained shear strength. Additionally, a data augmentation technique is employed to improve performance. The models' performance is assessed for various slope cases, characterised by varying spatial variability.
Two surrogate modeling approaches are employed: semi-surrogate modeling and full-surrogate modeling. In the semi-surrogate modeling approach, a small number of RFEM simulations are conducted for a specified case. The machine learning models are trained using the generated random fields as input data and the calculated factors of safety as output data. The mathematical models are then used to predict outcomes of FoS for a large number of random fields for the same specific slope case. In the full-surrogate modeling approach, many RFEM simulations are conducted for the training set, covering a range of spatial correlation lengths. Once trained, the full-surrogate models are ready for application to another different slope case without the need for any additional numerical simulation.
The results indicate that the prediction accuracy of the ML models typically decreases for slope cases with smaller scales of fluctuation. Nonetheless, the FoS predictions by the best-performing semi-surrogate model are highly consistent with the results from RFEM simulations for the whole range of considered slope cases. In terms of predicting the probability of failure for 2D-modeled slopes, the accuracy is high, with relative errors within 10% across the cases considered. This level of accuracy is achieved using no more than 13% of the total number of realisations needed for RFEM analysis. Consequently, the computational time for reliability analysis involving 4000 realisations reduces from 67 hours using the RFEM to between 4 and 8 hours using a semi-surrogate model, with the time increasing as the spatial correlation length decreases. Predicting the p_f for 3D slopes using a semi-surrogate model showed larger errors, indicating a need for improvement.
The full-surrogate models prove to be accurate for testing cases characterised by spatial correlation lengths within the training set's range. Notably, the best-performing full-surrogate model in 3D predicted the p_f within a relative error of 10% for two slope cases. This model performs a stochastic analysis of 4000 simulations within seconds, compared to 83 days of computational time required for RFEM reliability analysis.

Site characterization is indispensable in the design phase of geotechnical engineering projects. As a key factor in site characterization, the characterization of soil undrained shear strength (S_u) is always in the spotlight. Various methods, including laboratory and in-situ tests, have been developed to measure S_u. Nevertheless, these measurements are usually sparse at a specific site due to limited time and budget. To enhance S_u characterization, other relevant geotechnical investigation data (e.g., cone penetration test data), can be transformed into S_u through empirical correlations (referred to as transformation models) to provide more information on S_u. Considering this process introduces the transformation uncertainty and a developed transformation model may not be fully applicable to a local site, probabilistic transformation models (PTMs) have been developed to characterize soil parameters in a site-specific way and quantify the uncertainty to augment engineers’ judgement. However, few PTMs incorporate the spatial correlation of soil parameters, especially in the horizontal direction. This limitation hampers the ability to probabilistically characterize S_u in 2D/3D space, which is significant in practice. Moreover, estimating the horizontal spatial correlation from pure geotechnical data is challenging because they are typically sparse. In light of these circumstances, this thesis first proposes a PTM-based scheme to probabilistically characterize S_u in 2D. Then it is proposed to integrate geophysical data into the scheme. Compared to typical geotechnical investigations, geophysical surveys provide abundant 2D/3D measurement data, which are often correlated with geotechnical data. The fusion of these two data sources benefits characterizing geotechnical data including S_u. Particularly the horizontal spatial correlation of S_u 2D domain can be estimated from the abundant geophysical data. To be specific, a well-established PTM, MUSIC-X, by which measured S_u and other relevant soil parameters can be used to preliminarily characterize S_u, is first adopted. In this case, characterization specifically refers to simulating 1D vertical profiles of S_u. It is then combined with the intrinsic collocated co-kriging (ICCK) model, by which primary data (i.e., S_u) in 2D or theoretically 3D space can be estimated through linearly combining the preliminarily characterized S_u from MUSIC-X modelling and observed secondary data (i.e., geophysical data). The secondary parameter considered in this study is interval velocity (V_int). The scheme, to combine the MUSIC-X and ICCK model to estimate S_u in 2D space by the fusion of geotechnical and geophysical data, is applied to a real case study at Hollandse Kust (west) wind farm zone to demonstrate its effectiveness. The results indicate that such a scheme can robustly estimate a 2D cross section of S_u with quantified uncertainty. A comparative analysis is conducted between the proposed scheme and two alternatives, one lacking preliminary S_u characterization (i.e., without MUSIC-X modelling) and one lacking geophysical data, confirming that the proposed scheme has a relatively high accuracy in the estimated cross section. The research reveals it is sensible to combine MUSIC-X and ICCK for 2D S_u characterization and brings a new perspective that integrating geotechnical and geophysical data is promising to characterize soil parameters in higher dimensional space.

The Netherlands is prone to flooding as more than a quarter of the country lies under sea level. To combat flooding and ensure that the country remains dry structures such are levees and dikes have been installed. However, older water retaining structures are more than ever failing the stringent safety standard assessments. These older conventional reinforcement measures, including berm constructions, are not only costly but require an expanse of ground to ensure performability.

Backward erosion piping is an internal erosion mechanism during which shallow pipes are formed in the direction opposite to the flow underneath water-retain structures as a result of the gradual removal of low cohesive material by the action of water. This mechanism is an important failure mechanism in both levees and dams where a cohesive layer covers a sand layer. Although failure resulting from backward erosion piping is not common, several levee failures in the United States, China and the Netherlands have been attributed to this mechanism.

There are mitigation measures known to stop the backward erosion mechanism. One such measure is the placement of a seepage wall, to create a physical barrier directly in the flow path trying to reach the lowest region of the hydraulic head. A review of the literature showed that current design rules only consider groundwater flow calculations when determining the likelihood of hydraulic heave, one of the failure modes within the backward erosion process. Hydraulic heave in the backward erosion piping context is closely linked to the quicksand condition, essentially stating that once the effective stress is zero, the sand particles become suspended, liquifying a solid layer. The absence of an assessment of the effective stresses during the design process in conjunction with hydraulic heave has contributed to the main research question addressed by this thesis; How does a restricted exit for groundwater flow affect hydraulic heave compared to Terzaghi’s free exit situation?.

Flood protection infrastructure requires constant investments to cover the increasing flood risk. However, due to over-conservatism in (dyke) safety assessments, poorly targeted investments can be made. Over-conservatism can be avoided by understanding the entire failure process, from the initiation of failure until flooding. Dyke slope instability is one of the main initiation mechanisms evaluated during a safety assessment. Following an initial instability, a slope failure occurs, where large deformations may occur as the failure mass slides along the failure surface. A large initial failure mechanism may immediately trigger flooding, but in most cases secondary mechanisms, such as new slope failures, are required to flood the hinterland. The dyke may have enough resistance to prevent secondary mechanisms and thereby prevent flooding. Therefore, dyke assessments can be optimised by assessing the potential for secondary failures.

The standard methods for dyke slope stability assessment cannot model large deformations. This thesis therefore develops and applies the Material Point Method (MPM), a large deformation variant of the Finite Element Method, to investigate the residual (remaining) resistance of a dyke against flooding after an initial slope instability. The residual dyke resistance has been assessed within a risk-based framework using the Random MPM (RMPM), which accounts for the effects of soil heterogeneity on the failure process by combining random fields with MPM. From the realisations of an RMPM analysis, both the probability of initial failure as well as the probability of flooding may be determined. Moreover, with RMPM, the likelihood of failure processes can be evaluated such that the process between initial failure and flooding can be understood.

To model the external water level in the RMPM analysis, the application of boundary conditions in MPM has first been investigated. The thesis shows that the boundary conditions should systematically match the MPM discretisation. Improvements of MPM, such as the Generalized Interpolation Material Point Method (GIMP), often change the discretisation. Therefore, the accurate application of a boundary condition can therefore depend on the version of MPM being used. Consistent boundary conditions are described in this work for MPM and GIMP. For standard MPM, a consistent boundary condition is proposed for simple 1D problems. However, it is shown that this solution is not generally applicable for dyke slope failures or other higher dimensional problems. For GIMP, two generally applicable algorithms for (almost) consistent boundary conditions are proposed: one algorithm constructs the exact material boundary, while the other merges the support domains of all material points. The algorithms are shown to outperform other boundary condition methods presented in literature.

The residual (dyke) resistance has been investigated by modelling both a 2D dyke failure and 3D slope instability using RMPM. It is shown that secondary failures (required to trigger flooding) often do not occur or may not be large enough to trigger flooding. Therefore, the probability of flooding can be significantly lower than the probability of an initial failure due to residual dyke resistance. In the best case scenario for the problem analysed, a reduction of the probability of flooding compared to the probability of initial failure of more than 90% has been observed, while in the worst case only a 10% reduction was found. The reduction was high (90%) for a material without layering of the spatial variability of the strength properties and decreased when the spatial variability was more layered. However, note that, to reduce computational costs, the probability of initial failure was unrealistically high in these examples, i.e. the dyke was relatively weak. In stronger slopes, secondary failures are less likely and more residual dyke resistance is therefore expected. Additionally, secondary slope failures are less likely in 3D simulations compared to 2D simulations, generally due to the additional resistance of the sides of the failure surfaces (the so-called 3D-effect). A 2D simulation can therefore be seen as a conservative estimate of the residual dyke resistance. In 3D, the failure process more often spreads sideways rather than backwards. This is also beneficial for dyke slope stability assessments, where backward failures are required to trigger flooding.

The degree of anisotropy of the soil heterogeneity changes the expected failure process. For smaller horizontal scales of fluctuation, i.e. less layering of the soil, secondary failures are less likely to occur, since the initial and secondary failures are mostly uncorrelated. Additionally, in the 3D simulation, smaller horizontal scales of fluctuation triggered small failure blocks, again likely to reduce the risk of flooding. For larger horizontal scales of fluctuation, initial failure in a weaker layer can more easily trigger secondary failures through the same layer, thereby decreasing residual dyke resistance. A depth trend, i.e. a linear increase with depth, in the mean resistance of the material, typical due to compaction processes, also impacts the failure process. For a material without a depth trend, progressive failure occurs along approximately circular failure surfaces, whereas for a material with a depth trend, a steady flow like behaviour along a gentle ’straight’ slope occurs. Moreover, retrogressive failure can flow in any direction for a material with a depth trend while avoiding local strong zones.

This thesis highlights that RMPM can provide estimates of the residual dyke resistance, thereby more accurately estimating the probability of flooding due to dyke slope instability in many situations. This leads to more targeted and cost effective dyke reinforcements. RMPM also provides insight into the size and shape of the initial and subsequent failures. RMPM can therefore be used in future research to develop guidelines for practice to approximate the probability of flooding, for example based on the probability and the shape of the initial failure computed with a small deformation model.@en

The need of shear strength measurements of soil in the design phase of geotechnical engineering is almost indispensable. Many methods have been applied to estimate the shear strength of soil, including various laboratory test, in-situ test and analytical methods. As an in-situ test method, cone penetration test (CPT) is a powerful and cost-effective tool for the investigation of subsoil conditions. CPT data is usually complemented by the laboratory test data for verification. The laboratory-based studies of subsoil, however, can be not only a complex but also tedious and expensive task for large projects involving large amount of data. Therefore, new approaches for estimating the soil shear strength are demanded. Having demonstrated superior predictive ability for many material properties compared to traditional methods, machine learning methods have been increasingly popular and widely used. This thesis focus on the prediction of soil undrained shear strength through cone penetration test data. The major objectives of this master thesis include testing how machine learning could help us lower the need for laboratory test data. At first, the research starts with a literature review of various methods used to evaluate the soil shear strength. Comparing to the machine learning methods, the laboratory and in-situ test methods are relatively more time-consuming, costly and labour-intensive. And the analytical methods are considered lacking in precision. Then the training dataset which consists of 526 samples is introduced. In each sample, there are four input variables obtained from cone penetration test, namely the effective stress (σ′v ), cone tip resistance (qt − σv), effective cone tip resistance (qt − u2) and the excess pore pressure (u2 − u0). The undrained shear strength obtained from laboratory test is taken as the output variable. Next, the training dataset is fed to five machine learning techniques, namely the artificial neural network, support vector machine, Gaussian process regression, random forest and XGBoost, to train models. The hyperparameters are tuned with k-fold and group k-fold cross-validation strategies in the validation process. After that, the testing dataset which consists of 20 samples is established. Cone penetration test data that are in close vicinity to the location of the samples are processed by Gaussian process regression to obtain representative cone penetration test data at the sample location, which is taken as the inputs in the testing dataset. The undrained shear strengths of the samples are measured by Consolidated-Undrained shear test and are taken as the outputs of the testing dataset. Finally, the five machine learning models are tested on the testing dataset. The crossvalidation results, together with the prediction results of the models on the training and testing dataset are evaluated, gathered and compared by various statistic metrics to show the relative performance of the models. XGBoost appears to be the most accurate of all the tested algorithms on this dataset. And Gaussian process regression is chosen as the second option due to its ability to capture uncertainties. The robustness of these two models are then validated from a statistical point of view by applying Monte Carlo analysis. The importance of the input parameters in this study is evaluated by applying random forest for the sensitivity analysis. The results from random forest indicate that the excess pore pressure and the cone tip resistance - total vertical stress are the most influential inputs to the undrained shear strength

A random field generator based on Local Average Subdivision (LAS) method is proposed in this study in order to achieve probabilistic soil classification and quantify the uncertainty of the generated most probable geological cross section. CPT data and Robertson’s soil classification chart (1990) are adopted to classify the soil. The sole application of LAS makes the random field unconditional, which has been improved to conditional random field generator by using Kriging interpolation. Both unconditional and conditional generator are tested in an illustrative example and the results indicate that the improvement from unconditional to conditional random field reduces the uncertainty of the most probable result of classifications and the classifications in the unconditional random field will converge if there are enough realizations. Additionally, the conditional random field generator is further applied in a case with three conducted CPTs, which build up a domain with very large scale of fluctuations. It’s found that the uncertainty of the generated most probable result of classifications is pretty low so it’s speculated that the proposed generator can be best applied in a large scale of fluctuation scenario. Another finding in the case study is that the proposed random field generator can be used to verify the reliability of conducted CPTs.

The numerical modelling of a cone penetration test (CPT) has long been a challenging task due to the large deformations associated with the penetration of a CPT. Recent developments in advanced numerical methods have shown promising results in overcoming these difficulties by using the Material Point Method (MPM). In this thesis it is researched whether the MPM is able to reliably produce CPT results in dry sand by using a state-dependent constitutive model. Calibration chamber (CC) tests are modelled for dry sand and results are compared with experimentally performed CC tests in the laboratory. Features regarding the numerical setup and applied boundary conditions which quantitatively influence modelling results are identified and assessed before the model is validated to real CC test data. Validation results show that the model is able to accurately produce cone resistance values for different types of sand for soil states that can be categorised as moderately-dense to dense. Last, it is shown how parameters within the constitutive framework affect the model output and a quantification of the sensitivity of the parameters to model results is presented.

The Dutch Water board revised their guidelines for the safety assessment of dikes in 2017. A major change for the assessment of macro-stability is the use of the SHANSEP method to estimate the strength of impermeable cohesive layers. The failure probabilities are estimated with a deterministic analysis using limit equilibrium methods. The use of design values and safety factors to account for uncertainty is basic and proven to be conservative leading to over engineering.

In this thesis, it is investigated how the SHANSEP method can be incorporated to the more advanced Random Finite Element Method. It is found that three random fields for SHANSEP parameters S,m and POP are required. The random fields do not show particular trends in mean or standard deviation. A random field generator is coded in Python. a simple version of the in-house FEM is modified to read the generated random fields. This code is used to test various geotechnical assumptions.

A final version of the assumptions is coded into a more advanced version of the simulator to do the comparison. The output of the FEM code are the FOS and failure mechanism of a single evaluation with a combination of three random fields. A mean and standard deviation of the FOS results are calculated. The probability of failure is estimated by the area under the probability density function of a lognormal distribution for values below unity. The probability of failure of the deterministic case is estimated using the First Order Second Moment method.

The results show that the probability of failure is overestimated in a FOSM analysis by one order of magnitude compared to the most conservative RFEM simulation. It is expected that this difference is even higher for the more conservative deterministic approach the Dutch guidelines prescribe. The slip surfaces of RFEM were found to be similar to their deterministic counterpart. The RFEM slip surfaces went through local weak zones in random fields.

It is recommended to Dutch policy makers to investigate the use the random finite element method. Although conservatism is preferable in safety assessments, an conservatism of this significance compared to the RFEM approach is unnecessarily costly.

The material point method (MPM) is gaining increasing amounts of attention due to its capacity to solve geotechnical problems involving large-deformations. While some problems require dynamic analysis, simulating the (infinite) continuous domain using typical Dirichlet (fixed) boundary conditions induce spurious reflections causing (1) unrealistic stress increments at the domain oundary and (2) the appearance of multiple unnatural stress waves in the domain. Aiming to eliminate this numerical artifact in MPM, two solutions for absorbing boundary conditions found in FEM are implemented and investigated; these are (1) a viscous boundary condition and (2) a viscoelastic boundary condition. The use of such dynamic boundary conditions in MPM is scarce and no validation of them has yet been presented in the literature. In this work, these absorbing conditions are implemented alongside other recent developments, which improves the numerical stability (Double-Mapping, Generalized Material Point Method, Composite Material Point Method), using two approaches: (1) directly imposing at the external active boundary nodes and (2) imposing via shape function interpolation. The proposed solutions are then validated with a one-dimensional benchmark: a soil column under dynamic load in small-deformation and large-deformation, and a 2D symmetric plane strain model under one loading pulse. The benchmark results demonstrate that the numerical reflections that lead to inaccuracies of stress and velocity can be removed using GIMP interpolation and, together with the other numerical technics, render high-quality and realistic results. A study of shallow foundation failure under repeated loading is also presented, showing the potential of applications of the proposed solution for modelling extreme geotechnical events.

The deterministic approach for interpreting CPT soil profiles poses the serious limitation of not taking data uncertainty into account. Therefore, a Bayesian model was developed by Wang et al. (2013) that, for a given CPT profile, determines the most probable number of soil layers and most probable soil layer thicknesses by simulating and comparing multiple ‘model classes’ with different complexities. In this study, this proposed model is implemented into the Python coding environment after which the functionality is verified by conducting a case study on a 23 푚 CPT profile from the Groningen area (NE Netherlands). For the given CPT profile, the model distinguishes 6 separate soil layers from which the position and thickness are in agreement with the deterministic analysis and the available borehole data. However, the case study suggests that the model fails to correctly identify the most probable soil types for CPT measurements within the vicinity of the edges of the Robertson chart. This is most-likely related to a “cut-off”-effect of the joint Gaussian distribution describing the uncertainty of a single datapoint. A subsequent study on the integration of the statistical parameters within the model is therefore required. Additionally, the code includes several optimizing strategies, but remains time consuming for very complex model classes. Further optimization is suggested to achieve greater model precision and efficiency.

This thesis tried to quantify the strength reduction of the soft soil layers of the experiment through a back analysis of the failure of the dyke without sheet pile wall. Residual strength hypothesis were formulated based on a literature study and the current design norm of dyke design on soft soil layers in the Netherlands. A Limit Equilibrium Method analysis of the pre and post failure geometry served as a basis to determine the peak and residual strength of the soft soil layers following the residual strength hypotheses. The 3D effect was taken into account in these analyses, and determined to be around 20\%. A Material Point Method model of the experiment without sheet pile wall was then created to test the different residual strength hypotheses. Factors were applied on the strength properties calculated from the LEM model to match the MPM model. A factor of 1.23 was applied to account for 3D-effect and 1.16 to correct for water on passive side being absent from MPM discretization.

The Material Point Method model showed that during failure the behaviour of the clay layers could be expressed with an Undrained SHANSEP formulation. In this formulation the residual strength of clays was found to be independent of the Over-Consolidation Ratio. In a Mohr-Coulomb formulation this results in a complete loss of cohesion. Leaving the OCR out of the strength formulation of clayley layers resulted in horizontal displacement of 4.5 m, which is close to the 6 to 8 m found during the experiment. Further decrease in S-ratio of 30\% resulted in horizontal displacement going up to 7.5m in the MPM model. A reduction of 0 to 30\% of the S-ratio could therefore be concluded to be a range of friction softening. This was concluded to be in accordance with what was found in literature. Laboratory testing and correlations based on index properties effectuated prior to the experiment showed residual friction angle around 30\textdegree. The residual strength backcalculated are much lower, and therefore in contradiction with the laboratory testing results effectuated. The use of cyDSS and LDSS tests was therefore deemed inappropriate for the determination of residual strength.

The Limit Equilibrium Method analysis of the dyke with sheet pile wall was deemed inappropriate for the back analysis of the soft soil layers. The horizontal forces induced by the soil-structure interaction cannot be disregarded. It is recommended to back calculate the peak and residual strength of the peat layer using a Finite Element Method analysis. The displacement measurements of the failing dyke with sheetpile wall in the peat layer showed similarities with strain localisation in a DSS test.

Seismic design codes are currently moving from a force-based design approach to a performance-based design approach. For example, in a performance-based design approach it could be specified how many lanes must be available during the lifetime of a bridge given a certain earthquake intensity.The problem with this approach is that it is not specified what the probability must be that the performance criterion is satisfied. This raises the question whether the design codes are acceptably safe or not. Focus is laid on the Canadian Highway Bridge Design Code (CHBDC), in which a total resistance factor approach is used. Because the total resistance factor in the CHBDC is a multiplicative factor, lower resistance factors lead to stronger foundation designs. The goal of this thesis is to calibrate the design procedure in the CHBDC for geotechnical systems under seismic loading, by finding a relationship between resistance factors and the lifetime probabilities of failure of said systems. The resistance factor can then be fine-tuned to a lifetime probability of failure that is consistent with the lifetime probability of failure targeted in static design. As an example problem, the bearing capacity of a shallow foundation on a clay with a pseudo-dynamic earthquake load is tested. The research question that is answered in this thesis is: "What should the resistance factors for geotechnical seismic design be in order to achieve a target lifetime probability of failure that is consistent with static design targets?'' Not every possible combination of soil strengths and forces on the superstructure can be taken into account, and therefore the random finite element method is used in a Monte Carlo simulation. Thousands of realization sare performed for each resistance factor, design return period, and "actual''return period that the designed foundations are tested against. By seeing how many realizations of the Monte Carlo simulation fail given a certain earthquake intensity, the conditional probability of failure given that earthquake intensity can be estimated. The total lifetime probability of failure can then be estimated from the conditional probabilities of failure with the total probability theorem. As part of a parametric study, the lifetime probabilities of failure are estimated for six different scenarios, each of which has different sources of uncertainty.  The resulting lifetime probabilities of failure are interpolated in order to find a resistance factor that targets a lifetime probability of failure consistent with static design targets. Currently, the resistance factor that the CHBDC recommends for geotechnical systems under seismic loading are defined as the static resistance factor for that geotechnical system incremented with 0.20, meaning that compared to static design, weaker foundations are designed for seismic load cases. The resistance factor found in this thesis is closer to the resistance factor for static design than to the resistance factor for seismic design. It should therefore be considered to lower the seismic resistance factor to the value of the static resistance factor so that a sufficient lifetime reliability can be targeted.

Displacement control is of utmost importance in deep excavation design and is usually based on numerical modelling, e.g. Finite Element Method (FEM). Numerical methods tend to be more conservative when analysing soil behaviour during deep excavation, whereas for practical and economic reasons this is not favoured. The inverse analysis allows for the identification of the soil parameter set that can provide the measurements observed in the monitoring when it is applied in the model. When performed in a probabilistic concept, it reduces parameter uncertainty and enables the stochastic prediction of future soil behaviour. In this thesis, capabilities and limitations of difference advanced constitutive models are investigated. The Generalized Hardening Soil Small strain model (GHS) presents a positive aspect in modelling soil behaviour during deep excavation with its various stress/strain dependency settings. Because of the uncertainties originating from the size of the domain and limitations of site investigation, the soil parameters can only be shown as probability distributions. In order to make that distribution more accurate, comparative selection of several inverse analysis optimization algorithms is performed. Thereafter, choice of the relevant parameters is done based on the conducted sensitivity analysis and engineering judgement. Having the most competitive optimization approach selected, remote scripting with Python is used to utilise Finite Element (FE) modelling in the 2D Plaxis software. The input parameters are iteratively updated with response observation (diaphragm wall deflections) using the Ensemble Kalman filter optimisation algorithm based on a chosen excavation stage. The re-calibrated parameters are checked with the data, which was used to create synthetic measurements made using the same FE, to perform reliability assessment of the developed Python-based algorithm and investigate its capabilities and limitations. The further development of the presented optimisation method is expected to increase certainty in setting alarm thresholds in the applications of the Observational Method.

The earthquake-induced liquefaction is a high-risk phenomenon for dredging industries, which need to set strict requirements in order to avoid potential disastrous effects for the project. Different types of liquefaction exist which can be triggered over a wide range of soil types and for different loading conditions. The liquefaction triggering due to an earthquake event is dependent on the soil behaviour under undrained cyclic loading. The assessment of the liquefaction hazard during an earthquake is mainly based so far on empirical procedures. The most common used in practise is the NCEER method (Youd & Idriss, 2001) which is established according to empirical evaluation of field observations and in-situ testing. However, the NCEER method can be inaccurate for the design primarily due to its empirical nature as it is capturing different soil types and loading conditions. For that purpose, advanced constitutive models can provide more precise assessments as they can be calibrated for specific site conditions. Such a model is the PM4Sand, which is very attractive for practical applications because there are only a few model parameters to be determined in the calibration process.The first part of the current thesis project includes the validation of the PM4Sand model for both earthquake-induced and static liquefaction according to undrained Cyclic Direct Simple Shear (CDSS) tests and undrained Direct Simple Shear (DSS) respectively, performed on Ottawa F-65 Sand. The influence of the model parameters is examined throughout a parametric assessment analysis. It was observed, that the model approximates well the general features of both cyclic and static loading. Regarding cyclic loading, it produced similar responses in terms of excess pore pressures generation and stress paths even though it slightly overpredicts the cyclic resistance for small number of loading cycles and underpredicts the cyclic resistance for large number of loading cycles. Regarding static liquefaction, even if the model had initially overestimated the response, it was able to simulate successfully the static liquefaction behaviour after a recalibration process was established.The next part of the project includes the performance of the PM4Sand model for the prediction of earthquake-induced liquefaction in hydraulic fills, which are analysed for several different seismic motions. The fill is placed over a different range of relative densities and it is modelled in Plaxis software as a 1-D soil column. The fill layers that are prone to liquefy, are modelled with the PM4Sand model while the layers that are not susceptible to liquefaction are modelled with Hardening Soil Small (HSS) model. The PM4Sand layer is calibrated according to factors that are accounting for the in-situ state of the fill and the magnitude of the earthquake motions. The dynamic analyses are performed with and without consolidation and the lateral boundaries used are tied degrees of freedom. The results in terms of excess pore pressures generation are examined throughout the whole earthquake motion. Moreover, the onset of liquefaction in the hydraulic fill is captured when the excess pore pressure ratio has reached a value of around 1.0 (ru≈1). It is shown, that the PM4Sand model is indeed applicable for the prediction of earthquake-induced and static liquefaction in hydraulic fills. The effect of the in-situ state of the fill, in particular the relative density, has a critical role on the liquefaction susceptibility, which is a lot representative to what has been observed in reality. According to PM4Sand model, the loosely-packed fills DR=30% and DR=40%) are indeed more susceptible to liquefaction compared to the densely-packed fills (DR=50% and DR=60%) which showed less or even no liquefaction potential due to the earthquake events. On the other hand, the largest drawback of the NCEER method it its empirical nature, as for the current project it is proved to be conservative for the design. More specifically, it predicted liquefaction for almost all the hydraulic fills (DR=30% to DR=60%) analyzed for all different earthquake motions. Regarding the dynamic analyses with consolidation, the results related to the earthquake-induced liquefaction of the fills are more representative to realistic conditions as there is a better distribution of excess pore pressures along the soil column with respect to the dynamic analyses without consolidation. For the latter type of analysis, in the loosely-packed fills (DR=30% and DR=40%) there is a better diffusion of excess pore pressures more for the signals of low dominant frequencies regardless the peak ground acceleration values of the input signal. In the densely-packed fills (DR=50% and DR=60%) the same phenomenon takes pace more for the signals of high dominant frequencies. However, a localization of liquefied zones is observed in distinct parts along the fill layer for the rest of the signals.

This research has shed light on the capability of the HASP model in reproducing the dilatant behaviour of OC clays in drained and undrained conditions. The model is formulated by employing a combined hardening rule and uses the void ratio as a state variable while maintaining the simplicity of the MCC model. However, a sensitivity analysis has revealed that the model is sensitive to some input parameters which when varied slightly can largely affect the outcome of an analysis.

This has led to the formulation of the PLAXIS OC clay model while maintaining the framework of the HASP model but replacing the void ratio with volumetric strain as the state variable. Thus resulting in the use of the modified compression (λ*) and swelling (κ*) indexes which are used in obtaining the soil stiffness parameters used as model inputs. The PLAXIS OC clay model features the use of real soil stiffness parameters instead of soil indexes, the addition of small stain stiffness by T. Benz to improve model prediction in the small strain region and the elimination of the sensitivity issues noticed when using the HASP model.

The PLAXIS OC clay model is validated for boom clay (BC) at single stress points by simulating CU test and comparing with the available experimental data for the BC. Good agreement is found with experimental data as shown in the stress strain, pore water pressure and stress path plots obtained from the analysis.

Furthermore, the model is used to simulate boom clay in a trial excavation. Piezometers and extensometers are installed into the BC layer prior to the excavation to monitor the changes in porewater pressure and vertical displacement (heave) on the BC during the excavation. A comparison of the numerical and experimental data shows that good agreement is observed in porewater pressures and vertical displacement in the BC layer.

The thesis presents a numerical study on dredging induced undrained instability and subsequent static liquefaction of submarine landslides. For the study, a pre-existing hydro-dynamic uncoupled submarine slope numerical model, developed by Molenkamp (1999), has been modified to incorporate a fully hydro-dynamic coupled interaction between submerging water mass and submarine slope. The modified model is able to simulate transient quasi-static and dynamic phenomena up-till and including the immediate post-liquefaction behavior of submerged slopes of loose undrained homogenous fine sands in a 2 dimensional Updated Lagrangian (UL) finite element (FE) frame of work. To simulate soil behavior under dredge loading applications the model incorporates a Monot soil constitutive model and for submerging water behavior a Lagrangian expression of Navier stokes for nearly-incompressible visco-elastic, irrotational, fluid model.
The study primarily addresses the effect of dynamics of submerging water on the liquefying submerged slope. The research findings suggest that the dynamic motion of submerging water barely affects the occurrence of instability. However, it may decrease the rate of post-instability liquefied flow as compared to the commonly sorted uncoupled scenario, where dynamics of submerging water mass is ignored and only constant hydrostatic pressure heads due to water level is considered at the slope interface. Moreover, the findings suggest that about 50% of the loss in the potential energy of soil is consumed by the potential energy of the submerging water at the very initial stages of post-instability and that the contribution of kinetic energy of water amounts to mere 3.4%.
Next, as a secondary issue, the study also provides a valuable insight into the effect of the liquefying slope on the motion of the submerging water mass. The findings show a surface impulse wave formation post-instability, moving along the direction of landslide. Moreover, it shows a development of a distinct circular motion of fluid along the slope interface. Other than this, the thesis also attempts to provide some similarities and differences between the current findings and the published conventional research studies which make use of basic slide shapes such as viscous or rigid sliding wedge blocks.
Finally, the thesis also addresses some numerical shortcomings such as the hour-glass effect, the shake-down by the procedure to define the “initial state” effect etc., and thereby providing necessary recommendations useful for future computational modelling work.

The subsurface is often modelled as multiple homogeneous soil layers. In reality, the soil properties of these layers are heterogeneous and spatially variable, making it difficult to properly deduct the engineering properties that are required in the design of geotechnical structures. This leads to uncertainty in the design. The probabilistic framework known as the random finite element method (RFEM) can take the spatial variability of soil properties into account in geotechnical designs in an effort to reduce the uncertainty. The spatial variability is represented by the mean, the standard deviation and the scale of fluctuation. The scale of fluctuation is defined as the distance in which a significant correlation is found within the soil properties and is related to the geological deposition of the soil layer. The calculated reliability of man-made structures, such as dykes, is the output of the RFEM and is influenced by the scales of fluctuation of the subsurface. Since many man-made structures in the Netherlands are constructed on soft soils, they cause deformations in the subsurface. These deformations change the geological layering of the soil. It is therefore expected that the deformations caused by man-made structures influence the scales of fluctuation and that the changed scales of fluctuation influence the calculated reliability of the man-made structures. The objective of this thesis is to research how the spatial variability of soft soils is influenced by man-made structures. To research this, a case study is undertaken at the Leendert de Boerspolder in the Netherlands. One-hundred cone penetration tests (CPTs) have been taken in an area of 15 by 50 meters in the polder as well as in the dyke as part of the research programme Reliable Dykes. From the cone resistance data of the CPTs the spatial variability of the area can be estimated. The case study is accompanied by two theoretical studies. The first study aims to estimate the effect of deformations on the spatial variability of generated random field data for different types of deformations. The second study of this thesis will provide the deformations that have occurred in the area of the case study from a plane-strain finite element analysis of the dyke and polder.
The first part of the research focussed on four deformations that are often encountered in the proximity of a dyke: vertical compression, horizontal compression, simple shear and rotation. The effect of these deformations on the spatial variability has been studied by applying the deformations to three resolutions of random fields and then estimating the change in scales of fluctuation. Since no soil properties were assigned to the random fields, the results can only be used qualitatively. The research is then continued with the estimation of the spatial variability of the Leendert de Boerspolder. It was necessary to identify the soil layers first. These are a top peat layer, a clay layer which shows a transition from organic at the top to silty at the bottom, and a bottom peat layer. For each of the CPT lines parallel to the dyke the outliers of the cone resistance data were removed for the different soil layers and the depth trends were removed. This resulted in normalised cone resistance data from which the spatial variability in different directions could be estimated. Finally, a two-dimensional plane strain analysis of the soil layers at the Leendert de Boerspolder was made. In the first step of the analysis the depths of the soil layers are at the level before construction of the dyke. Then the dyke was constructed step-wise. This allowed the initial soil layers to deform as close as possible to the actual deformations under influence of the dyke. Finally, the following conclusions are made by combining the resulting deformations with the results of the other research. At the Leendert de Boerspolder, the spatial variability of the polder next to the dyke is not influenced by the dyke and can be used as the baseline spatial variability to which to compare the spatial variability in the proximity of the dyke. Underneath the crest a significant vertical compression causes both the vertical and horizontal scales of fluctuation to decrease, while a smaller amount causes a change in just the vertical scale of fluctuation. Underneath the slope and toe of the dyke, rotation increases the larger vertical scale of fluctuation and decreases the larger horizontal scale of fluctuation. Horizontal extension of the soft soil underneath the slope results in a larger horizontal scale of fluctuation. Horizontal compression underneath the toe results in a smaller horizontal scale of fluctuation.