A positive semi-definite problem of three-dimensional incompressible contamination treatment from nuclear waste in porous media is discussed in this paper. The mathematical model is defined by a nonlinear initial-boundary system consisting of partial differential equations. Four important equations (an elliptic equation, two convection-diffusion equations and a heat conductor equation) determine the physical features. Considering the physical natures and computational efficiency, the authors introduce the conservative mixed finite volume element, upwind approximation and multistep difference to solve this system. The pressure and Darcy velocity are computed by a mixed finite volume element. The concentrations and temperature are solved by a combination of upwind approximation, multistep difference and mixed finite volume element. A multistep difference is used for approximating the partial derivative with respect to time. Mixed finite volume element and upwind differences are given for solving the convection-diffusions equations. Numerical dispersion and nonphysical oscillations could be eliminated, and the computational efficiency is improved by using a large time step. Furthermore, a conservative law is preserved and error estimates in $L^2$-norm is obtained. Finally, two numerical experiments are given to show the efficiency and possible applications.