I am given this following PDE with the initial and boundary conditions with $0 < r < 1$, $t > 0$, and $v_0$ being a constant:
$u,_t,_t + 2bu,_t = u,_r,_r + \frac{1}{r} u,_r$
$u(t,r=0) = 0, \vert u(t,r=0)\vert < \infty$
$u(t=0,r)=0, u,_t(t=0,r) = v_0$
Given this information I am suppose to obtain the following solution::
$$u(t,r) = 2v_0e^{-bt}\sum_{i=0}^\infty \frac{J_0(z_{0n}r)}{z_{0n}J_1(z_{0n})} \frac{\sin(t\sqrt{z_{0n}^2 -b^2})}{z_{0n}^2 -b^2}$$
I am at a loss on how to start this problem. My best guess would be use to separation of variables in which case I would get the following:
$u(t,r)=h(t)*f(r)$
$\frac{h''}{h'} + \frac{2bh'}{h} = \frac{f''}{f} + \frac{f'}{rf} = -\lambda$
With the separated ODEs:
$h'' + 2bh' +h\lambda = 0$
$f'' + f' + fr\lambda$ = 0
Question: Am I going in the right direction to obtain the given solution and if not what am I missing?
You are on the right track, assume that $u(r, t) = R(r)T(t)$, so the equation
$$ \frac{\partial^2 u}{\partial t^2} + 2\frac{\partial u}{\partial t} = \frac{\partial^2 u}{\partial r^2} + \frac{1}{r}\frac{\partial u}{\partial r} $$
transforms to
$$ \frac{1}{T}\left[\frac{{\rm d}^2T}{{\rm d}t^2} + 2b\frac{{\rm d}T}{{\rm d}t} \right] = \frac{1}{R}\left[\frac{{\rm d}^2R}{{\rm d}r^2} + \frac{1}{r} \frac{{\rm d}R}{{\rm d}r}\right] $$
as usual, you have two equations
\begin{eqnarray} \frac{{\rm d}^2T}{{\rm d}t^2} + 2b\frac{{\rm d}T}{{\rm d}t} &=& \alpha T \\ \frac{{\rm d}^2R}{{\rm d}r^2} + \frac{1}{r} \frac{{\rm d}R}{{\rm d}r} &=& \alpha R \end{eqnarray}
for some constant $\alpha$
The solution for $T$ is the one of a damped harmonic oscillator
$$ T(t) = A e^{-\mu_+ t} + B e^{-\mu_+ t} $$
where
$$ \mu_\pm = b \pm \sqrt{b^2 + \alpha} $$
With the condition $u(r, t = 0) = 0 = R(r)T(0) \Rightarrow T(0) = 0$ we conclude that $A + B = 0$. Moreover, since $\alpha > 0$ leads to a purely exponential decay, we also take $\alpha < 0$ to ensure oscillations. $\color{blue}{\alpha = -\lambda^2}$ so the solution for the time-part is
$$ T(t) = A e^{-b t}[e^{+ i t\sqrt{\lambda^2 - b^2}} - e^{- i t\sqrt{\lambda^2 - b^2}}] = A'e^{-bt}\sin\left(t\sqrt{\lambda^2 - b^2}\right) \tag{1} $$
For the radial part you have
$$ R(r) = C J_0 (\lambda r) + D Y_0 (\lambda r) $$
The function $Y_0$ diverges at the origin so we require $D = 0$, the result is then
$$ R(r) = C J_0 (\lambda r) \tag{2} $$
With the boundary condition $u(r = 1, t) = 0 = T(t)R(1) \Rightarrow R(1) = 0$, we get $J_0(\lambda) = 0$, the solutions to this equation are the zeros of the Bessel function
$$ \lambda = z_{0n} ~~~ n = 0, 1, \cdots $$
That means that our solution so far is
$$ u(r, t) = e^{-bt}\sum_{n = 0}^{+\infty} A_n J_0(z_{0n} r) \sin\left(t\sqrt{z_{0n}^2 - b^2}\right) $$
Now we can impose the final condition $\partial u(0, r)/\partial t = v_0$, which leads to
$$ \frac{\partial u(r, t = 0)}{\partial t} = v_0 = \sum_{n = 0}^{+\infty} A_n\sqrt{z_{0n}^2 - b^2}J_0(z_{0n}r) \tag{3} $$
I'm not going to show it here, but it is possible to show that if $f(x) = \sum_n a_n J_\nu(z_{\nu n}x)$ then
$$ a_n = \frac{2}{J_{\nu + 1}^2(z_{\nu n})} \int_0^1 x f(x)J_{\nu}(z_{\nu n}x){\rm d}x $$
if you apply this to Eq. $(3)$ with $\nu = 0$ you will find then
$$ \sqrt{z_{0n}^2 - b^2} A_n = \frac{2v_0 }{z_{0n}J_{1}(z_{0n})} $$
Putting everything together you then arrive to the final solution
$$ \bbox[5px, border:2px solid blue] { u(r, t) = 2v_0e^{-bt}\sum_{n = 0}^{+\infty} \frac{J_0(z_{0n} r)}{z_{0n}J_1(z_{0n})} \frac{\sin\left(t\sqrt{z_{0n}^2 - b^2}\right)}{\sqrt{z_{0n}^2 - b^2}} } $$