3D hydrodynamics models of deep stellar convection exhibit turbulent entrainment at the convective-radiative boundary which follows the entrainment law, varying with boundary penetrability. We implement the entrainment law in the 1D Geneva stellar evolution code. We then calculate models between 1.5 and 60 M$_{odot}$ at solar metallicity ($Z=0.014$) and compare them to previous generations of models and observations on the main sequence. The boundary penetrability, quantified by the bulk Richardson number, $Ri_{mathrm{B}}$, varies with mass and to a smaller extent with time. The variation of $Ri_{mathrm{B}}$ with mass is due to the mass dependence of typical convective velocities in the core and hence the luminosity of the star. The chemical gradient above the convective core dominates the variation of $Ri_{mathrm{B}}$ with time. An entrainment law method can therefore explain the apparent mass dependence of convective boundary mixing through $Ri_{mathrm{B}}$. New models including entrainment can better reproduce the mass dependence of the main sequence width using entrainment law parameters $A sim 2 times 10^{-4}$ and $n=1$. We compare these empirically constrained values to the results of 3D hydrodynamics simulations and discuss implications.