To assess the present and future conditions of potential deep geological repository sites, understanding their evolution in the past is mandatory. Here, glaciation cycles strongly affected the long-term thermo-hydro-mechanical (THM) evolution of the geosystem. The AREHS project studied the effects of changing boundary conditions (BCs) on the long-term evolution of large-scale hydrogeological systems. The focus is on numerical modeling the far-field using the open-source multi-field finite element code OpenGeoSys with full THM coupling. The impact of the glacial THM loading is taken into account using complex time-dependent THM BCs. In the present study, a generic geological model for a clay host rock formation including predominantly sedimentary rock layers is applied. The elasto-plastic behavior of all the layers is described with the Modified Cam clay model. Thus, a range of relevant effects (dilatancy, contractancy, consolidation etc.) can be considered with few material parameters. Special emphasis is put on the specification of a suitable initial state: To this end, an initial simulation is carried out, where a reasonable plastic pre-consolidation is adjusted. Then, the thermodynamic state is transferred in full to the subsequent 2D simulation of two glacial cycles. As a main result, the glacial cycles lead to persistent deviations in the subsurface, e.g. long-term pressure anomalies. Large glacial over- and underpressure cells disturb the in-situ hydraulic gradients and alter the flow fields around the DGR. As the host rock horizons represent aquitard layers they prevent periglacial circulation flow. No dilatant deformation was observed in this study. Contractant plastic deformation in combination with the HM coupling plays an important role as it significantly increases pore pressure peaks during glacial transit.

This paper aims to investigate the wave-induced evolution of small-strain stiffness and its effects on seismic wave propagation. To this end, an advanced numerical framework based on the dynamic porous media theory was developed, in which the Iwan multi-surface constitutive model was adopted to model the soil behavior during cyclic loading. Moreover, the numerical framework integrates key parameters such as ocean wave characteristics and depth-dependence seabed conditions to model the intricate interactions between waves and the seabed. Following model verification via analytical solutions and previous experimental data, comprehensive parameter studies are conducted, from which the effects of different wave conditions and seabed properties on the dynamic response of the seabed were obtained, revealing the wave-induced small-strain stiffness spatial and temporal variation. Subsequently, simulations of geophysical monitoring instants are conducted, assessing the impact of evolving small-strain stiffness on seismic wave propagation. The findings highlight the implications of stiffness changes on seismic wave propagation characteristics. The study provides valuable insights into the challenges and opportunities associated with interpreting geophysical data in dynamic submarine environments, offering implications for subsurface characterization and monitoring applications.

This paper is intended to serve as a low-hurdle introduction to nonlocality for graduate students and researchers with an engineering mechanics or physics background who did not have a formal introduction to the underlying mathematical basis. We depart from simple examples motivated by structural mechanics to form a physical intuition and demonstrate nonlocality using concepts familiar to most engineers. We then show how concepts of nonlocality are at the core of one of the most active current research fields in applied mechanics, namely, in phase-field modeling of fracture. From a mathematical perspective, these developments rest on the concept of convolution in both its discrete and its continuous forms. The previous mechanical examples may thus serve as an intuitive explanation of what convolution implies from a physical perspective. In the supplementary material we highlight a broader range of applications of the concepts of nonlocality and convolution in other branches of science and engineering by generalizing from the examples explained in detail in the main body of the article.

Deep geological repositories (DGRs) are designed to isolate radioactive waste (RW) from the biosphere over extremely long-time scales (i.e., several hundred thousand years). In order to assess the robustness of a safety case for a DGR, it is therefore necessary to analyse probable, less probable and hypothetical future developments due to, for example, climate change. Climate models predict that several ice sheets will advance and retreat over the next several thousands years. Glacial isostatic adjustment (GIA), resulting from the large moving mass of an ice sheet, can alter the displacement and far-field stress field of a DGR.Due to their extremely low matrix permeability, crystalline rocks are suitable host rocks for the disposal of RW in DGRs. However, the mechanical properties of crystalline rocks often promote crack growth and faulting, which in turn compromise their barrier function with respect to groundwater flow. In the INFRA project (DFG NA1528/2-1 and MA4450/5-1), we quantify how faults prone to reactivation during glacial events can affect radionuclide migration around a DGR in crystalline rock.We apply boundary conditions derived from an established GIA model [1,2] to a finite element model [3,4] of coupled fluid flow and radionuclide transport to numerically solve the component transport problem before and after fault reactivation. The Coulomb failure stress criterion is used as an indicator of permeability changes. The simulations show that GIA can increase permeability in the upper 400m of the reactivated faults. There, groundwater flow enhances radionuclide migration along the fault. In contrast, groundwater flow is reduced in the direction perpendicular to the fault plane. Although the proposed numerical workflow has been applied to the case of GIA, it can also be adapted to study hydromechanical processes induced by seismic events or by hydrofracking in enhanced geothermal systems. [1] Argus, D.F., Peltier, W.R., Drummond, R., Moore, A.W.: The Antarctica component of postglacial rebound model ICE-6G C (VM5a) based on GPS positioning, exposure age dating of ice thicknesses, and relative sea level histories. Geophysical Journal International 198(1), 537–563 (2014)https://doi.org/10.1093/gji/ggu140[2] Peltier, W.R., Argus, D.F., Drummond, R.: Space geodesy constrains ice age terminal deglaciation: The global model. Journal of Geophysical Research: Solid Earth 120(1), 450–487 (2015) https://doi.org/10.1002/2014JB011176[3] Bilke, L., Flemisch, B., Kalbacher, T., Kolditz, O., Helmig, R., Nagel, T.: Development of Open-Source Porous Media Simulators: Principles and Experiences. Transport in Porous Media 130(1), 337–361 (2019)https://doi.org/10.1007/s11242-019-01310-1[4] Bilke, L., Fischer, T., Naumov, D., Lehmann, C., Wang, W., Lu, R., Meng, B., Rink, K., Grunwald, N., Buchwald, J., Silbermann, C., Habel, R., Günther, L., Mollaali, M., Meisel, T., Randow, J., Einspänner, S., Shao, H., Kurgyis, K., Kolditz, O., Garibay, J.: OpenGeoSys. Zenodo (2022)https://doi.org/10.5281/zenodo.709267

To properly assess the present and future conditions of potential nuclear waste repository sites, understanding their evolution in the past is mandatory. Here, glaciation cycles strongly affect the long-term thermo-hydro-mechanical (THM) evolution of the geosystem. The AREHS project studies the effects of time-dependent boundary conditions on the long-term evolution of large-scale hydrogeological systems. The focus is on numerical modeling using the open-source multi-field finite element code OpenGeoSys with THM couplings. The impact of the glacial THM loading is taken into account using appropriate time-dependent THM boundary conditions. The generic geological model for a clay host rock formation includes almost only sedimentary rock layers. Within the scales of the physical problem, it can be assumed that the plastic flow behavior of the different sedimentary rocks shares qualitative features. Therefore, the same generic material model is used for all layers: The elasto-plastic modified Cam clay (MCC) model can describe qualitatively a range of relevant effects (dilatancy, contractancy, consolidation, normal and overconsolidation effects, etc.) with a manageable number of material parameters. The chosen approach is considered sufficient to demonstrate the effect of inelastic deformations. Comparing the results to the elastic reference case, the role of plastic flow in the context of such a long-term simulation is elaborated. Special emphasis is put on the specification of a suitable initial state: To this end, an initial simulation is carried out, where plastic flow can occur under the gravitational load. The fields for temperature, pore water pressure and effective stress are then transferred as initial values for the glacial cycle. Additional internal variables of the material model are automatically adopted. This means that the thermodynamic state is transferred in full. The simulation results are analyzed with respect to potential safety-critical parameters, such as maximum temperature, hydraulic pressure, subsidence, equivalent effective stress and strain.

To assess the robustness of a safety case for a deep geological repository (DGR), it is necessary to analyze a range of scenarios covering likely, less likely, and hypothetical future developments. Crystalline rock can, under ideal conditions, provide a suitable hydrogeologic barrier due to its extremely low matrix permeability. However, this host rock is often fractured, which can compromise its hydro-mechanical (HM) barrier function. We quantify how faults that are prone to reactivation during glacial events can affect radionuclide migration around a DGR in a crystalline host rock. We extend a previously developed finite element model of coupled fluid flow and radionuclide transport to numerically solve the component transport problem before and after fault reactivation. Assuming that fault reactivation is triggered by changes in mechanical boundary conditions, we derive heterogeneous permeability distributions in the reactivated faults by evaluating the Coulomb failure stress criterion of finite element solutions of a complementary hydro-mechanical problem. Specifically, we evaluate the consequences of glacial isostatic adjustment (GIA) during a glacial cycle. We find that the increased permeability in the reactivated faults accelerates the migration of radionuclides along the fault by channeling the flow, while it is reduced in the direction perpendicular to the fault. The channeling observed is also a result of heterogeneous permeability enhancement, and the flow fields differ from those of the previous model which postulated a homogeneous permeability enhancement. Although the proposed numerical workflow has been applied to the case of GIA, it is adaptable to study hydro-mechanical processes induced by seismic events or by hydrofracking in enhanced geothermal systems.

This paper is intended to serve as a low-hurdle introduction to non-locality for graduate students and researchers with an engineering mechanics or physics background who did not have a formal introduction to the underlying mathematical basis. We depart from simple examples motivated by structural mechanics to form a physical intuition and demonstrate non-locality using concepts familiar to most engineers. We then show how concepts of non-locality are at the core of one of the moste active current research fields in applied mechanics, namely in phase-field modelling of fracture. From a mathematical perspective, these developments rest on the concept of convolution both in its discrete and in its continuous form. The previous mechanical examples may thus serve as an intuitive explanation of what convolution implies from a physical perspective. In the supplementary material we highlight a broader range of applications of the concepts of non-locality and convolution in other branches of science and engineering by generalizing from the examples explained in detail in the main body of the article.

A 3D numerical simulation of radionuclide migration from a potential repository of solid radioactive waste is carried out. The federal site is located in the Krasnoyarsk region (Siberia, Russia), at a distance of a few kilometers from the eastern bank of the river Yenisei. The models take into account the actual topography of the area and the heterogeneity of the hydraulic conductivity of the enclosing rocks. The migration of long-lived radionuclides is considered, e.g. 241 Am (half-life of 432.6 years). The simulated time intervals range from the waste disposal closure to the fault emergence. The model settings consider different distances between the fault and the repository. The highly permeable near-vertical faults run either from South to North (meridional) or from West to East (latitudinal). The simulation results show that the meridional faults do not exert a substantial influence on contaminant plume movement. By contrast, the latitudinal faults can represent a significant ecological hazard depending on the distance between the repository and the fault. While this study is not an assessment of the safety of the potential repository site, the results point out that further investigations of the stress state of the enclosing rocks are needed to better constrain potential faults and flow directions.

Extensive mine dumps consisting of loosely deposited sands have been created as a result of open-pit lignite mining, with a risk of soil liquefaction under high water saturation and a corresponding initiating event. Soil compaction is one of the feasible methods for reducing the probability of liquefaction. For the monitoring of liquefaction events and the evaluation of compaction work, seismic survey methods with sensitivity to changes in soil saturation and structure may thus complement other methods. Compared to exploration methods for deep systems, the shallow subsurface presents some unique challenges. To this end, an open-source, customizable code based on Biot's theory was developed in the FEniCS library, which takes into account partial saturation and porosity dependence of stiffness, permeability, and other quantities. Following code verification, a comprehensive investigation of parameter studies is conducted, from which the effects of different factors on wave propagation characteristics were obtained. The numerical model was applied to simulate the expected changes in seismic response following soil compaction. Furthermore, the position of the high saturation area could be detected from the reflection and refraction P waves. The goal of this work is to provide an analysis framework for the assessment of compaction works and monitoring liquefiable soils in mine dumps under conditions of variable saturation due to rising groundwater tables.

<p>Important aspects for subsurface installations, such as Deep Geological Repositories (DGRs), in crystalline rock are the presence and evolution of fractures and faults, since they control the subsurface flow regime. According to climate extrapolation, it is expected that cold and warm period will alternate, accompanied by ice sheet progression and regression. The large moving mass of an ice sheet causes a dynamic response of the earth's crust, referred to as glacial isostatic adjustment (GIA) [1]. GIA changes the displacement and stress field not only under and near the ice sheet but also in its far-field. In view of the long-term assessments, we apply boundary conditions derived from an established GIA model [2] in order to analyze induced far-field stress and pore pressure changes and their impacts on existing faults in a hydromechanical simulation. As indicator for permeability changes we apply the Coulomb failure stress criteria. To quantify the consequences on the subsurface flow we run a component transport simulation before and another after the fault reactivation, revealing how the faults canalize the radionuclid propagation. For both kinds of simulations, hydromechanical and component transport, we apply Finite-Element methods (FEM) [3].<br />The INFRA project is funded by the DFG under grants NA1528/2-1 and MA4450/5-1.<br /> <br />[1] Holger Steffen, Patrick Wu. "Glacial isostatic adjustment in Fennoscandia - a review of data and modeling". Journal of geodynamics 52.3-4, p. 169-204, 2011. https://doi.org/10.1016/j.jog.2011.03.002<br />[2] Georg Kaufmann. "Program package ICEAGE". Manuscript, Institut für Geophysik der Universität Göttingen, vol. 40. p. 840, 2004.<br />[3] OpenGeoSys 6.4.3. Lars Bilke, Thomas Fischer, Dmitri Naumov, Christoph Lehmann, Wenqing Wang, Renchao Lu, Boyan Meng, Karsten Rink, Norbert Grunwald, Jörg Buchwald, Christian Silbermann, Robert Habel, Linda Günther, Mostafa Mollaali, Tobias Meisel, Jakob Randow, Sophia Einspänner, Haibing Shao, Kata Kurgyis, Olaf Kolditz, Jaime Garibay. 2022. https://doi.org/10.5281/zenodo.7092676</p>

We revisit the classical and solved problem of the terrestrial brachistochrone, the fastest path between two points in earth's gravitation field, by an approach we refer to as experimental numerics. By this term we mean arriving at a qualified guess by deliberately taking inspiration from numerical results that are easily available.Since in many cases verification is easier than derivation, this approach may have some educational merits.Current software tools such as Jupyter Notebooks blend coding with documentation and allow leveraging this approach to enable new ways in modern teaching.The intended audience are graduate students with prior knowledge of multivariate calculus, ordinary differential equations (ODEs), calculus of variations and classical physics, particularly mechanics.

One of the most important aspects for Deep Geological Repositories (DGRs) in crystalline rock is the presence and evolution of fractures and faults, since they dominate the subsurface flow regime and thus the possible transport of contaminants. In the considered period of one million years, it is expected that cold and warm period alternate, accompanied by ice sheet progression and regression. The large moving mass of an ice sheet causes a dynamic response of the earth's crust, referred to as glacial isostatic adjustment (GIA). GIA changes the displacement and stress field not only under and near the ice sheet but also in its far-field. In view of the long-term safety assessments for DGRs, we analyze GIA-induced far-field stress and pore pressure changes and their impacts on existing faults.For that purpose, we use Finite-Element methods (FEM) to simulate the hydromechanical processes around an exemplaric DGR of the Yeniseiskiy Site, Russia, applying boundary conditions derived from established GIA models [1,2]. As result, we obtain the Coulomb failure stress for varying instances of assumed faults.The INFRA project is funded by the DFG-RFBR program:DFG funds: NA1528/2-1 and MA4450/5-1RFBR funds: 20-55-12009, АААА-А20-120012190168-5References[1] Patrick Wu. “Using commercial finite element packages for the study of earth deformations, sea levels and the state of stress”. In: Geophysical Journal International 158.2 (2004), pp. 401–408.[2] G. Spada et al. “A benchmark study for glacial isostatic adjustment codes”. In: Geophysical Journal International 185.1 (2011), pp. 106–132.

While guidelines for the location and licensing of a deep geological repository (DGR) for high-level radioactive waste depend both on national government policies and international regulations, it is mandatory to select a site where the hydrogeological setting provides sufficiently safe natural conditions for long-term waste isolation from groundwater flow. Therefore, safety assessments of a suitable location of a DGR require the evaluation of future external events and processes that may affect its long-term evolution.Here, glaciation cycles are of special importance: Ice sheets evoke crustal deflections (including deformation), and impose pronounced hydraulic heads, both of which change the large-scale hydrogeological conditions. To properly assess the present and future conditions of a DGR site, its evolution in the past should be understood. For this, a sedimentary basin [3] is considered here as a large-scale hydrogeological benchmark. The evolution during one glacial cycle is simulated using the open-source multi-field finite element code OpenGeoSys. The hydraulic-mechanical impact of the glacial loading is taken into account using appropriate time-dependent boundary conditions. For comparison with a previously published study [3], the same (heuristic) displacement field is prescribed and the groundwater evolution is regarded. Then, a more realistic displacement field obtained from large-scale GIA simulations [1,2] is prescribed. Using a one-sided mechanical-hydraulic coupling with a staggered solution scheme it is possible to consider not only the hydraulic head from the glacier and the crustal deflection but also the crustal compression. Especially in regions at the margin of the glacier this is could have an impact on the hydraulic behavior at the depth of a DGR.References & FundingThis research is funded by the Federal Office for the Safety of Nuclear Waste Management under Grant No. 4719F10402 (AREHS project)[1] Patrick Wu. “Using commercial finite element packages for the study ofearth deformations, sea levels and the state of stress”. In: Geophysical Journal International 158.2 (2004), pp. 401–408.[2] G. Spada et al. “A benchmark study for glacial isostatic adjustment codes”. In: Geophysical Journal International 185.1 (2011), pp. 106–132.[3] V.F. Bense and M.A. Person. “Transient hydrodynamics within intercratonic sedimentary basins during glacial cycles”. In: Journal of Geophysical Research: Earth Surface 113.F4 (2008).

Recent applications require the extension of stability theory to axially loaded viscoelastic bars. For elastic bars it is known that a pulsating load leads to a Mathieu equation. Meaning, there are some instable regions before the static buckling load and some stable regions beyond, depending on the excitation (amplitude, frequency, offset). To fill this gap we follow along the lines of Weidenhammer and adopt the classical Euler‐Bernoulli beam kinematics for coupled longitudinal and bending vibrations. Modeling viscoelasticity by the Standard Linear Solid model introduces an internal variable. Hence, there are auxiliary transversal and longitudinal displacements corresponding to the internal variable in addition to the physical transversal and longitudinal displacements of the bar centerline. The modeling further assumes only longitudinal vibrations in the stable regime, thus there are no bending vibrations and the description simplifies to a system of linear partial differential equations. Further, a one‐sided coupling is assumed for the stability analysis, i.e. the longitudinal vibrations are prescribed and induce bending vibrations. Still, the time‐variant coefficients of the transversal dynamics make it difficult to find an analytical solution. Consequently, an approximation is obtained by Ritz method from Hamilton's principle. Applying the eigenforms of elastic buckling as trial functions finally leads to a system of rheolinear ordinary differential equations. Its stability chart with respect to static and harmonic load components is calculated numerically by Floquet theory.

From previous studies it is evident that decoupled simulations lack the ability to capture certain coupled effects, such as the Noordbergum effect or the Mandel-Cryer effect in a hydraulic-mechanical context. Thus, for detailed simulations of geotechnical or geological system, coupled simulations are usually chosen. For example, thermal-hydraulic-mechanical (THM) coupled systems, and even chemical and biological couplings (THMCB), are considered in simulations used to assess barrier integrity over long time spans in the context of geological waste disposal.This paper is restricted to coupled hydraulic-mechanical (HM) systems. A monolithic approach is both stable and accurate for strongly coupled systems. However, as site-scale models of geological disposal facilities are also large in spatial dimensions, it is worth to investigate how staggered methods may cut down the computational costs. The fixed-stress split appears to be a promising approach for staggered schemes in terms of stability, consistency, accuracy, and efficiency.While adding another iteration level in comparison to monolithic schemes, staggered schemes allow for lower-order approximation spaces, whereas monolithic schemes require Taylor-Hood elements resulting in a larger number of degrees of freedom per element. Both coupling schemes are implemented in the the open-source finite-element (FE) software OpenGeoSys and used to simulate a large-scale model, which is oriented towards a real site in planning in Russia. Simulation results are compared in terms of accuracy, coupling effects and performance.

Abstract. Under ideal conditions, owing to its extremely low matrix permeability, crystalline rock can constitute a suitable hydro-geological barrier. Mechanically, its high strength and stiffness provide advantages when constructing a repository and for long-term stability. However, crystalline rock usually occurs in a fractured form, which can drastically alter hydromechanical (HM) barrier functions due to increased permeability and decreased strength. Seismic events have the potential to alter these HM properties by activating faults, increasing their transmissibility, creating new fractures or altering network connectivity (De Rubeis et al., 2010). Therefore, it is of high importance to build computational models to allow assessment of the HM effects of seismic events in a Deep Geologic Repository (DGR) in crystalline rock, as illustrated in Fig. 1. For this purpose, we consider a DGR in Russia (Yeniseysky site) for high-level waste in crystalline rock (Proterozoic and Archaean gneiss complexes) that is located close to a potentially seismically active area (Jobmann, 2016). Here, we present a coupled HM simulation, using OpenGeoSys (Kolditz et al., 2012), of a large-scale, three-dimensional finite-element model of the Yeniseysky site to assess the consequences of seismically induced stress-field changes on the local stress field and the fluid flow. This research also provides an outlook of current model development geared towards a more detailed assessment of seismically induced hydro-mechanical processes in porous and fractured rocks.

In real machines and mechanisms, almost all joints are limited. A common "hard limiter" in mechanical systems is a collision between bodies, or between bodies and system boundaries. Often these collisions are modelled by impulse-like contact forces. During such a collision of two or more bodies, their positions go on continuously, but their velocities jump. We consider this non-smoothness and its consequences for simulation and control in a variational setting, from which a discretization follows in a systematic way. The resulting variational integrator (VI) for collisions is known from the literature, but the incorporation of the collision equations into the Discrete Mechanics and Optimal Control (DMOC) approach follows a new idea. As many technical systems are underactuated, we take this into account too. The feed-forward control for the swing-up of a pendulum-on-cart-system, a.k.a. inverted pendulum, proves the concept. The sequence of collisions is not optimized yet, but for a given sequence, we find the cart force to optimally steer the pendulum into its final state.

This tutorial is a ready-to-run LaTeX example that prospective authors of GAMMAS may substitute with their own content. Moreover, it contains information about some journal policies.

In terms of simulation and control holonomic constraints are well documented and thus termed standard. As non-standard constraints, we understand non-holonomic and unilateral constraints. We limit this survey to mechanical systems with a finite number of degrees of freedom. The long-term behavior of non- holonomic integrators as compared to structure-preserving integrators for holonomically constrained systems is briefly discussed. Some recent research regarding the treatment of unilaterally constrained systems by event- driven or time-stepping schemes for time integration and in the context of optimal control problems is outlined.

AbstractFor the special case of incompressible and highly viscous fluids, the interaction with a rigid body can be collected in a damping matrix, relating the velocities and angular velocities of the body with the fluid force and torque. This damping matrix (a.k.a. viscous resistance matrix) depends exclusively on the geometry and needs to be computed only once. We consider a rigid body moving in an unbound fluid. The generalized dissipative forces from the fluid onto the body enter the time discretization via the discrete D'Alembert Principle. As generally large rotations may occur, we chose quaternions for a singularity‐free description of the body orientation. The corresponding holonomic constraint of a unit quaternion is enforced on the position and momentum level by the RATTLE algorithm. The problem of Stokes drag on a sedimenting particle serves as an example.