We present a first calculation of the heavy flavor contribution to the longitudinally polarized deep-inelastic scattering structure function $g_1^{Q}$, differential in the transverse momentum or the rapidity of the observed heavy quark $Q$ or antiquark $overline Q$. All results are obtained at next-to-leading order accuracy in QCD within the framework of a newly developed parton-level Monte Carlo generator that also allows one to study observables associated with the produced heavy quark pair such as its invariant mass distribution or its correlation in azimuthal angle. First phenomenological studies are carried out for various heavy quark distributions in a kinematic regime relevant for a future Electron-Ion Collider with a particular emphasis on the expected size of the corresponding double-spin asymmetries and their sensitivity to the still poorly constrained helicity gluon distribution. Theoretical uncertainties associated with the choice of the factorization scale are discussed for selected observables.