  Many applications involving porous media--notably reservoir engineering and
geologic applications--involve tight coupling between multiphase fluid flow,
transport, and poromechanical deformation. While numerical models for these
processes have become commonplace in research and industry, the poor
scalability of existing solution algorithms has limited the size and resolution
of models that may be practically solved. In this work, we propose a two-stage
Newton-Krylov solution algorithm to address this shortfall. The proposed solver
exhibits rapid convergence, good parallel scalability, and is robust in the
presence of highly heterogeneous material properties. The key to success of the
solver is a block-preconditioning strategy that breaks the fully-coupled system
of mass and momentum balance equations into simpler sub-problems that may be
readily addressed using targeted algebraic methods. Numerical results are
presented to illustrate the performance of the solver on challenging benchmark
problems.
