The hadroproduction of heavy-flavoured mesons has recently attracted a growing interest e.g. within the people involved in global analysis of proton and nuclear parton distribution functions, saturation physics, and physics of cosmic rays. In particular, the D- and B-meson measurements of LHCb at forward direction are sensitive to gluon dynamics at small $x$ and are one of the few perturbative small-$x$ probes before the next generation deep-inelastic-scattering experiments. In this talk, we will concentrate on the collinear-factorization approach to inclusive D-meson production and describe a novel implementation --- SACOT-$m_{\rm T}$ --- of the general-mass variable flavour number scheme (GM-VFNS). In the GM-VFNS framework the cross sections retain the full heavy-quark mass dependence at $p_{\rm T}=0$, but gradually reduce to the ordinary zero-mass results towards asymptotically high $p_{\rm T}$. However, the region of small (but non-zero) $p_{\rm T}$ has been somewhat problematic in the previous implementations of GM-VFNS, leading to divergent cross sections towards $p_{\rm T} \rightarrow 0$, unless the QCD scales are set in a particular way. Here, we provide a solution to this problem. In essence, the idea is to consistently account for the underlying energy-momentum conservation in the presence of a final-state heavy quark-antiquark pair. This automatically leads to a well-behaved GM-VFNS description of the cross sections across all $p_{\rm T}$ without a need to fine tune the QCD scales. The results are compared with the LHCb data and a very good agreement is found. We also compare to fixed-order based calculations and explain why they lead to approximately a factor of two lower D-meson production cross sections than the GM-VFNS approach.