A numerical method for the simulation of buoyancy-induced macrosegregation during solidification processes is presented. The physical model is based on volume averaged conservation equations for energy, mass, momentum, and solute. The resulting partial differential equations are solved by a finite element method, which considers two different length-scales: on the one hand, the scale of the overall process, on the other, a small critical zone near the solidification front where solutal inhomogeneities are initiated and the fluid velocity is non-zero. A domain decomposition method using two adapted grids has been developed. The overall computational domain is discretized using a 'coarse' finite element mesh adapted to the process scale. At each time step, the energy conservation equation is solved using this discretization and new values of temperature and solid fraction are calculated at each finite element node. Based on these values, the computational domain is subdivided into three subdomains: the so-called solid, mushy, and liquid regions. The mushy subdomain corresponds to the critical zone near the solidification front and is adaptively discretized with a finer finite element mesh, whereas the liquid domain uses the initial coarse grid and the solid is no longer considered. The fluid flow and solute transport equations are then solved on the different meshes using a Dirichlet-Neumann substructuring iterative method, together with the mortar technique to deal with the non-conforming discretizations at the subdomain interfaces. The method is applied to several test problems such as the Hebditch-Hunt macrosegregation experiment or the prediction of freckles, and the performances and limits of the approach are pointed out and discussed.