A microscopic understanding of the complex solute-defect interaction is pivotal for optimizing the alloy’s macroscopic mechanical properties. Simulating solute segregation in a plastically deformed crystalline system at atomic resolution remains challenging. The objective is to efficiently model and predict a physically informed segregated solute distribution rather than simulating a series of diffusion kinetics. To address this objective, we coupled molecular dynamics (MD) and Monte Carlo (MC) methods using a novel method based on virtual atoms technique. We applied our MD-MC coupling approach to model off-lattice carbon (C) solute segregation in nanoindented Fe-C samples containing complex dislocation networks. Our coupling framework yielded the final configuration through efficient parallelization and localized energy computations, showing C Cottrell atmospheres near dislocations. Different initial C concentrations resulted in a consistent trend of C atoms migrating from less crystalline distortion to high crystalline distortion regions. Besides unraveling the strong spatial correlation between local C concentration and defect regions, our results revealed two crucial aspects of solute segregation preferences: 1. defect energetics hierarchy and 2. tensile strain fields near dislocations. The proposed approach is generic and can be applied to other material systems as well.