Circular Image

P.J. Vardon

info

Please Note

54 records found

Accurate estimation of the minimum horizontal stress is essential for well design, drilling safety, and predicting induced seismicity. In many regions of the Netherlands, however, only formation integrity tests (FIT) and leak-off tests (LOT) are available, raising the question of whether these drilling tests can be used to reliably infer the minimum horizontal stress. This study evaluates the reliability of stress estimation from hydraulic fracturing data, with a particular focus on the differences between fracture initiation and closure pressures.A comprehensive dataset of extended leak-off tests (XLOT), micro-fracture tests, and conventional LOT/FIT data is analysed. Multiple closure interpretation methods are compared, and their variability, bias, and applicability are quantified. In addition, initiation-based pressures are evaluated within a mechanical framework, including sensitivity analyses and stress consistency checks.The results show that closure-derived pressures provide the most reliable estimate of the minimum horizontal stress, whereas initiation-based pressures are strongly influenced by near-wellbore effects and operational conditions. Closure interpretation is inherently method-dependent, with inter-method differences up to 25~bar, significantly exceeding the repeatability within a single method. Among the evaluated techniques, the semilogarithmic derivative method demonstrates the best balance between robustness and applicability.Pressures from FIT and LOT tests often exceed closure-derived minimum horizontal stress and exhibit large variability, including non-physical results in a significant fraction of cases. While these pressures can be interpreted as upper or lower bounds within a mechanical framework, their quantitative reliability is limited. Furthermore, an offset data case study from nearby wells did not provide a reliable predictor the minimum horizontal stress.This study establishes a hierarchical interpretation framework in which closure-derived pressures form the primary estimate of the minimum horizontal stress, while initiation-based pressures provide only supplementary constraints.
...
Doctoral thesis (2026) - P. Kolah Kaj, H.A. Abels, A. Barnhoorn, P.J. Vardon
In the Netherlands, geothermal energy is considered a major contributor to achieving climate and energy goals. The success of geothermal projects depends strongly on understanding the reservoir. Knowing the thermo physical and mechanical properties of reservoir rocks, which govern heat transfer and mechanical stability, is therefore essential. However, these properties are often poorly constrained due to high measurement costs, the substantial time required for laboratory testing, and limited availability of suitable rock material. In addition, reliable rock property prediction is challenging because of strong heterogeneity in lithology, mineralogical composition, and diagenetic history.
This thesis generates a comprehensive database of thermal, acoustic, and mechanical properties for key Dutch geothermal formations. Based on measured data and their integration with downhole petrophysical logs, several predictive equations and models were developed, including machine learning approaches. These models improve property prediction tailored to the Dutch subsurface and enhance geothermal reservoir characterisation in general.
The research begins with a comprehensive study of Permian Rotliegend sandstones, a key geothermal reservoir in the Netherlands. More than 1100 core plugs were analysed to determine porosity, density, acoustic velocities, thermal properties, and mineralogy. The results confirm that porosity is the primary control on most rock properties. Higher porosity corresponds to lower density, acoustic velocity, thermal conductivity, and diffusivity. Systematic deviations from porosity trends were linked to mineralogical and diagenetic factors. For example, nacrite and other kaolinite group minerals enhanced thermal conductivity beyond porosity based predictions, whereas other clay types reduced it. Porosity dominates, but mineralogy and texture impose measurable secondary effects.
The analysis was extended to the Triassic Main Buntsandstein Subgroup, with more than 700 core plugs studied and compared directly to the Rotliegend dataset. Similar porosity dependent trends were observed, but systematic inter formation differences emerged. At equal porosity, Buntsandstein samples show lower thermal conductivity than Rotliegend samples. This difference is attributed to variations in clay type and distribution, as well as mineralogical features such as dolomite cementation and replacive clays. The lower Cretaceous Delft Sandstone Member was investigated to assess coupled mechanical and thermal behaviour. Laboratory tests included ultrasonic velocity measurements, thermal properties, and mechanical loading. Dynamic elastic moduli derived from ultrasonic data were systematically higher than static moduli measured during loading. A lithology specific workflow was developed to convert dynamic to static Young modulus, enabling continuous static modulus logs. Sandstones follow trends comparable to Permian samples, while clay rich intervals exhibit distinct but explainable behaviour due to low porosity.
The final part focuses on machine learning based prediction of thermal properties using laboratory and well log data. Ensemble models and regularised regression achieved promising results for thermal conductivity prediction, even in unseen wells. Thermal diffusivity remained poorly predictable, reflecting its sensitivity to mineralogical and microstructural factors. Density and acoustic features dominate conductivity prediction, whereas no single parameter controls diffusivity.
Overall, this thesis establishes a coherent framework for predicting thermo physical and mechanical properties of Dutch geothermal sandstones. It combines laboratory measurements, petrophysical analysis, and machine learning to improve reservoir characterisation and support reliable geothermal resource assessment.
...

From static to dynamic processes by means of Material Point Method analyses

Doctoral thesis (2026) - M. Bolognin, S.N. Jonkman, P.J. Vardon, M.A. Hicks
This dissertation focuses on validating the Material Point Method (MPM) for analyzing underwater flow slides and slope stability issues to enhance understanding and predictive capabilities. By addressing the complexities of multi-phase, multi-physics, and multi-scale problems, the research aims to establish a reliable numerical solution for simulating these phenomena, which are critical for flood defense systems. The research bridges gaps in knowledge by integrating experimental data, advanced numerical modeling, and validation studies. ...
Doctoral thesis (2026) - W. Luo, P.J. Vardon, Florian Amann, G.G. Drijkoningen, A.A.M. Dieudonné
Near-borehole coupled thermo-hydro-mechanical (THM) processes are critical in defining the success or failure of geothermal projects. Injectivity decline due to clogging or the increased viscosity of cold water reduces the safe operational window of geothermal projects. However, near-borehole fracturing, a result of near-borehole coupled thermo-hydro-mechanical (THM) processes, is thought to possibly be able to contribute to maintaining or improving the re-injection performance. Understanding and being able to simulate the fracturing processes is therefore of great significance for both analysing re-injection and designing stimulation operations. This thesis presents a geo-mechanical tool, based on the finite element method and cohesive zone model, to simulate the fracture initiation and propagation under THM loadings. The geo-mechanical tool is used to study several stimulation scenarios, including monotonic, stepwise, cyclic, and stepwise combined with cyclic stimulation. To monitor such stimulation processes, this thesis proposes a dual-cable distributed acoustic sensing (DAS) system in a single vertical well and investigates its feasibility to localise and understand near-borehole cracking events.

In the numerical method, possible discontinuities are represented by zero-thickness triple-nodded interface elements, which allow solid elements to separate with mechanical damage and the simulation of longitudinal and transversal fluid/heat flow in the discontinuity. The cubic law is used to simulate the fracture transmissivity changes, while an elasto-damage law is used to characterise the mechanical response of the discontinuity. To simulate the fracture initiation and propagation from high-permeability intact rock, interface elements are inserted in-between all the solid elements, with high stiffness and transversal hydraulic coefficient assigned to reduce artificial compliance. An artificial heat conductivity is introduced to stabilise the numerical solution, in which high Peclet numbers lead to numerical divergence. Substantial verifications and validation are implemented to demonstrate the performance of the developed method.

A new elasto-damage law is developed by incorporating a fatigue damage variable into the tensile branch, in order to account for the fatigue effects during the simulation of cyclic (thermal) stimulation to geothermal reservoirs. The fatigue damage variable is calibrated using the number of loading cycles and fatigue life at different load intensities, with Palmgren-Miner’s rule used to account for varying-amplitude cyclic loading. The proposed model is validated against extensive laboratory tests, including cyclic Brazilian test, cyclic hydraulic fracturing test and cyclic thermo-hydraulic fracturing test. The validation results show good agreement with the experimental data, demonstrating that the proposed model is capable of handling fatigue damage under cyclic and coupled THM loadings.

The developed tool is then used to study stimulation to a synthetic sedimentary reservoir, which, according to regional experience, is assumed to be clogged in the near-borehole region. THM simulations of various stimulation strategies - monotonic, stepwise, cyclic, and stepwise combined with cyclic - demonstrate that the stepwise stimulation yields the most favourable outcomes. Specifically, it enables a significantly lower peak injection pressure with more near-borehole damage. This performance is not achievable using either monotonic or cyclic strategies (assuming same Qinj and Tinj). Conversely, cyclic-injection-rate stimulation slightly underperforms (under high injection rate) or slightly outperforms (under low injection rate) the monotonic stimulation. A combined approach incorporating both cyclic and stepwise strategies may lead to slightly better stimulation performance, showing lower peak pressure, compared to corresponding monotnic stimulation, but is inferior to the stepwise stimulation alone.

The feasibility of using a single-well dual-cable DAS to fully localise and understand the near-borehole micro-seismic events is investigated based on synthetic signals, assuming homogenous and isotropic media. A localisation method is introduced to determine the source depth, epicentral distance and azimuth. Sensitivity analysis shows that the localisation accuracy is not sensitive to source with frequency varying from 50 Hz to 200 Hz. But a low signal-to-noise ratio and/or source-to-receiver azimuth close to 0◦ can lead to decreasing accuracy. Moreover, resolvability analysis suggest that non double-couple moment tensor components Mxx,Myy and Mzz can be reliably resolved with an epicentral distance within 20 meters, showing improvement on the case of only one cable in a well. A discussion based on the geo-mechanical simulation demonstrates that the single-well dual-cable DAS can be used to understand near-borehole tensile fractures induced during thermal stimulation, with a limited epicentral distance, which implies it is well suited to monitoring stimulation operations.

This thesis contributes to the energy transition by developing a geo-mechanical model to simulate cyclic and coupled THM processes, including the development of fractures, around the near field of the wellbore which can allow the design of novel cyclic thermal stimulation and by proposing a single-well dual-cable DAS configuration that is demonstrated to be feasible to localise and understand near-borehole micro-seismic events to monitor thermal stimulation operations.
...

How can the efficiency and reliability of energy pile system be enhanced ?

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. ...
Heating and cooling accounts for nearly half of the total energy consumption in the European Union,making geothermal systems a vital component for the energy transition. However, their operational lifetime is limited. Injected cold water forms a cold thermal front that spreads outwards. When this front reaches the production well, thermal breakthrough occurs and production temperatures drop. This report investigates whether seasonally reversing the flow direction by injecting hot water during periods of low heat demand can extend the reservoir’s lifetime or allow for reduced well spacing. This is important because it can optimize the use of geothermal systems in the future.To do this, annual energy demand was studied to find when heat demand is low. The available solar power in this period was then calculated for three different scenario’s to determine what reversed flow rate can be delivered. After establishing an analytical model and baseline simulation of the Delft geothermal reservoir, the three scenario’s were tested and the new breakthrough times computed. The effects of reservoir thickness and well-spacing were then studied.It was found that with the solar power currently available on the TU Delft campus, the lifetime of the geothermal reservoir can be extended by 7 years from the original 21 years. With a limited expansion of solar power on campus, this extension could be 11 years and if the system is run at maximum capacity,21 years is possible. It was found that if the lifetime were to be kept constant, the well spacing could be reduced from the original 1100 metres down to between 934.8 and 789.2 metres depending on the scenario.It was concluded that seasonally reversing the flow direction offers several opportunities to optimize geothermal systems. Either by extending the lifetime of the geothermal system or by allowing closer wel spacing to minimize costs whilst keeping a reasonable reservoir lifetime. ...

Integrating Renewable Sources with Thermal Storage for a Balanced Grid

As Amsterdam moves toward its 2040 carbon-neutral and 2050 renewable-heating targets, infrastructure must increasingly support the city’s energy transition. The Oostbrug, a planned cycling and pedestrian bridge across the IJ River, provides an opportunity to demonstrate an energy-autonomous bridge concept. Its combined electrical and thermal loads for lighting, operation, and winter de-icing are traditionally supplied by grid electricity and salt spreading. Integrating renewable energy sources with hydronic bridge-deck heating can reduce grid dependence, eliminate salt use, and enhance the bridge’s sustainability.
This thesis investigates the feasibility of a self-sufficient bridge energy system combining solar, wind, and thermal sources with on-site battery storage. Five configurations were modeled: three direct heating systems using energy piles or aquathermal heat pumps, and two indirect systems with Aquifer Thermal Energy Storage (ATES). A MATLAB-based hourly simulation using KNMI 2024 weather data and modeled IJ-water temperatures evaluated energy and peak coverage, embodied CO2, and battery requirements.
Among the direct systems, only the aquathermal heat-pump configuration (System 3) meets both annual and peak de-icing demands, though it requires substantial battery capacity and is sensitive to cold conditions. ATES-based systems (Systems 4 and 5) achieve full annual coverage with strong winter resilience. System 4, combining aquathermal recharge with ATES, provides the highest energy surplus of 2373 MWh yr−1 but with greater complexity and embodied emissions, while System 5, using bridgedeck regeneration, offers a simpler, more material-efficient solution with 1048 MWh yr−1 surplus and higher robustness. Both storage-led systems eliminate salt-based de-icing, avoiding 7.9 t CO2 yr−1 , salinization in the IJ river, and deterioration of the bridge construction. System 5 is identified as the preferred option for full-deck, energy-autonomous operation. The Oostbrug demonstrates how bridges can serve as energy-positive infrastructure, integrating structural, thermal, and electrical systems. Future work should extend model validation across multiple climate years and explore integration with Amsterdam’s district heating network. ...
Doctoral thesis (2025) - Josselin Ouf, Florian Amann, P.J. Vardon
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.... ...
Doctoral thesis (2024) - M. Mohsan, P.J. Vardon, Femke (F. C.) Vossepoel
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.
...

The Effectiveness of Generative Adversarial Networks in Subsoil Schematization

Master thesis (2023) - F.A. Campos Montero, P.J. Vardon, R. Taormina, B.E. Zuada Coelho
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. ...

Investigating the influence of soil heterogeneity on slope failure processes

Doctoral thesis (2023) - G. Remmerswaal, M.A. Hicks, P.J. Vardon
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. ...

Reservoir characterization, dynamic simulation and legislation in a medium deep aquifer

Master thesis (2022) - F.J.J. Smit, P.J. Vardon, D.F. Bruhn, Hans Veldkamp, Lies Peters
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. ...
Master thesis (2022) - B.T.M. van Esser, P.J. Vardon, E.A.C. Neeft, A.A.M. Dieudonné, H.A. Abels
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. ...
Master thesis (2022) - L.A. Kamphuis, J.P. Aguilar Lopez, Raymond van der Meij, P.J. Vardon, A.P. van den Eijnden
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. ...
Master thesis (2022) - L. Vrielink, P.J. Vardon, A. Daniilidis, M. Murali
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.
...
Master thesis (2022) - M. Novell Morell, P.J. Vardon, Siefko Slob, J.E.A. Storms, W.S.J. Uijttewaal, Marcela Busnelli
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. ...
Master thesis (2021) - C.K.L. Man, P.J. Vardon, G. Remmerswaal, K.A. Hildebrandt, M.J. Billeter, Patricia Ammerlaan
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.
...
Master thesis (2021) - C.P.B. de Beus, D.V. Voskov, X. Lyu, H. Hajibeygi, M. O. Saar, P.J. Vardon
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. ...