The computational burden associated with multidimensional, multicompoment, numerical solute transport models can be prohibitive. To solve these CPU intensive problems efficiently and accurately, we have implemented and extended a local adaptive grid refinement method of Berger and Oliger (1984) to track error-prone regions and supply high-resolution subgrids where they are locally needed while maintaining relatively few nodes elsewhere on a coarse, base grid. Novel features include a unique method for detecting a priori where the numerical error is unacceptable, variable time step control which allows smaller time steps on subgrids than on the base grid, and a modular framework which allows easy exchange of partial differential equation solvers to accommodate different problem formulations. Three examples show how the subgrids track the error-prone regions, how multiple subgrids are used for geometrically complex front shapes, and how increased resolution is provided for self-sharpening fronts when multiple reacting species are considered. For all three examples, solutions with accuracies comparable to those achieved with a uniform fine grid are obtained at between 20 and 30% of the computational cost. ¿ American Geophysical Union 1994