In this paper, we propose a mass conservative pseudostress-based finite element method for solving the Stokes problem with both Dirichlet and mixed boundary conditions. We decompose the velocity by means of a Helmholtz decomposition and derive a three-field mixed variational formulation, where the pseudostress, the velocity, both in \({\textbf {H}}(\text {div}\,)\) , and an additional unknown representing the null function, are the main unknowns of the system. By employing suitable finite element spaces, the velocity is approximated using \({\textbf {H}}(\text {div}\,)\) -conforming finite elements, ensuring the desired mass conservation property. The proposed method offers several advantages, including simplicity of implementation and compatibility with existing software packages for partial differential equation solvers. In addition, other variables of interest can be approximated through a simple postprocessing of the pseudostress without applying any numerical differentiation, thus, avoiding further sources of error. We extend the study to incorporate mixed boundary conditions for the Stokes problem and complement the analysis with the introduction of a reliable and efficient residual-based a posteriori error estimator. Numerical examples are provided to validate the theoretical results, demonstrating the effectiveness and accuracy of the proposed method.
