  This work proposes an original preconditioner that couples the Constrained
Pressure Residual (CPR) method with block preconditioning for the efficient
solution of the linearized systems of equations arising from fully implicit
multiphase flow models. This preconditioner, denoted as Block CPR (BCPR), is
specifically designed for Lagrange multipliers-based flow models, such as those
generated by Mixed Hybrid Finite Element (MHFE) approximations. An original
MHFE-based formulation of the two-phase flow model is taken as a reference for
the development of the BCPR preconditioner, in which the set of system unknowns
comprises both element and face pressures, in addition to the cell saturations,
resulting in a $3 \times 3$ block-structured Jacobian matrix with a $2 \times
2$ inner pressure problem. The CPR method is one of the most established
techniques for reservoir simulations, but most research focused on solutions
for Two-Point Flux Approximation (TPFA)-based discretizations that do not
readily extend to our problem formulation. Therefore, we designed a dedicated
two-stage strategy, inspired by the CPR algorithm, where a block preconditioner
is used for the pressure part with the aim at exploiting the inner $2 \times 2$
structure. The proposed preconditioning framework is tested by an extensive
experimentation, comprising both synthetic and realistic applications in
Cartesian and non-Cartesian domains.
