In this paper a numerical analysis of in situ-biorestoration is presented. Improved stability and error estimates are derived for a P 1 finite-element method applied to coupled system of non linear partial differential equations modeling flow in porous media. These results improve upon previously derived stability and error estimates in two respects: first, a stability of the approximate solution is demonstrated in the case of a non linear diffusive flux λ ( u ) ∇ u with weaker norm assumptions than before and depends only on the initial conditions and time interval, and second, error estimates are optimal as in linear case. Extensions include the finite element approximation of flow field, described by Darcy’s law under the stream function–vorticity formulation because of its multiple advantages. Finally, numerical simulations for a correlation exponential plume are presented.