This paper focuses on multigrid methods for flow in heterogeneous porous media. We consider Darcy flow and the local permeability K(x) being a stationary random field of lognormal distribution. We apply the recently developed coarse graining method for the numerical upscaling of permeability, and develop a new multigrid method which applies this technique to obtain the coarse grid operators. The coarse grid operators are adjusted to the scale-dependent behaviour of the system as it incorporates only fluctuations of K on larger scales. This kind of action is essential for an efficient interplay with simple smoothers. We investigate important properties of the new multigrid method such as dependence on the boundary conditions and on grid refinement for the coarse graining and dependence on the mesh size. We compare the resulting method with the algebraic method of Ruge and Stüben, a Schurcomplement method, and matrix-dependent multigrid methods by solving the flow equation with K being random realizations as well as periodic media. The numerical convergence rates show that the new method is as efficient as the algebraic methods for variances σ2f ≤ 3 of K.