We introduce a block Lanczos (BL) recursive technique to construct quasi-one-dimensional models, suitable for density-matrix renormalization group (DMRG) calculations, from single- as well as multiple-impurity Anderson models in any spatial dimensions. This new scheme, named BL-DMRG method, allows us to calculate not only local but also spatially dependent static and dynamical quantities of the ground state for general Anderson impurity models without losing elaborate geometrical information of the lattice. We show that the BL-DMRG method can be easily extended to treat a multi-orbital Anderson impurity model. We also show that the symmetry adapted BL bases can be utilized, when it is appropriate, to reduce the computational cost. As a demonstration, we apply the BL-DMRG method to three different models for graphene: (i) a single adatom on the honeycomb lattice, (ii) a substitutional impurity in the honeycomb lattice, and (iii) an effective model for a single carbon vacancy in graphene. Our analysis reveals that, for the particle-hole symmetric case at half filling of electron density, the ground state of model (i) behaves as an isolated magnetic impurity with no Kondo screening while the ground state of the other two models forms a spin singlet state. We also calculate the real-space dependence of the spin-spin correlation functions between the impurity site and the conduction sites for these three models. Our results clearly show that, reflecting the presence of absence of unscreened magnetic moment at the impurity site, the spin-spin correlation functions decay as $r^{-3}$, differently from the non-interacting limit ($r^{-2}$), for model (i) and as $ r^{-4}$, exactly the same as the non-interacting limit, for models (ii) and (iii) in the asymptotic $r$, where $r$ is the distance between the impurity site and the conduction site.