Multiscale simulation frameworks are essential to quantify the CO2 trapping and migration in large-scale saline aquifers, which entail highly-resolved fine-scale heterogeneous properties. However, classical upscaling approaches which aim to define effective properties on larger grid sizes can lead to significant and systematic overestimation of the solubility and residual trapping mechanisms. Reliable assessment of these two trapping mechanisms is crucial to ensure the integrity of the storage process and properly mitigate the leakage risks. Therefore, it is essential to develop advanced simulation technologies that are both accurate and efficient (i.e., scalable) for simulation of complex CO2 plume dynamics within large-scale heterogeneous reservoir models. To overcome this challenge, in this work three advanced strategies are developed and investigated: Effective Values (EV) for parameters, Local Grid Refinement (LGR) and Algebraic Dynamic Multilevel (ADM). The numerical investigations specially include a set of consistent models in the Ponta Aguda saline aquifer, with a total area of 40,000 km2[jls-end-space/], located offshore the Brazilian coast. The results indicate that the ADM is a promising method, delivering stable and robust results in a representative section of the field. This encourages further extensions of this method for real-field deployment. Specially, LGR and EV are found to be limited in their scopes for field simulations, since they depend on a matching pre-procedure (against a reference solution) for their upscaled parameters before any new simulations can be run. In addition, their tuned parameters cannot be transferred from one model to another. ADM, on the other hand, does not require any upscaling procedure, as the multiscale basis functions allow for consistent mapping across resolutions.

Salt cavern simulations involve many numerical challenges that need to be addressed in order to ensure accurate and meaningful results. Firstly, lithological structures and solution-mined salt caverns always present fairly complex shapes, which favors the use of tetrahedral meshes with local refinements for adequate domain discretization. Secondly, salt rocks are known to creep under deviatoric stresses, meaning that deformations take place at constant volume (isochoric). The combination of isochoric deformations with tetrahedral meshes is particularly problematic for low-order finite element formulations. This work presents a stabilized mixed finite element (FE) formulation for linear tetrahedrons, where the mean stress is a primary variable, incorporating all the relevant deformation mechanisms for salt rocks. The stabilization consists of enriching the displacement FE approximation in the mean stress equation by obtaining an approximation for the Laplacian of the displacement that accounts for inelastic strains. This is achieved by using the Physical Influence Scheme (PIS) with the concept of secant Young’s modulus, which promotes local stabilizations where necessary. When combined with a proper calculation of a geometric parameter h[jls-end-space/], this stabilization technique is shown to produce oscillation-free and physically consistent results without any sort of tuning parameter. The proposed technique is analyzed in relevant test cases for salt cavern simulations and the results show the effectiveness of the proposed stabilization to eliminate spurious numerical oscillations with low-order tetrahedral meshes.

Underground hydrogen storage (UHS) is a potential technology that can resolve renewable energy supply-demand challenge at seasonal (terawatt-hours) scales. Enabling this technology and optimizing its performance require a wide range of analyses from hydrodynamics to geomechanics and biogeochemistry, among which understanding the transport (and trapping) of hydrogen in porous rocks stands out. A key parameter in quantification of hydrogen transport in partially brine-saturated geological formations is its relative permeability (k_rh). In this study, we develop a theoretical k_rh model using upscaling concepts from effective medium approximation and percolation theory. Our theoretical model, developed based on pore-scale characteristics, estimates k_rh from pore size distribution, capillary pressure curve or mercury intrusion capillary pressure curve, and critical hydrogen saturation, S_hc, at which k_rh approaches zero. We evaluate the proposed model using eight experimental datasets and eleven pore-network simulations. Discrepancies are observed for some of the carbonate samples, likely due to secondary porosity effects (e.g., presence of vugs and/or fractures), and in some of the sandstone rocks, possibly due to imprecise S_hc estimation. These observations highlight the importance of improving pore structure characterization to better account for such heterogeneities and enhance model accuracy for reliable quantification of the k_rh relevant to UHS applications. These findings also highlight the critical role of accurate parameter estimation, such as determining the S_hc in estimating k_rh. Overall, the study demonstrates that the proposed approach provides a cost-effective and practical alternative to extensive experiments and simulations, offering a promising tool for quantifying k_rh relevant to UHS applications.

This study extends the Algebraic Dynamic Multilevel (ADM) method for simulating contaminant transport and retention in vadose zones. Building upon a fully implicit scheme that couples variably saturated flow and contaminant transport, the developed ADM framework effectively predicts contaminant plume migration across both unsaturated and saturated media under heterogeneous conditions. During the simulation, ADM dynamically adjusts grid resolution based on the spatial gradients of primary variables, applying fine-scale grids in regions with steep gradients and coarsening the mesh where fields remain smooth. These dynamic adjustments are achieved through prolongation and restriction operators that transfer solutions across multilevel grid systems. As both water content and contaminant concentration evolve spatiotemporally, dual coarsening criteria are introduced to simultaneously capture flow and transport dynamics. Results show that the developed model reproduces the contaminant migration obtained from the fully resolved solution using substantially fewer grids. Moreover, it offers the flexibility to trade off numerical accuracy against computational cost by selecting an appropriate coarsening criterion.

This study introduces a multiscale simulation framework, termed Projection-based Embedded Discrete Fracture Modeling with Algebraic Dynamic Multilevel method (pEDFM-ADM), which integrates an embedded discrete fracture network representation with a fully algebraic, front-tracking-based mesh adaptation strategy. Incorporating a fully implicit scheme, compositional thermodynamics, and algebraic multilevel operators, the framework captures essential subsurface processes such as buoyancy-driven migration, convective dissolution, phase partitioning, and fracture-matrix interactions under geologically realistic conditions. The method constructs a hierarchy of multilevel grids and localized multiscale basis functions that introduce fine-scale heterogeneities at each coarse level. Adaptive mesh refinement and coarsening are driven by local variations in CO₂ mass fraction and executed through algebraic prolongation and restriction operators, enabling efficient projection between grid levels. The framework is systematically evaluated across a sequence of test cases with increasing complexity, including systems with low-permeability flow barriers, highly conductive fractures, striking a trade-off between computational resource and detailed simulation accuracy. Overall, the pEDFM-ADM framework provides a scalable, fully algebraic, and physically adaptive modeling tool for large-scale CO₂ storage simulations in fractured porous media, supporting predictive simulation and risk assessment for long-term carbon sequestration.

This work addresses numerical instabilities that can appear when computing the mean stress in linear elasticity and coupled poroelasticity problems discretized with low-order finite elements. The linear elasticity and coupled poroelasticity models are solved using both primal and mixed finite element formulations. Stabilization is obtained by enriching the finite element approximation with an approximation of the Laplacian of displacements. This Laplacian is then evaluated with the Physical Influence Scheme (PIS) by leveraging the underlying governing equation. A key step in the proposed stabilization is the calculation of a parameter h, often computed in the literature as a characteristic length of the element. In this work, we calculate h by solving an optimization problem at the element level. To avoid the high computational cost associated with this procedure, a machine learning model is proposed to predict the optimal h. The benefit of combining PIS with an appropriate computation of h is that the resulting stabilization scheme does not rely on any type of heuristic or user-specified tuning parameter, as often required in other stabilization methods. The results show that the proposed stabilization strategy can effectively remove both saddle-point and Gibbs mean stress oscillations in linear elasticity. We also report, for the first time, that mean stress oscillations can also appear when solving coupled poroelasticity problems, and, differently from pore pressure oscillations (which naturally vanish with time), mean stress instabilities are persistent throughout the whole simulation time, unless deliberately removed. The proposed stabilized mixed formulation is able to remove both pore pressure and mean stress oscillations in coupled poroelasticity problems. Finally, the calculation of h is shown to be critical for the quality of the stabilization, with the machine learning-based approach providing the best compromise between numerical diffusion and accuracy.

Large-scale geological storages of hydrogen (H₂) and carbon dioxide (CO₂) in saline aquifers present feasible options for a sustainable energy future. We compared the plume migration of CO₂ and H₂ in aquifers using the FluidFlower benchmark, incorporating the state-of-the-art thermophysical and petrophysical properties. The H₂ plume, with its higher buoyancy and mobility compared to CO₂, remains predominantly in the gas phase due to its lower solubility, increasing the chances of escaping through fractures or migration to distant regions. This additionally leads to a higher pressurized reservoir, which, along with higher buoyancy, increases the chance of caprock penetration. Dissolution trapping of CO₂ into brine increases over time due to its fingering, while H₂ does not show fingering. Our findings show that while geological carbon storage (GCS) benefits significantly from all structural, dissolution, and residual trapping, underground hydrogen storage (UHS) relies mainly on structural trapping, making the integrity of sealing elements of the system a key factor in its performance.

Residual trapping is a critical mechanism influencing the efficiency of Underground Hydrogen Storage (UHS). This study investigates the underlying processes of residual trapping by bypassing, through bifurcating geometries, focusing on how geometrical parameters and flow characteristics affect the trapping process. We develop a dynamic simulation framework based on the lattice Boltzmann method (LBM) to simulate full drainage/imbibition cycles. Various geometries, based on the pore doublet model, were investigated and supported by theoretical analysis. In addition, trapping behavior of hydrogen was compared to that of CO₂ and CH₄. It is found that the channel width ratio, specially across the local bifurcating geometries, and the roundness of the grains, are among the key factors which control hydrogen trapping. Results indicate that the suited reservoirs for underground hydrogen storage have narrower channel-size ratios and smoother edges at micro-scale. Operational conditions also play a significant role. Lower flow rates enhance bypassing, which increases trapping.

Finite-size effects of transport properties computed from molecular dynamics simulations are investigated for Weeks-Chandler-Andersen systems at reduced densities of 0.05 (dilute gas), 0.45 (dense gas), and 0.85 (fluid close to the solid-liquid transition). Viscosities, self-diffusivities, Onsager coefficients, and electrical conductivities are computed for various system sizes ranging from 64 to 8192 WCA particles at each density. At dilute and intermediate densities, finite-size corrections to the transport properties significantly deviate from the widely used Yeh–Hummer correction, which was originally developed for the liquid phase.

This work introduces a novel application of the Algebraic Dynamic Multilevel (ADM) method for simulating CO₂ storage in deep saline aquifers. By integrating a fully implicit coupling strategy, fully compositional thermodynamics, and adaptive mesh refinement, the ADM framework effectively models phenomena such as buoyancy-driven migration, convective dissolution, and phase partitioning under various subsurface conditions. The method starts with the construction of a hierarchy of multilevel grids and the generation of localized multiscale basis functions, which account for heterogeneities at each coarse level. During the simulation, the ADM method dynamically refines areas with significant overall CO₂ mass fraction gradients while coarsening smooth regions, thus optimizing computational resources without compromising the accuracy required to capture essential flow and transport characteristics. This dynamic grid adjustment is facilitated by algebraic prolongation and restriction operators, which map the fine-scale system onto a coarser grid suited to the evolving distribution of the CO₂ plume. This feature allows the ADM to navigate various coarsening thresholds efficiently, striking a trade-off between computational economy and detailed simulation accuracy. Even at relatively higher thresholds, key trapping mechanisms are captured with sufficient detail for quantification. These capabilities make the ADM framework well suited for long-term CO₂ sequestration in highly heterogeneous reservoirs, where large-scale models may otherwise become impractically expensive, offering a practical balance between the need for detailed simulations and manageable computational requirements. Overall, the ADM framework proves to be a robust tool for designing, monitoring, and analyzing large-scale CO₂ storage operations, supporting reliable and cost-effective solutions in carbon management.