We propose and analyze a new mixed finite element method for the nonlinear problem given by the coupling of the steady Brinkman–Forchheimer and double-diffusion equations. Besides the velocity, temperature, and concentration, our approach introduces the velocity gradient, the pseudostress tensor, and a pair of vectors involving the temperature/concentration, its gradient and the velocity, as further unknowns. As a consequence, we obtain a fully mixed variational formulation presenting a Banach spaces framework in each set of equations. In this way, and differently from the techniques previously developed for this and related coupled problems, no augmentation procedure needs to be incorporated now into the formulation nor into the solvability analysis. The resulting non-augmented scheme is then written equivalently as a fixed-point equation, so that the well-known Banach theorem, combined with classical results on nonlinear monotone operators and Babuška–Brezzi’s theory in Banach spaces, are applied to prove the unique solvability of the continuous and discrete systems. Appropriate finite element subspaces satisfying the required discrete inf-sup conditions are specified, and optimal a priori error estimates are derived. Several numerical examples confirm the theoretical rates of convergence and illustrate the performance and flexibility of the method.