The accurate prediction of shrinkage defects remains a major challenge in the numerical modeling of casting processes. Many existing FEM-based casting models rely on node-based linear expansion coefficients or density approximations to simulate volumetric changes during solidification. In contrast, widely used commercial casting packages often employ FDM/FVM- or VOF-based strategies. The present paper therefore focuses on an element-based FEM formulation that calculates volume changes directly from temperature-dependent density variations and extends the analysis to a multi-domain casting–mold system. This paper presents an element-based numerical methodology that calculates volume changes directly from temperature-dependent density variations. The proposed FEM framework adopts a multi-domain approach, encompassing both the casting and the molds. The solidification model employs a multi-step iterative correction algorithm: an initial heat conduction analysis is followed by an enthalpy-based latent heat correction, with element statuses identified through a tagging system. By updating density values at the element level, the method provides a physically consistent representation of the mass-volume balance throughout the cooling process. Heat transfer between the casting and the mold domains is modeled using the fourth-kind boundary condition. The validity of this approach is demonstrated through a comparative study of a lead (Pb) ingot, where the predicted shrinkage cavity geometry is benchmarked against experimental results, supported by indirect temperature measurements. The FEM approach described herein has been fully implemented in the author’s proprietary software.