The energy crisis and climate change, is creating an urgent need for sustainable and energy-efficient solutions. A significant portion of global energy consumption comes from households. Within this sector, the largest share is used for space heating and cooling. Therefore, efforts to reduce energy usage and cut carbon emissions should primarily target heating and cooling systems. A promising solution to this challenge lies in harvesting shallow geothermal energy through technologies such as energy piles, which have a multifunctional role to the building.

Energy piles are a specialized form of Borehole Thermal Energy Storage (BTES) systems that utilize shallow geothermal energy by taking advantage of the ground's stable temperature throughout the year. They are gaining popularity as an efficient solution for both heating and cooling, primarily because they serve a dual function within a building. On the one hand, they act as heat exchangers, enabling the transfer of thermal energy to and from the ground. On the other hand, they provide structural support, as they are typically constructed from reinforced concrete. This multifunctional role makes them a cost-effective choice and reduces the initial costs investments. These systems are commonly integrated with Ground Source Heat Pumps (GSHPs), which facilitate the exchange of thermal energy between the energy piles and the circulating fluid within the system. Their high Coefficient of Performance (COP) and environmentally friendly operation contribute positively to the transition towards more sustainable energy solutions.

This research focuses on an existing energy pile system located beneath a building on the TU Delft campus, which is responsible for covering the building's heating and cooling demands. The main objective of this research is to assess how the efficiency and reliability of this system can be improved through the implementation of advanced control strategies. To achieve this, a precise simulation model was developed to capture the underlying physical processes and monitor key performance indicators. Energy balance is a key metric and evaluation point for the system that secures the efficiency and sustainability of it. Additionally, various scenarios with different operational parameters were designed to evaluate the system's capabilities and performance.

The main findings initially indicated that the heating load is higher than the cooling load, revealing a significant imbalance. By creating various scenarios with different temperature setpoints, heating and cooling months, and durations of heating and cooling modes, an energy balance was achieved. However, the system was still unable to adequately cover the energy needs of the building during winter. Scenarios that utilized solar gains —by opening the building’s sunblinds— enabled the heating and cooling system to supply the majority of the heating load during winter. This aligns with one of the main goal of the system while maintaining energy balance throughout the year. Among the scenarios evaluated, scenario 5 demonstrated the highest heating and cooling energy delivery. According to performance evaluations, it also consumed the least electricity among the compared scenarios. Additionally, Scenario 5 provided the highest levels of visual and thermal comfort for occupants, as the sunblinds remained open during non-operational months. The system was found capable of operating under loads 50% higher than normal, although these levels push its operational limits. Moreover, it was observed that the energy piles system can store surplus energy during summer and completely cover the heating demand of a neighboring apartment.

The thermal plume interaction is investigated thoroughly and shows the rate and the magnitutude of the expansion or contraction of those. It was observed that energy is lost due to the open boundary with the ambient air temperature but varies between the scenarios. More details will be discussed in the following sections.

The doctoral research has been carried out in the context of an agreement on joint doctoral supervision between RWTH Aachen University, Germany and Delft University of Technology, the Netherlands.

Geothermal energy stands as a promising avenue for low-carbon energy production. Traditional methods involve extracting hot water into one well and the injection of cold water in another, with hydrothermal systems relying on fluid flow through either pores or natural fractures in rock. Enhanced Geothermal Systems (EGS) are deployed in fractured/faulted rock settings where fluid flow is insufficient, necessitating enhanced permeability of fractures for improved heat transfer. Challenges persist in ensuring productivity, sustainability, and safety, with fault and fracture reactivation presenting significant concerns....

Slope stability applications are vital assets for a country. These slope stability systems include dikes, dams, levees, embankments and enable applications such as open-pit mining. The failure of these systems pose huge impacts on society and the economy and hence the accurate stability assessment of these systems are of primary concern. Existing methods such as limit equilibrium methods, numerical methods (e.g. the finite element method, FEM), empirical methods and probabilistic methods all provide an approximate estimate of the factor of safety (FoS), and are often observed to have inaccuracies when failures occur. This lack of ability to make accurate predictions is due to many reasons, such as missing physical processes incorporated into the methods, inaccurate boundary and initial conditions, constitutive model selection, uncertainty in model parameters and limited mechanism understanding. This thesis suggests using data assimilation to combine monitoring data with a finite element model to improve the predictive capabilities of the FEM model. These days, geotechnical systems are equipped with measurement devices to monitor their response to external changes. These measurements can be in the formof surface displacements, porewater pressures, strains, etc. These measurements can be obtained from in-situ devices (such as inclinometers, strain gauges, etc.) or can be measured remotely (with Light Detection and Ranging (LIDAR), Interferometric Synthetic Aperture Radar (InSAR), etc.). These measurements can be assimilated into the popular ensemble-based well-established data assimilation methods, e.g., the ensemble Kalman filter (EnKF), ensemble smoother (ES) and ensemble smoother with multiple data assimilation (ESMDA) to improve the predictability of FEM models.

In the first stage, an FEM model of slope stability has been integrated with EnKF. Based upon the slope deformation measurements, this approach estimates the key material parameters (strength and stiffness parameters), the state (displacement), and the FoS of a slope. The effect of two different constitutive models (Mohr-Coulomb (MC) and Hardening Soil (HS) model) on the FoS was studied via a synthetic twin experiment. The HS model was able to estimate the FoS with a narrow posterior distribution, starting from a wide prior distribution of material parameters, including those not encompassing the actual parameters, demonstrating the advantage of using advanced constitutive models when combining with data assimilation.

In the second stage, the constitutive model which produced relatively more accurate results (the HS model) was selected from the first stage has been tested with three data assimilation schemes, i.e., EnKF, ES and ESMDA. Each of these schemes was integrated with the FEM to assimilate measurements of deformation of the slope and the crest of the slope stability system. The accuracy of these schemes was evaluated by comparing their FoS to the synthetic true FoS and evaluating their computation time in a synthetic twin experiment. The results of the synthetic twin experiment show that EnKF estimated an FoS that was close to the true FoS with a small standard deviation. ESMDA, when using four iterative assimilation steps, was also able to estimate an FoS close to the truth, yet had a higher standard deviation compared to EnKF. The ES and ESMDA (with two iterative assimilation steps) were not able to reconstruct the true FoS as well as the other schemes, most likely due to the mostly linear updates of these schemes. The theoretical computation time required by the ES was the smallest, followed by ESMDA with two iterative assimilation steps, ESMDA with four assimilation steps, and finally the EnKF.

In the third stage, a data assimilation scheme was implemented on a case study of an open pit mine in Cottbus, Germany. The LIDAR measurements of the vertical displacements were assimilated into a FEM model of slope stability. Model parameters, displacement ensemble and FoS are estimated from this analysis. The posterior estimation of FoS is compared with slope failure observed in the field. The data assimilation results provide better results than only using FEM models when comparing the ground truth of slope failure. However, it was clear that not all physical processes were included in the model, resulting in a considerable mismatch of the modeled and observed deformations, although a considerable improvement was observed. This initial observation led to the choice of a data assimilation method, which is able to update the parameters to generally improve the results, as opposed to those which incrementally improved parameters.

Furthermore, as the data assimilation approach developed involved multiple FEM analyses, it is computationally expensive and therefore developing a real-time assessment system is likely to be impractical. Therefore, an effort was made to reduce the required computational resources by developing a surrogate model. The surrogate model was trained and tested based on the output of the FEM model ensemble. Specifically, it used the displacements at different locations as input and the FoS as output. The output of the surrogate model in the validation stage was compared with the observed FoS from the case study. It was found that the prediction made by the surrogate model was not reliable. This is probably due to the mismatch between the training/testing dataset (from FEM) and the validation dataset (i.e., the measurements from LIDAR). This mismatch was identified to be due to the identified missing physical processes in the model, and the fact that the on-ground measurements had a different nature than training and testing data. It is further suggested that a surrogate model can only be used provided the training testing and validation datasets are compatible - and as the FoS is rarely identifiable in reality leads to challenges using surrogate models to predict slope failure.

This thesis introduces a novel Generative Adversarial Network application called SchemaGAN, which has been adapted from the Pix2Pix architecture to take Cone Penetration Test (CPT) data as a conditional input and generate subsoil schematizations. For training, validation and testing, a database of 24,000 synthetic schematizations of size 32x512 pixels was created, representing a broad spectrum of stratigraphical complexity in the layered models. Each synthetic cross-section was additionally transformed into a CPT-like image with less than 1% of the original data remaining at random locations along the model. After training for 200 epochs, the best-performing SchemaGAN Generator was chosen from the validation, and the effectiveness of SchemaGAN in generating subsoil schematizations was tested against traditional interpolation methods (Nearest Neighbour, Inverse Distance Weight, Kriging, Natural Neighbour) and newer methods such as Inpainting. The evaluation metrics obtained reveal that SchemaGAN outperforms all other methods, with results characterized by clearer layer boundaries and accurate anisotropy within the layers. In contrast, Nearest Neighbour and Kriging are characterized by a lack of continuity and blurry layer boundaries respectively. Inverse Distance, Natural Neighbour and Inpainting fail to come close to the performance of the other methods. The superior performance of SchemaGAN is confirmed through a blind survey, in which SchemaGAN ranked as the top-performing method in 78% of cases according to experts in the field. Results also suggest that SchemaGAN is the least affected method by the location of CPT data along the cross-section. In a real case study, SchemaGAN demonstrates better predictive accuracy for known CPT data than both Nearest Neighbour and Kriging interpolation methods. The future potential lies in refining its performance by considering enhancements such as training with real CPT data, incorporating additional conditional inputs, and exploring larger inputs or specialized databases. All the code related to the project has been made publicly accessible.

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.

Geothermal energy in the Netherlands is growing in general, however, medium deep geothermal (200 - 1500 m depth) is lagging behind while the potential at these depths could be very high. Development is held back by little geologic knowledge of this interval as well as legislative constraints. Better understanding of the reservoir characteristics and well design in this depth range is necessary in order to accelerate the heat transition away from fossil fuels. The Breda Formation is deemed as one of the most promising formations within the medium depth range. Therefore, this formation, particularly the part found in the Zuiderzee Low (ZZL), is the focus point of this thesis. The reservoir characteristics are determined using newly obtained cutting information, novel biostratigraphic analyses, old and new seismic interpretations and petrophysical log data, eventually leading to a static reservoir model. Consequently, the static model is upscaled and used as input for dynamic simulations in the Harderwijk case study area. The static model provides new thickness, depth, porosity, permeability and transmissivity maps of the Breda Formation in the ZZL. These maps were generated using more data (types) than the current ones present in the REGIS model, and should therefore contain less uncertainty. The dynamic model provides insights into doublet lifetimes and flow rates for three different well trajectories in two possible reservoir scenarios. The horizontal trajectory proved superior in terms of doublet lifetime and flow rate for both reservoir scenarios. Finally, the current Dutch legislative framework covering the medium depth is limiting medium deep geothermal development in the Netherlands. A revision of the 500 m depth boundary separating the Mining and Water Law jurisdictions should be considered if the full potential of the medium depth is to be accessed.

New solutions for future flooding problems are becoming more important than ever before. Emergency measures against failure of flood defenses may play an important role in flood risk management and is therefore currently being investigated. One of these measures is the BresDefender, which is a submersible floating pontoon that can be placed on the outer slope of a weakened dike section to prevent or delay failure. The application of the BresDefender can be divided into two scenarios: (1) to stop breach formation and (2) to limit water infiltration into the dike. The latter application scenario is expected to influence the development of the phreatic surface inside the dike, which affects dike safety. This study focuses on this scenario by assessing the effect of this seal on the development of the phreatic surface.

A numerical modelling approach is developed to describe this effect using the finite element method program PLAXIS, with add-on module PlaxFlow, to describe transient groundwater flow problems. This numerical modelling process starts off with the development of a model for a laboratory scale dike. In a previous study at the WaterLab, the phreatic surface level of the dike was measured over time for different degrees of sealing using a steel plate. These measurements are used to calibrate the numerical model, where the connection between seal and dike slope is described with a transmissive interface layer. It is concluded that this interface layer approximates the effect of the seal on the phreatic surface adequately.
The next step is to extent the numerical model to a larger scale dike consisting of heterogeneous soils. This model represents the dike located in the Flood Proof Holland test facility, which consists of a permeable core covered by a low-permeability layer. By modelling this transition from simple to complex, influences of permeability, heterogeneity, and damage of the cover on the phreatic surface are identified.

The effect of a seal on the phreatic surface is also studied using physical model tests of the dike at Flood Proof Holland. The position of the phreatic surface is measured by pressure sensors in standpipes, which are spread over the crest and inner slope of the dike. Different scenarios for the seal are examined: stiff plate, flexible textile, and no emergency measure (reference case). For every scenario, experiments are carried out with and without a damaged location in the outer slope of the dike, which lead to a total of six test cases.
The results of the physical model tests show multiple effects of the seal on the phreatic surface. First, the seal shows a delaying effect on the position of the phreatic surface, which implies that the seal delays the rise of the phreatic surface over time. The time until the phreatic surface reached a steady-state over the entire dike was increased with around 15% for plate cases and around 25% for textile cases when compared to the corresponding reference cases. Second, no decreasing effect on the phreatic surface can be observed, which means that the phreatic surface level in its steady-state condition is not affected by the placement of a seal on the outer slope. Third, the delaying effect is larger for a seal that consists of a flexible textile rather than a stiff plate. The connection between seal and dike cover is proved to be important since leakages underneath the seal influence its performance, especially when the dike is locally damaged.
Last, the textile seal has a three-dimensional effect on the development of the phreatic surface. A certain area of influence can be identified where the phreatic surface is affected. The effect of the phreatic surface is seen to be stronger for locations near the textile and this effect diminishes over time. For the textile case on a damaged dike, an initial decrease of the phreatic level is observed directly behind the textile ranging up to 30 cm.

Overall, the effect of a seal on the development of the phreatic surface is concluded to be relatively small. The effect is only of a time-varying nature and, in three dimensions, the effect diminishes for larger distances. The application of a textile in terms of dike safety shows only a marginal improvement, which occurs only in a limited time period.

High-level waste can be radiotoxic for thousands of years and should be carefully handled to prevent accidents. The research for geological disposal of radioactive waste focuses on the construction and durability of a geologic disposal facility in Dutch clays or salts. One of the challenges is assessing these host rocks’ capability to retard radionuclides and prevent them from entering the biosphere. This research focuses on disposal in clay. Because of clay’s low permeability, water movements are slow, and radionuclides transport is expected to occur predominantly by diffusion. The Peize and Waalre formations, situated on the interface between brackish and salt groundwater, serve as a natural analogue for targeted deeper, poorly indurated host rocks. The known disparity in chlorine levels between the aquifers adjacent to these clays qualify these formations for NaCl-tracer research. A cutting sampling study determined the conductivity of these clays’ porewater. Combining this with a one-dimensional modelling study, the saline history of this formation has been simulated. NaCl gradients were demonstrated at different drilling locations. The most manifested gradient is attempted to fit in the one-dimensional model. The model results suggest that this gradient originates from adjacent aquifer salinity, the clays’ physical properties and the difference in hydraulic head. The observed salinity discrepancy between aquifers and ion concentration gradient in the firstWaalre Clay confirm the assumption that the member can be a natural analogue. Uncertainty on the continuity of the total system prohibits concluding that diffusion-dominated transport in Dutch poorly indurated clays can be assumed. The best fitting scenario this research found fitting the empirical Waalre clay salinity curve is combined transport by diffusion and advection. The finding of diffusion-advection transport implies that the Waalre clay shows more complexity than initially expected.

The extraction and processing of mineral and metal ores in the mining industry comes paired with large amounts of mine waste, also known as tailings. This waste consists of various chemicals, acids, and heavy metals which are used during these extraction processes. The tailings are usually stored off in tailings storage facilities (TSF) in a loose state, gradually consolidating over time. TSFs founded in seismically active areas are susceptible to liquefaction due to earthquake loading. Historic data show that about 35% of dam failures are due to liquefaction of the tailings, thus increasing the need to study the liquefaction responses of tailings under cyclic loading. Extensive studies have already been
conducted using cyclic direct simple shear (DSS) and cyclic triaxial (CTX) tests as these tests under constant-volume conditions can evaluate the change in pore pressure within a soil accurately. Previous study shows what importance the relative density and sloping ground conditions, known as drained shear bias, have on the cyclic resistance to liquefaction of the tailings. However, in practice the pore water is not bounded within the material and excess pore water can flow out through installed drains. A round-robin program, issued by the University of Western Australia (UWA), requested a study on the liquefaction response of a particular fine-sand tailings material. Inspired by this round-robin program, an interest raised in studying the cyclic shear response of tailings by using a direct shear box to investigate the cyclic behaviour under partially drained conditions. With use of the direct-shear apparatus of Wille Geotechnik, a test program has
been set up to study the influences of relative density and drained shear bias under stress-controlled cyclic shearing and under constant normal load conditions.
Results of the experiments met the expectations that denser soils have a 9.6% higher cyclic resistance ratio (CRR), and samples with applied drained shear bias have a 32% lower CRR compared to samples tested with level ground conditions. Furthermore, samples which underwent post-cyclic shearing showed strain-hardening responses and yielded higher shear stresses compared to the monotonic test, indicating that the constant normal load further densified the samples during cyclic shearing. However, during the experiments, it was
quickly found out that the loading frequency was not being applied optimally making it not possible to analyse influence of partially drained conditions . This study showed promise on its capabilities to study cyclic shear loading on a soil. For future work, it is suggested to perform similar tests under a uniform loading frequency with the use of a shear box to evaluate its capabilities to study on partially drained conditions. It is also recommended to conduct tests under constant volume conditions, to evaluate the shear-box apparatus’ capabilities to
study the liquefaction response of a soil due to excess pore pressure generation.

Climatic conditions in uence peak discharges in rivers and change sea levels; therefore, attention to the safety of dikes is of ever growing importance. Macro instability is one of the dike failure mechanisms that can inundate the hinterland. Soil heterogeneity plays an important role in assessing dike safety, especially for slope stability, because it is a major source of uncertainty. To assess a dike network for safety, numerical simulations for a full probabilistic analysis can be computationally expensive. Therefore, this study investigates how to build a state-of-the-art data-driven framework from a numerical model to predict the safety margins from the macro stability of dikes. Inputs and outputs of tens of thousands D-Stability simulations were used to create a training dataset. The most relevant features were selected based on global sensitivity analysis and the representation of soil heterogeneity in the framework. The maximisation of Shannon's information entropy and the generation of the training dataset was achieved by employing a smart sampling strategy for the input parameters. The sampling strategy consists of a Latin hypercube optimised uncorrelated uniform distributed dataset combined with a correlated dataset for optimal training eciency. The uncertainty due to soil heterogeneity is represented by a Gaussian random eld with a trend. This trend is commonly determined from a geotechnical cone penetration test. With a CPT, it also is possible to nd the vertical scale of uctuation, which is parametrised by the correlation length of the uctuations in soil strength. The second-order Markov correlation function is used to represent the correlation of the random elds. The Gaussian random eld is later mapped onto 16 stacked horizontal layers to model the heterogeneous soil properties. The surrogate model consists of an ensemble of thirteen machine learning models. The most important model is a multi-layer perceptron feed forward articial neural network. The other models are histogram based gradient boosting regression trees. Random search and Bayesian optimisation are used as hyperparametrisation techniques to optimise the prediction capability is of the individual ML algorithms. Weights for each model are determined based on optimisation for error reduction for maximum performance. The surrogate predicts the factor of safety (FOS) as well as the coordinates of the slop failure circles and line of depth from the Uplift-Van method. The surrogate model ensemble that predicted FOS is quite accurate with respect to the numerical FOS of D-Stability, and yet the prediction of the failure plane is still slightly worse. A case study was used to demonstrate the performance of the framework. Despite the uncertainty of the subsoil, due to the soil heterogeneity, the surrogate was able to accurately predict the failure probability. However, the prediction of the far end circle coordinates showed lower performance due to propagating errors. Concluding, application of the framework is possible for dike reinforcement optimisation, risk-based dike safety assessment, length effect, and effcient Monte Carlo simulations.

Interpretation of the thermal properties of soils is an important challenge in the field of geo-engineering, for example the development of geothermal energy solutions and for the design of electricity cable routes used for offshore wind farms. Of the thermal properties, the thermal conductivity is of most interest to find, as this determines the long-term thermal response of the soil. The soil volumetric heat capacity is of secondary interest, as this mainly influences the short-term thermal response.

To find the thermal properties of offshore soils, a new in-situ test is being developed, called the heat flow cone penetration test (HF-CPT). This test uses a module that can be attached to a cone penetration test (CPT) which contains a heating element and temperature sensors. In this test, the penetration trough the soil is stopped at a required depth, the heating element is then activated, and the thermal response of the probe is measured. This thesis presents an interpretation method that can predict the thermal conductivity of soils based on the thermal response of the HF-CPT. The interpretation method is validated by conducting laboratory tests in four different materials: moist sand, saturated sand, kaolin clay and a water-agar mixture. With the interpretation method, excellent results are found with the laboratory tests conducted in saturated sand, kaolin clay and the water-agar mixture.

The interpretation method is suitable for offshore testing, as the runtime of the method is short and the storage space is low. The interpretation method gives an accurate prediction for testing duration of about 300 seconds, which is fast when compared to other in-situ tests to measure the thermal conductivity of the soil. With this interpretation method, the HF-CPT can become a successful new in-situ test to determine the thermal conductivity of offshore soils. This way, the thesis contributes to the implementation of geothermal energy solutions and offshore cable routes for wind farms.

Tailings storage facilities (TSF) are engineered structures that retain mixed waste material (known as mine tailings material) from mining processes in liquid or slurry form. One of the issues regarding TSF management is the lack of site investigation and sample analyses, which cause the uncertainty of geotechnical properties of mine tailings materials. For a proper operation of tailings storage facilities, the owner should study the material involved and understand the physical and chemical properties associated with it. In case of mine tailings dam failure, the resulting flood wave supposes an environmental, social and economical disaster for society.

Mine tailings dam break studies use numerical models to predict the flooding area and assess the possible damaged area. Historically, these studies were carried out according to Newtonian modelling principles, but the presence of solids within the fluid suggests that the resulting flood wave of a TSF failure should be treated as a non-Newtonian fluid. Absence of laboratory data regarding the geotechnical properties of mine tailings materials make difficult the prediction of such flood wave, since the composition of the mixture is unknown. Therefore, the aim of this research is to study the flow behaviour of mine tailings materials in case of failure of tailings storage facilities. Understanding the flow behaviour of the non-Newtonian fluid is essential to analyse the possible failure event for an existing structure, in order to plan and organise emergency procedures that anticipate and mitigate downstream damages.

To protect the Netherlands better from flooding, and with an eye on sea-level rise in the rest of the world, more accurate assessments are needed for dykes. The calculation for the most occurring failure mechanism in dykes, i.e. macro-instability, is limited by not being able to calculate large deformation when sliding occurs within a dyke. The random material point method (RMPM) is able to capture the complete failure path, including the residual dyke strength, while taking the heterogeneity of the soil into account, thereby improving the assessment of failure processes. However, as the method is not (yet) used in current practice, clear communication of the results is essential for convincing a wider public of the contribution of this method. Visuals are an important tool in communication, as it increases comprehension of the subject matter if the visual is designed efficiently. This research aims to investigate the available software and which techniques are suitable to make realistic and informative visualizations for a given RMPM dataset of slope failure problems.
A method to create visualization with a certain graphical realism in three-dimensional space is developed. The technique uses a computer graphic software (Blender) combined with an add-on, i.e. an extension plugin. The add-on allows to work with VTK software to process scientific data for visualization, thereby maintaining the scientific correctness of the visualizations. Moreover, a rendering pipeline in Blender is created, which transforms the properties from scientific colors into realistic materials, making the visualizations more intuitive.
However, the dataset is too large to summarize in a straightforward illustration. Therefore, a data analysis is obtained to classify each realization into five pre-defined failure profiles, which are determined based on a literature study. Four failure profiles are classified based on the number of retrogressive failures and whether or not the realization resulted in flooding, while the fifth class describes horizontal failures. A technique has been developed to separate the horizontal failures from the other classes based on the plastic deviatoric strain attribute. Additionally, the data analysis aims to characterize the behavior of each failure profile from an early start, such that the findings could be used for current methods, which could not calculate the full failure profile.
Therefore, this thesis needs to investigate the reduction of the dataset to make it more time efficient when doing a data analysis. It is extended on the clustering algorithm, which has the function to detect failure blocks based on the displacement per dyke profile. The reduction method replaces an amount of data by one representative point per cluster. It not only reduced the size of the dataset significantly, from 3000 GB to 6 GB, it also made the comparison of attributes between realizations, and therefore the data analysis, easier.
The data analysis shows that it is hard to distinguish different failure profiles using only data of the initial failure, which shows the importance of using RMPM to account for post-failure behavior instead of using the current assessment i.e. FEM and LEM. One finding is that equilibrium of the initial failure block is often reached before a vertical crest displacement equal to 0.5 times the height of the dyke. This indicates that the crude estimation in the current assessment is highly conservative. Moreover, within the assumption, it is hypothesized that the secondary failure block will only form after the initial failure block has reached its equilibrium, which is shown otherwise within the data analysis of this thesis.
This work proposes a method for data analysis of RMPM using parallel coordinates, which can be extended to other RMPM datasets for macro-instability and can help to improve the prediction of the probability of flooding. Moreover, it proposes a method to visualize the prominent features, determined using parallel coordinates, in Blender-VTK. This work can, in future research, be extended to other geotechnical problems, such as 3-dimensional dyke slope failure.

Carbon dioxide (CO2) capture and sequestration (CCS) can play a significant role in reducing anthropogenic CO2 emissions while allowing society to slowly phase out traditional energy sources. One of the main challenge in CCS is the current absence of a clear revenue model. One idea is to use the sequestered supercritical CO2 as a working fluid for geothermal energy extraction from sedimentary reservoirs (Carbon-dioxide Plume Geothermal or CPG). Due to the large variability in density and mobility of CO2 under different temperatures, reservoir heterogeneities can give a rise to a combination of convective and conductive heat transfer. In this work, qualitative and quantitative descriptions
of the effect of reservoir properties on the performance of CPG in depleted
gas fields are provided using an example realistic depleted gas field. The primary focus is on the behaviour of the CO2 plume with regards to different reservoir properties such as porosity, permeability and thermal properties. The effect of large-scale reservoir structure, such as a presence of an aquifer, net-to-gross ratio and layering is also studied. In order to accurately model these effects, a thermal multi-component multiphase model based on a fugacity-activity Equation of State is built and validated for for pressures from 50-400 bar and temperatures from 35°C to 130°C. The developed thermodynamic model
is implemented into the Delft Advanced Research Terra Simulator. Numerous studies of 2D and 3D ensembles and sensitivity studies are carried out to examine the effects of isolated parameters on CPG performance. Results reveal that increased net-to-gross (N/G) ratio is associated with increased recovery factor. In addition, layering architecture becomes an important factor for the importance of conductive flux only at low N/G. Variations in the required pressure to sustain a production rate is associated with fluctuations
in production temperature and density due to expansive cooling. Varying reservoir properties and state also have a significant effect on brine upconing, which is detrimental to CPG performance. It appears that an increase in injection rate have a positive effect on the performance of CPG, but this should be studied in conjunction with a coupled wellbore and power plant model. Heterogeneous porosity-permeability realizations show a strong decrease in reservoir lifetime compared to their upscaled homogeneous counterparts, which is caused by a combination of preferential flow, reduced conductive flux and lower production BHP associated with the upscaled realizations. It was also found that reducing production rate delays the time of thermal breakthrough due to the combined effect of these factors.

Finite Element Analysis is an application used in many engineering fields, including geotechnical engineering. With increasing complexities of structures and material behavior, classical approach for non-linear finite element analysis has become computationally expensive (Gulikers, 2018). In order to improve the efficiency and accuracy of the model, Machine Learning (ML) is introduced. The ever-increasing accessibility and magnitude of computational power has made machine learning gained popularity since 1990s. Artificial Neural Networks have especially been proven to be a versatile machine learning framework, capable of solving complex tasks such as voice recognition to solving failure analysis of structures. This thesis aims at implementing soil constitutive models such as Linear Elastic and Elasto-plastic models and integrating them into a Finite Element Engine. Two different approaches have been attempted to implement the ANN into the constitutive models. The first approach involves the implementation of a "Generic" Neural network to simulate the natural anisotropic behavior of soils. Datasets have been generated with the help of theoretical soil tests classified as ‘Standard’ and ‘Non-standard’ tests. The standard tests depict an isotropic behavior, where the strains and stresses in different directions are proportionate. Non-standard tests consist of data representing anisotropic behavior where the strain increments in every direction are randomized. The neural network has first been implemented for a linear elastic soil model. Itis observed that when the network is trained and tested with similar kind of data, predictions are accurate. However, when the neural network is trained with a certain type of data (such as standard test), it is unable to predict the response for the other data (non-standard test). When the neural network (trained with standard data) is subsequently implemented in Finite Element analysis by deriving and substituting the material stiffness matrix from the neural network in a constitutive model, the difference between the isotropic linear elastic matrix (for linear elastic model) and the matrix derived from the neural network lead to significant errors, despite the accurate predictions of the neural network. The second approach is a component based approach in which the neural network is modeled to learn the typical behaviors of soil such as linear elasticity and plasticity. The data for training and testing the neural network are generated using standardized laboratory tests (CRS, DSS, Biaxial etc.). The neural network is modeled for linear elastic soil, by applying certain weight constraints. It is observed that the neural network has accurate predictions. Further implementation in constitutive model shows that the neural network is able to mimic the behavior of Linear elastic soil model. The neural network model, is subsequently implemented to predict Mohr-Coulomb behavior, and the predictions are significantly accurate. Further the generalization ability of the network is analyzed to evaluate the practical applicability of the network, followed by a sensitivity analysis to get an understanding about the influencing parameters in the network. The devised neural network, when implemented in a constitutive model displays a behavior similar to the theoretical constitutive model of Mohr-Coulomb soil. Finally, this thesis concludes with the implementation of the neural network constitutive model in Finite Element Analysis where it leads to convergence errors in the model.

A new approach is proposed which uses concepts of both the limit equilibrium method (LEM) and finite element method (FEM) to analyse slope stability. The main limitation of the LEM is that it can give an incorrect display of the initial stresses. The FEM is able to overcome this limitation, however, it uses the shear strength reduction method (SSRM) to determine the safety factor. This SSRM has its limitations as well. It’s a mathematical approach to obtain the factor of safety, without using the strength definitions given by advanced soil models. The proposed method is able to overcome the limitations of both LEM and FEM. It uses a genetic algorithm to search for the most critical slip surface and the safety factor is defined as the as the ratio of the total maximum available shear strength along the slip surface to the total mobilised shear stress. Multiple case studies are performed in this research and the genetic algorithm is able to find a reasonable critical slip surface and safety factor in every example. The proposed method has proven to give reliable results compared to other analysis methods: LEM, FEM and other optimisation methods proposed in literature.

Dikes protect people and lands all around the globe. With rising sea levels, the importance of well-designed dikes has never been more essential. One of the main factors that can compromise the stability of dikes are macro-instabilities. Macro-instabilities can cause dikes to lose their water bearing potential, leading to floods. Methods are developed to as accurately as possible determine the probability that a macro-instability takes place in order to prevent it. Therefore, recent advancements include the remaining strength after macro-instability, which may be able to prevent flooding. The foremost methods to determine this remaining strength after macro-instability use D-Stability or the Material Point Method (MPM). Both D-Stability and MPM have advantages and disadvantages. D-Stability can quickly determine the probability that a macro-instability takes place, but cannot model the process of failure, and must therefore simplify this process to estimate the remaining strength. MPM on the other hand can accurately model what happens after a macro-instability, but has a much larger computational cost, especially for probabilistic computations. By supplementing D-Stability and MPM, this thesis proposes a method that exploits the advantages of both methods and mitigates the disadvantages. For a dike with a single clay layer, a connection was made between D-Stability and MPM, allowing dike profiles to be transferred back and forth. In order to make this connection, the SHANSEP undrained shear strength model was successfully implemented in MPM. By using the quick probabilistic D-Stability calculation and the post-failure modeling option of MPM, the probability of failure and the effect of failure can be quickly determined. By taking into account the effect of failure, the method can determine the probability of flooding, without simplifying the failure process. The method can also be used to determine if flooding via retrogressive failure or a larger single instability is more likely. The method was tested via a case study and verified via a RMPM Monte Carlo analysis. For the case study, the probability of flooding was 5.189*10-3 compared to an initial probability of failure of 7.22*10-1, a reduction in the order of 139. This probability of flooding compared well to the probability of flooding of 5.308*10-3computed using the more accurate and computationally expensive RMPM. Based on these first results, the proposed method is a viable method to assess the probability of flooding after macro-instability for clay dikes.

The transition to renewable energy has lead to many new ways of energy production, such as wind and solar energy. For solar energy power, there is an excess of energy during summer and an deficiency of it during winter, which is reversed for the energy consumption. In order to amend this problem, there is an introduction of new ways of energy storage in the form of heat for domestic and industrial use. One of these is the use of Underground Thermal Energy Storage (UTES). These systems make use of the underground to store heat during summer and extract heat during winter by using a fluid as the heat transferring agent. There are many different ways of UTES installation, however there are two which are widely used in the Netherlands: 1) Aquifer Thermal Energy Storage (open systems) which are based on heat-transfer by convection and 2) Borehole Heat Exchangers (closed systems) which are based on heat-transfer by conduction. While the way of transferring heat differs, they both need vertically drilled boreholes to reach a certain storage depth. These boreholes are drilled close to, or underneath, pile-founded buildings. The drilling of boreholes influence in-situ soil stresses around them, which might influence the bearing capacity of the piles. The aim of this thesis is to investigate the influence of installing UTES systems on the bearing capacity of piles. This is done by: 1) investigating the stress changes in sand due to the drilling of boreholes (stress analysis) and 2) the influence of these stress changes on the bearing capacity of piles (bearing capacity analysis). All of the modelling was done in PLAXIS, by using the 2D-axisymmetric model for the first part and the 3D plane strain model for the second part. An advanced soil model for sand was used, with soil mechanical behavior according to the Hardening Soil small strain model. The stress analysis consists out of a parameter sensitivity analysis, including the key parameters: back-fill grout shrinkage and expansion, Over-Consolidation Ratio and Relative Density. After this, the influence of the borehole fluid pressure and borehole diameter on the stress states were investigated. As not all boreholes are drilled perfectly, the last part of the stress-analysis takes into account drilling complications. The bearing capacity analysis includes simulations of static pile load tests of non-displacement piles. Load-displacements curves are produced for several scenario's: 1) perfectly drilled boreholes, 2) boreholes with drilling complications, 3) influence of a soft soil layer and 4) displacement piles and the influence of a varying installation phasing. The results indicate no large effects of perfectly drilled boreholes on bearing capacity of piles located outside a zone of 1.5 times the borehole diameter. For boreholes with drilling complications, the bearing capacity is influenced for piles located in a zone of influence of only several meters (2-6m), depending on the severity of the complication. At last, loaded piles of existing buildings are extremely sensitive for stress changes due to the drilling process, which will cause severe loss of bearing capacity.

The stability of breakwaters in seismically active areas is not always guaranteed. A new approach to model a breakwater subjected to an earthquake is with PM4Sand. The goal of this research is to find if breakwaters subjected to earthquakes can be correctly modelled with PM4Sand.
To investigate PM4Sand a breakwater subjected to an earthquake on centrifuge scale is modelled in Plaxis. The soils of the centrifuge test are modelled with PM4Sand and UBCSand. After calibrating the soil parameters and incorporating the proper earthquake signal results are generated and compared. The investigated results focus on the settlements of the caisson, deformations of the breakwater and generated Excess Pore Water Pressures underneath the breakwater due to the earthquake. Comparing the results from the numerical models with the centrifuge test result show that both the UBCSand model and PM4Sand can give comparable results for the settlements and deformations. However, both UBCSand and PM4Sand were not able to give the correct EPWP development underneath the breakwater.
Due to the incorrect behaviour of the EPWP underneath the breakwater resulting from the numerical models, this research is not able to conclude that PM4Sand can be used for modelling breakwaters subjected to earthquakes. Further research is needed to investigate the development of EPWP underneath the breakwater during an earthquake. Focus points of future research can be: the influence of the amplitude of an earthquake signal on the EPWP, the influence of modelling a centrifuge test on the behaviour of the EPWP and the influence of the initial static shear stress on the EPWP development.

Geothermal energy can be a great solution for the downscaling fossil fuel society, but it can potentially lead to seismic hazards. A doublet system, with a cold water injection well and a hot water production well, alters the stress situation in the subsurface, which can result in (micro)fracturing and fault reactivation. Even in water filled reservoirs, aquifers, with relatively good permeabilities, the acting in-situ stress on already existing fault can be changed such that there can be a seismic hazard. The three dominant phenomena that influence the fault reactivation and are triggered by geothermal water injection and production are the direct pore pressure change, poro-elastic stress change and thermo-elastic stress change. To predict and subsequently diminish or limit the seismic hazards in geothermal operations, Seismic Risk Analysis (SRA) are to be completed before such operations can take place in an often seismic risky location, like densely populated areas. In the current (Dutch) geothermal environment mainly three SRA’s are used; “Methodiek voor risicoanalyse ontrent geinduceerde beving door gaswinning” by the Staatstoezicht op de Mijnen (SodM), “Defining the Framework for Seismic Hazard Assessment in Geothermal Projects V0.1” by Q-con/IF-technology [6] and an Excel-model created by TNO/Geomech. By investigating and reviewing these three SRA’s in this thesis their shortcoming and limitations are exposed, for example their lack of physical foundation and explanatory results. From the foundation of the currently excising SRA’s a new alternative SRA, which corresponds in some steps with the older SRA’s, is created in this thesis. In order to successfully finish the new SRA one of the three steps should be completed, starting with SRA Step 1. In this first step of the new SRA a new Physical Screening Model (PSM) is created. When completing the SRA an indication of what type of seismic monitoring there should be done during production. This PSM is a fairly quick and simple in its use but provides sufficient informative data to investigate the seismic hazard for most geothermal operations in the Netherlands. In four different steps in the PSM, the potential reactivation of faults over the whole reservoir during production will be evaluated. With the spatio-temporal evolution of ΔP, Δσporo, ΔT (PSM Step 1) this model can predict fault reactivation at any place and time inside the reservoir, while it can also look at which parameter dominated this reactivation. In this thesis the physical background and results of the PSM will be explained step by step. Eventually there are three final results from the PSM; a Mohr plot that predicts if certain faults (at certain locations) are stable or not and the maximum Moment magnitude (Mw) in combinations with the Peak Ground Velocity (PGV), which predict the severity of a possible event. Sensitivity analyses and case studies done with the PSM in this thesis show the influence of dominating parameters, like permeability and injection rate, and what results can be expected when using this model.