Numerical modeling of subsidence in saturated porous media: A mass conservative method

In this paper, a second order accurate cell-centered finite volume method (FVM) is coupled with a finite element method (FEM) to solve the deformation of a saturated porous layer based on Biot's consolidation model. The proposed numerical technique is applied to the fully unstructured triangular grids to simulate actual geological formations. To reconstruct the pressure gradient at control volume faces, the diamond scheme is implemented as a multipoint flux approximation method. Also the least square algorithm is used to interpolate pressure at the vertices from the cell-center values. The stability of this numerical model is studied in comparison to the different FEMs through various examples. It is shown that, although the Taylor-Hood FEM has been introduced as a remedy for violation of the inf-sup condition, it does not entirely remove the non-physical oscillations. Contrary to the linear and Taylor-Hood FEMs, the proposed discretization model provides monotonic solution without imposing any restriction on the mesh or time step size. Compared to the mixed FEM, the method achieves local mass balance with fewer degrees of freedom. To couple the flow and mechanical sub-problems, the fixed-stress operator split is implemented as an iterative sequential method, due to its unconditional stability, accuracy and high rate of convergence. The accuracy of the proposed model is verified via a range of examples including analytical and numerical solutions. The performance of this methodology is assessed through modeling of subsidence in an aquifer-interbed system. This problem illustrates the capability of the model in providing stable solution in heterogeneous domains with complicated shapes.