In this contribution, we make use of the QCD perturbative fragmentation function formalism to describe the one-particle inclusive fragmentation of a heavy quark produced in \(e^+e^-\) annihilation at \(\mathcal{O}(\alpha_S^2)\). We perform the computation analytically in Mellin-space. We resum soft-gluons effects in initial conditions and coefficient functions and perform evolution up to NNLL accuracy, obtaining the first NNLO + NNLL prediction for charm quark production. We study the impact of different matching schemes and Landau pole prescriptions in soft-gluon resummation, finding significant differences. We extract simple non-perturbative fragmentation functions for \(B\) and \(D^*\) mesons by comparing the perturbative prediction with the data from CLEO, BELLE and LEP experiments. We find that for charm mesons the experimental results from CLEO/BELLE and from LEP are not reconcilable with the standard DGLAP evolution.