Challenges in simulation of geological CO$_2$ sequestration and supercritical geothermal reservoirs

Büsing, Henrik; Clauser, Christoph (Thesis advisor); Bücker, H. Martin (Thesis advisor)

Aachen : RWTH Aachen University (2021)
Dissertation / PhD Thesis

Dissertation, Rheinisch-Westfälische Technische Hochschule Aachen, 2021, Kumulative Dissertation

Abstract

Multi-phase flow processes occur in many different scientific fields: During enhanced oil recovery, gas injection increases the recovery factor of the oil produced; in radioactive waste management, hydrogen is produced due to barrel corrosion; in groundwater management, non-aqueous phase liquids are remediated by injection of hot steam; in CO$_2$ geo-sequestration, CO$_2$ is injected into deep saline aquifers to mitigate its anthropogenic effects on climate; and during the production of geothermal energy, water and steam may occur in supercritical geothermal reservoirs. All these processes involve multiple phases, which, in addition, might react with one another. These reactions lead to the appearance or disappearance of different phases. Any multi-phase flow simulator needs to address these particular challenges. The objective of this thesis is to better understand the complexities associated with the numerical simulation of multi-phase, multi-component flow and provide efficient solution methods. Generally, fully coupled formulations are preferred for the simulation of multi-phase flow since they account for the strong coupling of the different equations and their highly non-linear behavior, while allowing for large time steps given an implicit time discretization. A first implementation of a fully coupled formulation for two-phase flow was performed using Matlab. This Matlab prototype used ADiMat to calculate exact derivatives via automatic differentiation (AD) for the Jacobian, avoiding the error-prone calculation of the Jacobian with finite difference approximations. AD derivatives, while being more accurate, also led to a decrease in the number of Newton iterations needed in each time step. Later, the Matlab prototype was implemented into SHEMAT-Suite using either constant fluid properties or pressure- and temperature-dependent properties of CO$_2$ and water. This implementation makes use of Tapenade for the calculation of the exact Jacobian and PETSc for the linear and non-linear solver parts. The efficiency of a Schur Complement Reduction method using algebraic multigrid (SCR-AMG) and a Constrained Pressure Residual method (CPR-AMG) was demonstrated on different test cases: advection- and diffusion-dominated heterogeneous two-dimensional test examples, the Sleipner reservoir and the SPE10 benchmark. CPR-AMG is particularly suited for the advection-dominated test case, while SCR-AMG is preferred for the diffusion-dominated case. The Sleipner test case also favors SCR-AMG. The SPE10 benchmark is the most challenging test case. Here, only simulations using AD finished in the desired simulation time of 72 hours. Finally, a weak scaling test on JUQUEEN, a BlueGene/Q supercomputer, was conducted showing efficiencies of over 95 % for both SCR-AMG and CPR-AMG. AD reduces the number of Newton iterations and linear iterations for the iterative solver compared to Jacobians computed via finite differences. In a next step, the simulation code for two-phase flow was coupled to a module calculating streaming potentials. Saturation dependence of the coupling coefficient and of electrical conductivity are taken into account when implementing this feature into SHEMAT-Suite. The objective of this new module was to examine the possibility of detection of CO$_2$ movement in the overburden of a storage reservoir due to CO$_2$ leakage through an abandoned well via low-cost self-potential (SP) measurements at the surface. In a typical leakage scenario including a leaky well with conductive metal casing, SP signals originating from injection can be identified at the surface. While SP monitoring can be applied for detecting abandoned wells, the smaller leakage signals, however, are masked by the injection signals. Three methods are proposed to detect the small leakage signals: i) simulate a non-leaky scenario and substract the simulated signals from the measured ones, ii) exploit the symmetry of the injection signal by analyzing potential differences of dipoles with a dipole center at the injection well, iii) measure SP signals when injection is interrupted. A prerequisite for method i) is a well-known reservoir structure, as, e.g., in a depleted oil or gas field. Best, use a combination of method i) and ii). In order to analyze the quality of these methods, three-dimensional numerical (SP) modeling of two-phase flow and the electrokinetic coupling between flow and streaming potential was performed.The SP signal is influenced by CO$_2$ injection, coupling coefficient, rock conductivity and the electrical conductivity of the well casing and is approximately ten times higher for conductive metal casings as compared to the presence of no metal casing.A major advantage of SP monitoring is the fact that SP signals at the leakage location arise shortly after injection started, even if CO$_2$ is still far away from the leaky well. The brine being displaced through the leaky well causes an SP signal. Thus, the measurement of SP signals is the only method which is able to detect a leakage before any fluid from injection actually leaks out of the reservoir. Hence, SP monitoring is recommended for subsurface CO$_2$ or fluid storage monitoring, but will benefit from a lower detection limit of the sensors. In a next step, the two-phase flow simulator has been enhanced, such that CO$_2$ may dissolve into brine and H$_2$O into gas making the (dis)appearance of phases possible. From the numerous different possibilities to deal with phase (dis)appearance, I chose the method of extended saturations. This method is relatively straight-forward to implement into an existing two-phase flow simulator and avoids primary variable switching (PVS) during phase changes. Switching of variables is associated with deteriorating convergence in Newton's method during phase change. The new method has been benchmarked against existing test cases. A primary feature of the two-phase, two-component module is the appearance of gravity fingering due to CO$_2$ dissolving into brine. To this end, I examined the three relevant trapping mechanisms, namely, stratigraphic, residual and solubility trapping for different synthetic reservoirs featuring homogeneous and heterogeneous permeability fields. In a last part, I implemented a pressure-enthalpy formulation for supercritical water/steam geothermal reservoirs. Relying on pressure and enthalpy no PVS is needed. During phase appearance a meaningful enthalpy can be calculated directly without the need to instantiate saturation with a small value (as needed for PVS methods). The implementation has been verified against existing test cases. In addition, the local geothermal reservoir around the well Venelle-2 has been simulated in a three-dimensional study. Further, Monte Carlo simulations of Venelle-2 try to explain the newly found elevated temperatures within the well.