This paper presents an adaptive multigrid approach, combining multigrid methods and adaptive local grid refinement, to solve 3D density-dependent flow and transport problems in the subsurface. This approach is incorporated with the Galerkin finite-element method to solve the modified Richards' equation in the flow module and the Lagrangian-Eulerian finite-element method to solve the advection-dispersion equation in the transport module, which are coupled through density effects. With multigrid methods, the linear/linearized matrix equations can be solved with O(n) manipulations to save computer time. With adaptive local grid refinement, computational accuracy is improved by refining only rough regions of the problem domain and the computation is efficiently achieved because computer efforts are focused on those regions. In this work, rough regions are determined by examining the mesh Peclet number during each nonlinear iteration of the flow module and by checking the smoothness of the Lagrangian concentra...