Sequential implicit nonlinear solver for geothermal simulation

Abstract Modeling of thermal multiphase fluid flow in subsurface porous media poses significant difficulties for nonlinear solvers. These difficulties are a result of the strong coupling between the mass and energy conservation equations through the thermodynamic relationships. In particular, for geothermal problems, this tight coupling results in an apparent “negative compressibility” for cells (control volumes) that have two fluid phases (liquid and vapor) during the time interval of interest. This behavior is associated with the flow of cold water into a cell that is at saturated conditions. As steam condensation occurs, the reduction in fluid volume due to the phase change is larger than the expansion of the liquid and vapor due to injection. As a result, the cell pressure decreases as condensation proceeds. This pressure decrease leads to an increase in the inflow of cold water until the condensation is complete. During the condensation period, the overall compressibility of the system is negative, and this poses serious convergence problems even when a fully implicit (i.e., backward Euler) time discretization is used. In order to deal with this problem, nonlinear solvers must be modified. In this work, we present a nonlinear solution strategy that overcomes the negative compressibility associated with the condensation front. We first consider a single-cell problem that demonstrates the negative compressibility phenomenon. Then, we describe a solution strategy that overcomes the nonlinear convergence difficulties. We then derive a criterion that can be used as a nonlinear preconditioning strategy to solve multicell problems. The effectiveness of the proposed nonlinear solution scheme is then demonstrated for more complex models. It was found that for one-dimensional problems, this nonlinear solution strategy converges for all timesteps sizes, while the fully coupled strategy only converges for a limited sized timestep. Furthermore, for a two-dimensional heterogeneous problem, the largest Courant–Friedrichs–Lewy (CFL) number for this method is at least an order of magnitude larger than the CFL number that can be used with the fully coupled approach.