We present numerical experiments aimed to study the correlation between the magnetic field strength, $B$, and the density, $n$, in the cold atomic interstellar medium (CNM). We analyze 24 magneto-hydrodynamic models with different initial magnetic field intensities ($B_0=$0.4, 2.1, 4.2, and 8.3 $μ$G) and/or mean densities (2, 3, and 4 cm$^{-3}$), in the presence of driven and decaying turbulence, with and without self-gravity, in a cubic computational domain with 100 pc by side. Our main findings are: i) For forced simulations, which reproduce the main observed physical conditions of the CNM in the Solar neighborhood, a positive correlation between $B$ and $n$ develops for all the $B_0$ values. ii) The density at which this correlation becomes significant ($\lesssim 30$ cm$^{-3}$) depends on $B_0$ but is not sensitive to the presence of self-gravity. iii) The effect of self-gravity, when noticeable, consists of producing a shallower correlation at high densities, suggesting that, in the studied regime, self-gravity induces motions along the field lines. iv) Self-gravitating decaying models where the CNM is subsonic and sub-Alfvénic with $β\lesssim 1$ develop a high density positive correlation whose slopes are consistent with a constant $β(n)$. v) Decaying models where the low density CNM is subsonic and sub-Alfvénic with $β>1$ show a negative correlation at intermediate densities, followed by a high density positive correlation.

Accepted for publication in MNRAS, 10 pages, 11 figures