We use direct numerical simulations to calculate the joint probability density function of the relative distance $R$ and relative radial velocity component $V_R$ for a pair of heavy inertial particles suspended in homogeneous and isotropic turbulent flows. At small scales the distribution is scale invariant, with a scaling exponent that is related to the particle-particle correlation dimension in phase space, $D_2$. It was argued [1, 2] that the scale invariant part of the distribution has two asymptotic regimes: (1) $|V_R| ll R$ where the distribution depends solely on $R$; and (2) $|V_R| gg R$ where the distribution is a function of $|V_R|$ alone. The probability distributions in these two regimes are matched along a straight line $|V_R| = z^ast R$. Our simulations confirm that this is indeed correct. We further obtain $D_2$ and $z^ast$ as a function of the Stokes number, ${rm St}$. The former depends non-monotonically on ${rm St}$ with a minimum at about ${rm St} approx 0.7$ and the latter has only a weak dependence on ${rm St}$.