While $z=\frac{H}{c} d$ is more or less the relation that Hubble originally found, it doesn't hold out to arbitrarily large distances in FLRW cosmology. $v=HD$ does hold out to arbitrarily large distances, provided $v$ and $D$ are interpreted correctly as FLRW recessional velocity and FLRW spatial distance.
The real problem here is that the $v$ in $v=HD$ and the $v$ in $1+z=\sqrt{\frac{1+v/c}{1-v/c}}$ are not the same physical quantity. They shouldn't have used the same letter for them in a single paper, and they certainly shouldn't have punned on it by substituting one quantity into a formula that expects the other.
The $v$ in the special-relativistic formula is defined in terms of global inertial/Minkowski coordinates. In most FLRW cosmologies, you can't define global Minkowski coordinates because spacetime isn't flat, so you can't apply the special-relativistic formula because the $v$ in it is just meaningless.
However, in the zero-density ($\Omega=0$) limit of FLRW cosmology, spacetime is flat, you can define Minkowski coordinates on it, and the SR formula does work, out to arbitrary distances.
There are two zero-density FLRW cosmologies. One is boring: the scale factor $a$ is constant and $v_{\small\text{FLRW}} = v_{\small\text{SR}} = 0$ and both formulas give $z=0$. The other is much more interesting; it's called the Milne model and it describes a linearly expanding universe. In this case it turns out that $$v_{\small\text{SR}} / c = \tanh (v_{\small\text{FLRW}} / c)$$ which means that $v_{\small\text{FLRW}}$ is in SR terms the rapidity. If you plug this into the SR redshift formula you get, after a bit of manipulation, $$1+z = \exp(v_{\small\text{FLRW}} / c).$$
Meanwhile, in FLRW coordinates we have $a(t) = \dot at$ (for some constant $\dot a$) and, for objects moving with the Hubble flow, $$v_{\small\text{FLRW}} = a'(t)x = \dot ax$$ $$1+z = a(t_r)/a(t_e) = t_r/t_e$$ $$x = \int_{t_e}^{t_r} \frac{c\,\mathrm dt}{a(t)} = \frac{c}{\dot a} \int_{t_e}^{t_r} \frac{\mathrm dt}{t}$$ (using e and r subscripts for emission and reception respectively), and the integral of $1/t$ is $\ln t$, so $1+z = \exp(v_{\small\text{FLRW}} / c)$ as before.
The results agree because there is only one kind of redshift in general relativity, and the special-relativistic and cosmological formulas are special cases of it. Since they're different special cases, usually at most one of them is applicable to any given problem. But in the overlap of their zones of applicability, they're different coordinate descriptions of the same phenomenon, so they must agree.
I like Davis and Lineweaver's paper a lot, but they didn't catch every misconception about cosmology, not even all of their own misconceptions, and when they talk about intrinsic expansion of space they're just spreading another misconception. In reality there's no difference in GR between the relative motion of galactic superclusters and any other relative motion.