Rigorous derivation of relativistic energy-momentum relation

Since $P = Fv$ we have $$\frac{dE}{dt} = \frac{dp}{dt} v$$ by Newton's second law. Integrating both sides with respect to $t$ gives $$\int \frac{dE}{dt} \, dt = \int v \frac{dp}{dt} \, dt = \int v \, dp$$ by the chain rule, aka ordinary $u$-substitution. We have $$p = \gamma m v = \frac{m v}{\sqrt{1-v^2}} \quad \Rightarrow \quad dp = \frac{m \, dv}{(1-v^2)^{3/2}}$$ where I set $c = 1$ for convenience and used the quotient rule. Integrating with initial and final velocities zero and $v_0$ gives $$E(v_0) - E(0) = \int_0^{v_0} \frac{mv}{(1-v^2)^{3/2}} \, dv = \frac{m}{\sqrt{1 - v_0^2}} - m.$$ At this point we cannot proceed further since we don't know the constant of integration. One can show by physical arguments that $E(0) = m$. Thus $$E(v) = \frac{m}{\sqrt{1-v^2}}$$ as desired. This isn't a hard derivation, but you're right: a lot of textbooks botch it.

For completeness, here's an arguably cleaner and simpler formulation of @knzhou 's answer:

We obtain

$$E = \int_{0}^{x_0} (\frac{d}{dt} p) \space dx = \int_{0}^{t_0} (\frac{d}{dt} p) \space v \space dt = \int_{0}^{p_0} v \space dp = \int_{0}^{v_0} v \space (\frac{d}{dv} p) \space dv$$

by applying a sequence of reparametrizations $dx = v \space dt$, $dp = (\frac{d}{dt} p) \space dt$ and $dp = (\frac{d}{dv} p) \space dv$ to the integral for $E$. Since $ \frac{d}{dv} p = m \space (1 - \frac{v^2}{c^2})^{-3/2}$, it follows that

$$ E = \int_{0}^{v} \dfrac{m v}{(1-\frac{v^2}{c^2})^{3/2}} dv = \frac{mc^2}{(1 - \frac{v^2}{c^2})^{1/2}} - mc^2.$$

Defining the total energy $\Sigma = E + mc^2$, since $\Sigma = \gamma m c^2$ and $p = \gamma m v$, it is easy to see by direct computation that $\Sigma^2 - c^2 p^2 = m^2 c^4$, hence

$$\Sigma^2 = m^2 c^4 + c^2 p^2 \space .$$