Abstract A hydrochemical model comprising a set of coupled partial differential equations has been proposed for the bioremediation of contaminated soil with a heterogeneous formation. The soil bed is a layer of water-saturated porous medium with dominant bypassing flow. Inside the layer, the groundwater flow is relatively slow such that the magnitudes of convection and dispersion are of the same order. This model describes, both spatially and temporally, the contaminant fate, microorganism growth and oxygen consumption. The results of simulation indicate that the velocity of the groundwater flow affects significantly the rate of contaminant depletion; the higher the velocity, the greater the availability of oxygen and the faster the substrate consumption. When the velocity is low, an anaerobic zone is generated in the clay layer, thereby rendering the biodegradation rate oxygen limited and, in turn, prolonging the remediation time.