  Strong coupling between geomechanical deformation and multiphase fluid flow
appears in a variety of geoscience applications. A common discretization
strategy for these problems is a continuous Galerkin finite element scheme for
the momentum balance equations and a finite volume scheme for the mass balance
equations. When applied within a fully-implicit solution strategy, however,
this discretization is not intrinsically stable. In the limit of small time
steps or low permeabilities, spurious oscillations in the pressure field, i.e.
checkerboarding, may be observed. Further, eigenvalues associated with the
spurious modes will control the conditioning of the matrices and can
dramatically degrade the convergence rate of iterative linear solvers. Here, we
propose a stabilization technique in which the balance of mass equations are
supplemented with stabilizing flux terms on a macroelement basis. The
additional stabilization terms are dependent on a stabilization parameter. We
identify an optimal value for this parameter using an analysis of the
eigenvalue distribution of the macroelement Schur complement matrix. The
resulting method is simple to implement and preserves the underlying sparsity
pattern of the original discretization. Another appealing feature of the method
is that mass is exactly conserved on macroelements, despite the addition of
artificial fluxes. The efficacy of the proposed technique is demonstrated with
several numerical examples.
