We present a theoretical framework to study open-charm production in relativistic heavy-ion collisions. Charm quarks are regarded as an effective probes of the deconfined medium formed in these collisions, the quark-gluon plasma (QGP). The initial conditions of such collisions are simulated with a transverse profile described via a Glauber-based model, and a longitudinal behavior modeled by a data-inspired parametrization. The space-time evolution of the temperature and the flow velocity field of the medium are quantified by means of (3+1)-dimensional relativistic viscous hydrodynamics. The Brownian motion of charm quarks propagating through the QGP is described by utilizing the Langevin transport equation. The subsequent hadronization is implemented via a dual model, including fragmentation and heavy-light coalescence mechanisms. In particular, in the coalescence also the contribution from higher hadronic states components is considered. The parameters of the model are tuned based on comparison to data. The coupling strength between the charm quarks and the QGP constituents, quantified by the spatial diffusion coefficient 2πT D s , is obtained by performing a phenomenological fit analysis to the lattice QCD calculations, resulting in 2πT D s = const. (model A) and 2πT D s = 1.3 + (T /T c ) 2 (model B). We find that the relative azimuthal distribution of the initially back-to-back generated cc pairs presents a broadening behavior, which is more pronounced for cc pairs with small initial p T and when the model B approach is adopted. The competition between the initial drag and the subsequent collective effects tends to restrict the time dependence of charm quark R AA . Concerning the theoretical uncertainty on final D-meson nuclear modification, the nuclear shadowing and pp baseline components are dominant at high and low p T (p T 3 GeV/c), respectively. The measured D-meson R AA (p T ) favors model A's assumption for the diffusion coefficient both at the Relativistic Heavy Ion Collider and the Large Hadron Collider, while v 2 (p T ) prefers model B at moderate p T . These results confirm the necessity to consider the temperature and/or momentum dependence of 2πT D s to simultaneously describe the D meson R AA and v 2 .