DARSim VI: Geochemistry

“New framework for modeling of complex phase behavior involving aqueous electrolytes with the potential for precipitation and/or dissolution of multiple solid phases with application to CO2 sequestration and other subsurface flow and transport processes.”

Main ingredients of simulation framework

1. GHC equation of state

The GHC EOS (Lucia, 2010) is a recent modification of the SRK EOS that constrains the energy parameter, a, to satisfy the Gibbs-Helmholtz equation and uses Monte Carlo simulations to incorporate molecular length scale information. The novel features of the GHC equation include the use of molecular information in the energy parameter expression and the estimation of co-volume parameter b from pure component density data. The GHC EOS only uses parameters based on pure component properties (e.g., mixture critical properties from Kay's rules) and pure component Monte Carlo molecular simulation (mixture internal energies of departure are computed from a linear mixing rule), even in systems containing electrolytes, and is truly predictive. The GHC EOS implemented as a part of GFLASH library (Lucia, et al., 2015) which contains a robust stability/flash algorithm and the capability to handle simultaneous phase and homogeneous/heterogeneous chemical reaction equilibrium with mineral deposition/dissolution.

2 Chemical reaction equilibrium model

A full detailed description of the numerical methods used for solving the equilibrium reactions for molecular salt formation included in this work can be found in the paper by Lucia et al. (2014). For the examples below, the formation/dissolution of solid salts was limited to NaCl, Na2CO3, CaCl2, and CaCO3, described by the following reactions:

In addition, the following reaction of dissolved carbon dioxide with water to generate carbonate ion was included:

The equilibrium constants for the reactions were calculated from tabulated standard Gibbs free energies of formation data and corrected for temperature using tabulated standard enthalpies of formation and the van’t Hoff equation.    

3 Reservoir-scale simulator

For the flow and transport simulation we used Automatic Differentiation General Purpose Research Simulator (ADGPRS) (Voskov, 2011, Voskov and Tchelepi, 2012, Zhou 2014) which have a wide range of capabilities including:

1. Flexible treatment of all nonlinear physics.

2. A fully thermal-compositional formulation for any number of phases.

3. Multi-phase CSAT for efficient and robust computation of phase behavior.

4. A variety of spatial and temporal discretization schemes.

5. Thermal geo-mechanical modeling including the effects of fractures.

6. A fully coupled, thermal, multi-segmented well model with drift-flux.

7. An adjoint-based optimization module.

Architecture of the framework

Figure 1 illustrates the architecture of the ADGPRS/GFLASH modeling framework, with a focus on the flow of information in the GFLASH library.   Fluid densities and fugacities, as well as chemical reaction equilibrium constraints and their derivatives for the fluids in each grid block in the reservoir model, are calculated at given conditions (temperature, pressure and composition) and returned to the simulator. The molecular level information required to use the GHC equation of state is obtained a priori and stored in pure component look-up tables for the components of interest.

Figure 1: Schematic diagram illustrating information flow in ADGPRS/GFLASH framework

Small-scale CO2 sequestration model

It is believed that due to the time-scale of chemical reactions in brine-CO2 systems, mineral trapping does not affect the early stages of the CO2 sequestration process. However, small-scale precipitation may affect dissolution trapping due to the change in the dynamics of the instabilities. Precipitation and dissolution can change the porosity, which in turn can affect diffusion and finger formation. In the presence of a capillary transient zone, this process will be coupled to phase behavior and become quite challenging to predict without a reliable simulation tool. Here we present unique simulation results where all of these mechanistic details are taken into account.

For the simulation of CO2 fingering in the presence of chemical precipitation, the second type of small-scale model from Elenius et al. (2015) was used. In this model, the entire two-phase region was located in the high volume area, which maintains the original CO2 profile. This model represents the leading part of the CO2 plume. The higher concentration of CO2, dissolving in the area below the plum, triggers the carbonate reaction, which increases the concentration of CO2 and, in turn, initiates precipitation of CaCO3.

The influence of chemical precipitation on the generation of fingers is shown in Fig. 2. The top figures correspond to the concentration of CO2 in brine at different simulation times. This distribution is shown for the reference case where the precipitated mineral is changing (decreasing) the porosity based on a constant mineral volume and permeability of the reservoir using the Kozeny-Carman equation. The resulting porosity changes at different times are shown in the middle row of Fig. 2. The lower row of Fig. 2 shows the absolute difference in concentrations between the reference simulation shown in the upper row and the simulation without an update in porosity and permeability.

Figure 2: Small-Scale Simulation Results with CaCO3 Precipitation: (upper row) overall CO2 composition; (middle row) porosity; (lower row) difference in CO2 composition for cases with and without porosity and permeability update.

Note that the difference in CO2 distribution is insignificant, which clearly suggests that the influence of early-time mineralization can be ignored in the accurate estimation of the CO2 dissolution rate for the system studied. However, this does not guarantee that at different thermodynamic or chemical reaction  conditions, the same conclusion will hold.

In the next simulation we keep the same system, but amplify the molar volume of minerals to update the porosity 1000 fold. It can be seen in Fig.3 that the porosity variation in this case is more significant due to the larger pore volume occupied by minerals resulting from precipitation. In this case, the fingering process and the subsequent macroscopic dissolution rate have changed.

Figure 3: Small-Scale Simulation Results with CaCO3 Precipitation at Magnified Molar Volume:  (upper row) overall CO2 composition; (middle row) porosity; (lower row) difference in CO2 composition for the previous case and the case with magnified molar volume

Obviously, the last case is hypothetical since the molar volume of CaCO3 for the porosity update is unrealistically large.

In the next simulation, we assume a presence of NaCl with an initial concentration CNaCl = 1 kmol/m3 and the same under-saturated brine. The dynamics of the instability completely changes in this case, and no fingers are formed in the reservoir formation. See Fig. 10. Here, the increase in porosity due to the dissolution of NaCl in the diffusion zone and the reduction of solubility due to the presence of Na+ and Cl- ions stabilizes the CO2-brine interface and prevents fingers from growing.

Figure 4: Small-Scale Simulation Results with CaCO3 Precipitation and NaCl Dissolution: (upper row) overall CO2 composition; (lower row) porosity changes.


Elenius, M. T., Voskov, D. V., & Tchelepi, H. A. (2015). Interactions between gravity currents and convective dissolution. Advances in Water Resources, 83, 77–88. doi:10.1016/j.advwatres.2015.05.006

Lucia, A. (2010). A multiscale Gibbs-Helmholtz Constrained cubic equation of state. Journal of Thermodynamics, 2010, 1–10. doi:10.1155/2010/238365

Lucia, A., Henley, H., & Thomas, E. (2014). Multiphase equilibrium flash with salt precipitation in systems with multiple salts. Chemical Engineering Research and Design, 93(May), 662–674. doi:10.1016/j.cherd.2014.04.034

Voskov, D.V. (2011). An extended natural variable formulation for compositional simulation based on tie-line parameterization, Transport in Porous Media, 92(3), 541-557. doi:10.1007/s11242-011-9919-2

Voskov, D. V., & Tchelepi, H. A. (2012). Comparison of nonlinear formulations for two-phase multi-component EoS based simulation. Journal of Petroleum Science and Engineering, 82-83, 101–111. doi:10.1016/j.petrol.2011.10.012

Zhou, Y. [2012] Parallel General-Purpose Reservoir Simulation with Coupled Reservoir Models and Multi-Segment Wells. PhD Thesis, Stanford University.

Figure 1: Physics-based operators of mass accumulation for gas component (a), mass flux for gas component (b), energy accumulation (c) and energy flux (d) terms in a high enthalpy geothermal model with methane co-production. Click to enlarge.