This is a similar question to [this](https://mathoverflow.net/questions/470070/estimating-a-solution-to-an-euler-type-ode/470111?noredirect=1#comment1220446_470111) but with a different ODE. 

Let $f$ be a continuous function in $L^2([1,\infty)$ satisfying $\sup_{r\geq 1} r|f(r)| <\infty$. Let $\ell$ be a positive integer, $R>1$, and $a$ be a real number. 

Let $u(r)$ be a function on $[R,\infty)$ solving the following IVP: 
$$(r^2-1)u''(r)+2ru'(r) - \ell(\ell+1) u(r) = f(r), \quad r\in [R,\infty)$$
$$u(R)= a, \quad \lim_{r\to \infty} u(r) = 0$$

We should be able to estimate the solution $u$ by $f$ and $a$. I want to understand precisely how $\ell$ will appear in this estimate.

I am speculating the following: there exists a constant $C$ independent of $\ell$, $a$, $u$ and $f$ satisfying 
$$\ell \sup_{r\geq 1}r|u(r)| + \sup_{r\geq 1}r^2|u'(r)| \leq C(\ell^{-1}\sup_{r\geq 1} r|f(r)| + \ell^{-1}\lVert f\rVert_{L^2} + \ell|a|)$$

Can anyone verify if this is true? How would I prove this estimate (or the correct one in case the above is incorrect)? 

The solutions to the homogeneous problem are the Legendre polynomials of the first and second kind $P_{\ell}$, $Q_{\ell}$. Supposing $a=0$ for simplicity, we can use the method variation of parameters to get

$$u(r) = -P_{\ell}(r) \int_{r}^{\infty} Q_{\ell}(t) f(t)W(t)^{-1}dt + Q_{\ell}(r) \int_{R}^r P_{\ell}(t) f(t) W(t)^{-1} dt$$

where $W = P_{\ell}Q_{\ell}'-P_{\ell}'Q_{\ell}$. I have tried to use this formula to get the estimate but couldn't. 

Any help is appreciated.