What concise method can solve the ODE, $y'' + \left(\frac{1}{x}+2x \right)y' + 4y =0$?

Question

While solving the vorticity transport equation in cylindrical coordinates, $\frac{\partial \omega}{\partial t}=\nu \left(\frac{\partial^2 \omega}{\partial r^2}+\frac{1}{r}\frac{\partial \omega}{\partial r}\right)$ and performing a transformation by a similarity solution, $\omega (r,t)=\frac{f(\eta)}{t}$ (where $\eta = \frac{r}{\sqrt{4\nu t}}$), I reduced the PDE to the ordinary Sturm-Liouville form, \begin{align*} f'' + \left(\frac{1}{\eta}+2\eta \right)f' +4f &= 0 \end{align*} with a solution of the form, $f(\eta)=c_1 f_1(\eta) +c_2 f_2(\eta)$.
As far as I know, Frobenius's method is the only approach that I've been able to obtain the general solution, \begin{align*} f(\eta)=c_1 e^{-\eta^2}+c_2 \operatorname{Ei}(\eta^2) e^{-\eta^2} \end{align*}

$\textbf{Attempt:}$ Assuming a solution (and its derivatives) to the DE, $y'' + \left(\frac{1}{x}+2x \right)y' +4y = 0$ is of the form, \begin{align*} y(x) &= \sum_{k=0}^{\infty}c_k x^{k+r}\\ y'(x) &= \sum_{k=0}^{\infty}c_k (k+r)x^{k+r-1}\\ y''(x) &= \sum_{k=0}^{\infty}c_k (k+r)(k+r-1)x^{k+r-2} \end{align*} By substitution, \begin{align*} 0&= \sum_{k=0}^{\infty}c_k (k+r)(k+r-1)x^{k+r-2} + \left( \frac{1}{x}+2x\right)\sum_{k=0}^{\infty}c_k (k+r)x^{k+r-1} + 4\sum_{k=0}^{\infty}c_k x^{k+r}\\ &= \sum_{k=0}^{\infty}c_k (k+r)(k+r-1)x^{k+r-2} + \left\{\frac{1}{x}\sum_{k=0}^{\infty}c_k (k+r)x^{k+r-1} + 2x\sum_{k=0}^{\infty}c_k (k+r)x^{k+r-1} \right\} + 4\sum_{k=0}^{\infty}c_k x^{k+r}\\ &= \sum_{k=0}^{\infty}c_k (k+r)(k+r-1)x^{k+r-2} + \sum_{k=0}^{\infty}c_k (k+r)x^{k+r-2} + 2\sum_{k=0}^{\infty}c_k (k+r)x^{k+r} + 4\sum_{k=0}^{\infty}c_k x^{k+r}\\ &= \sum_{k=0}^{\infty}c_k \left[(k+r)(k+r-1) + (k+r) \right]x^{k+r-2} + \sum_{k=0}^{\infty} 2c_k \left(k+r + 2\right)x^{k+r} \\ &= \sum_{k=0}^{\infty}c_k \left(k+r\right)^2x^{k+r-2} + \sum_{k=0}^{\infty} 2c_k \left(k+r + 2\right)x^{k+r} \\ &= \sum_{k=0}^{\infty}c_k \left(k+r\right)^2x^{k+r-2} + \sum_{k-2=0}^{\infty} 2c_{k-2} \left(k-2+r + 2\right)x^{k-2+r} \\ &= \sum_{k=0}^{\infty}c_k \left(k+r\right)^2x^{k+r-2} + \sum_{k=2}^{\infty} 2c_{k-2} \left(k+r\right)x^{k+r-2} \\ &= c_0 r^2 x^{r-2}+c_1 \left(1+r\right)^2 x^{r-1}+\sum_{k=2}^{\infty} \left[ c_k \left(k+r\right)^2 + 2c_{k-2} \left(k+r\right) \right]x^{k+r-2} \end{align*} With $c_1=0$, we have a double root in which $r^2 =0$ and $r_1=r_2=0$, thus the indicial equation becomes, \begin{align*} c_k k^2 + 2c_{k-2} k=0\\ \therefore c_k = -\frac{2c_{k-2}}{k} \:\:,\:\: k\geq 2 \end{align*} Testing the first few values of $c_k$, we observe the pattern,

$c_1=0$

$c_2 =-c_0$

$c_3=0$

$c_4=-\frac{2}{4}c_2 = \frac{1}{2}c_0$

$c_5=0$

$c_6=-\frac{1}{3}c_4 =-\frac{1}{3}\left(-\frac{1}{2}c_2 \right)=-\frac{1}{3}\left(\frac{1}{2}c_0 \right)$

$c_7=0$

$c_8 = -\frac{2}{8}c_6 =-\frac{1}{4}\left( -\frac{1}{3}\left(\frac{1}{2}c_0 \right)\right) $

$c_9=0$

Therefore, $\left\{ -c_0 , \frac{1}{2!}c_0,-\frac{1}{3!}c_0, \frac{1}{4!}c_0, -\frac{1}{5!}c_0 ,\dots \right\}$ and the first solution, $y_1 (x)$, becomes, \begin{align*} y_1 (x) &= \left | x\right |^{r_{1}} \sum_{k=0}^{\infty}c_k x^{k} \:\:\:,\:\:\: c_{2m-1}=0\:\:,m\in \mathbb{N}\\ &= (1)\left( c_0 -\frac{c_0}{1!}x^2 + \frac{c_0}{2!}x^4 - \frac{c_0}{3!}x^6 +\dots \right)\\ &= \sum_{k=0}^{\infty}c_0\frac{(-x^2)^k}{k!}\\ &=e^{-x^2} \:\:,\:\: c_0=1 \end{align*} The second solution, $y_2(x)$, can be found by substituting $y_2 (x)=u(x)e^{-x^2}$ into the original DE: \begin{align*} 0 &= y'' + \left(\frac{1}{x}+2x \right)y' +4y \\ 0 &= \left\{ u''e^{-x^2} - 4xu'e^{-x^2} -2ue^{-x^2}+4x^2 ue^{-x^2} \right\} + \left(\frac{1}{x}+2x \right)\left\{ u'e^{-x^2} -2xue^{-x^2} \right\} +4ue^{-x^2} \\ &\dots \textit{lots of stuff}\dots\\ u''&= -\left(\frac{1}{x}-2x \right)u' \:\:\:,\:\:\: \textit{solve with,}\:\: w(x)=u'(x)\\ u(x)&= \frac{C_1}{2} \operatorname{Ei}(x^2) \:\:\:,\:\:\: C_2 =C_1 /2\\ \therefore y_2 (x)&= C_2 \operatorname{Ei}(x^2) e^{-x^2} \end{align*} The general solution is then, \begin{align*} y(x) =c_1 y_1 (x) + c_2 y_2 (x)\\ \boxed{y(x)= c_1 e^{-x^2} + c_2 \operatorname{Ei}(x^2) e^{-x^2} } \end{align*}

$\textbf{Question:}$ is there another (faster) approach?

Answer 1

Multiplying by $x$ gives \begin{align*} 2 y^{\prime} x^{2}+x y^{\prime \prime}+4 x y+y^{\prime}&=0 \end{align*} Integrating both sides of the ODE w.r.t $x$ gives \begin{align*} \int \left(2 y^{\prime} x^{2}+x y^{\prime \prime}+4 x y+y^{\prime}\right)d x &= 0 \\ x y^{\prime}+2 x^{2} y = c_1 \end{align*} Which is now solved for $y$. This is linear first order ode whose solution can easily be found to be (using an integrating factor). $$ y = {\mathrm e}^{-x^{2}} \left(-\frac{c_1 \,\operatorname{Ei}_{1}\left(-x^{2}\right)}{2}+c_2 \right) $$

Another method is to use Kovacic algorithm. But the above is shorter.

Answer 2

For notation's sake, I work with the ODE in the title, $$ y'' + \left( \frac 1 x + 2x \right) y' + 4y = 0 \tag{1} $$ I then try a substitution $y(x) := u(x) e^{-x^2}$, in the interest of killing of the $2x$ in the parenthetical. (Something about this form of an ODE is just begging it for me; I have no real intuition for why, I guess just old textbook exercises screaming to me?)

Then $$ \begin{align*} y'(x) &= e^{-x^2} \Big( - 2x u(x) + u'(x) \Big) \\ y''(x) &= e^{-x^2} \Big( u''(x) - 4 x u'(x) + 4 x^2 u(x) - 2 u(x) \Big) \end{align*}$$ This turns $(1)$ into, with some minor algebra, $$ u''(x)e^{-x^2} - 2x u'(x)e^{-x^2} + \frac 1 x u'(x) e^{-x^2} =0 $$ Divide through by $e^{-x^2}$, and factor the $u'$ terms to get $$ u'' + \left( \frac 1 x - 2x \right) u' = 0 $$ and perform a reduction of order substitution, $v'=u$, to get $$ v' + \left( \frac 1 x - 2x \right) v = 0 \tag{2} $$ This ODE is separable, with solution $$ v(x) = C \frac{e^{x^2}}{x} \tag{3} $$ where $C$ is a constant dependent on initial data.

Then one only need integrate $(3)$ to get $u$ (as $u = \int u' = \int v$): $$ u(x) = C \int \frac{e^{x^2}}{x} \, \mathrm{d}x $$ In the integral, perform the substitution $w=x^2$ to get $$ \int \frac{e^{x^2}}{x} \, \mathrm{d}x = \frac 1 2 \int \frac{e^{w}}{w} \, \mathrm{d}w = \frac 1 2 \operatorname{Ei}(w) = \frac 1 2 \operatorname{Ei}(x^2) $$ so (abusing notation with the constants of integration $C$ and now $K$) $$ u(x) = C \operatorname{Ei}(x^2) + D \tag{5} $$ Since we defined $y(x) = u(x) e^{-x^2}$, we now have $$ y(x) = C e^{-x^2} \operatorname{Ei}(x^2) + D e^{-x^2} $$ as desired.

Answer 3

In accordance to Nasser's answer, I shall elaborate on the nature of this ODE. It is a special Sturm-Liouville DE in which, \begin{align}\tag{1} p(x)y'' + \frac{d}{dx}\left[q(x)y \right]&= 0 \end{align} for some differentiable functions, $p(x)$, $q(x)$, and $y(x)$; $y: \mathbb{R}\to \mathbb{R}$, $x\in \mathbb{R}$. Using integrations by parts, we have \begin{align}\tag{2} \int \left[ p(x)y'' + q(x)y' + q'(x)y \right] dx &= C \\ \left\{p(x)y' -\int p'(x)y' dx\right\} + \left\{q(x)y -\int q'(x)y dx\right\} + \int q'(x)y dx &= C \\ p(x)y' -\int p'(x)y' dx + q(x)y &= C \end{align} An easily obtainable solution to this DE is if $p'(x)$ is a constant. Let $q(x)=w(x)+p'(x)$. Then, \begin{align}\tag{3} p(x)y' -p'(x)\int y' dx + \left[w(x)+p'(x) \right]y &= C \\ p(x)y' -p'(x)y + w(x)y +p'(x)y &= C \\ p(x)y' + w(x)y &= C\\ y' + \frac{q(x)-p'(x)}{p(x)}y &= \frac{C}{p(x)} \end{align} If $p'(x)=c_1$ for some $c_1 \in \mathbb{R}$, then $p(x)=c_1 x$. Equation (1) becomes, \begin{align}\tag{4} c_1xy'' + \left[w(x)+c_1 \right]y' + w'(x)y &=0 \end{align} and Equation (3) likewise becomes, \begin{align*} y'+\frac{w(x)}{p(x)}y &=\frac{C}{p(x)} \end{align*} Equations (1) and (4) have the solution, \begin{align*} & y'e^{\int \frac{w(x)}{p(x)}dx}+\frac{w(x)}{p(x)}e^{\int \frac{w(x)}{p(x)}dx}y =\frac{C}{p(x)} e^{\int \frac{w(x)}{p(x)}dx}\\ &\int \frac{d}{dx}\left[ ye^{\int \frac{w(x)}{p(x)}dx} \right]dx =C \int \frac{e^{\int \frac{w(x)}{p(x)}dx}}{p(x)}dx\\ &y(x) =C e^{-\int \frac{w(x)}{p(x)}} \int \frac{e^{\int \frac{w(x)}{p(x)}dx}}{p(x)}dx \\ &y(x) = C e^{-\int \frac{q(x)-p'(x)}{p(x)}} \int \frac{e^{\int \frac{q(x)-p'(x)}{p(x)}dx}}{p(x)}dx \end{align*} More generally, suppose that in, \begin{align}\tag{5} y'' +f(x)y' +g(x)y =0 \end{align} there exists some $p(x)$ such that, \begin{align*} p(x)y'' + p(x)f(x)y' +p(x)g(x)y &=0 \\ \equiv \:\: p(x)y'' + \frac{d}{dx}\left[ q(x)y \right] &=0 \end{align*} for which, \begin{align}\label{A5} \left\{\begin{matrix} q(x)=p(x)f(x) \\ q'(x)=p(x)g(x) \end{matrix}\right. \end{align} Then solving for $p(x)$, \begin{align}\tag{6} \frac{d}{dx}\left[ p(x)f(x) \right] &= p(x)g(x)\\ p'(x)f(x) &=p(x) \left[ g(x)-f'(x)\right] \\ \frac{dp}{p} &= \int \frac{g(x)-f'(x)}{f(x)}dx \\ & \boxed{p(x) = e^{\int \frac{g(x)-f'(x)}{f(x)}dx}} \end{align} Using the integration by parts from Equation (2), \begin{align*} p(x)y' -\int p'(x)y' dx + q(x)y &= C\\ p(x)y' -\int p'(x)y' dx + p(x)f(x)y &= C \end{align*} and necessitating $p'(x)=\textit{constant}$, we obtain \begin{align*} p(x)y' - p'(x)y + p(x)f(x)y &= C\\ p(x)y' + \left[ p(x)f(x)-p'(x) \right]y &= C\\ y' + \left[ f(x)- \frac{p'(x)}{p(x)} \right]y &= \frac{C}{p(x)}\\ y' + \left[ f(x)- \frac{\frac{g(x)-f'(x)}{f(x)}e^{\int \frac{g(x)-f'(x)}{f(x)}dx}}{e^{\int \frac{g(x)-f'(x)}{f(x)}dx}} \right]y &= C e^{-\int \frac{g(x)-f'(x)}{f(x)}dx}\\ y' + \left[ f(x) - \frac{g(x)-f'(x)}{f(x)} \right]y &= C e^{-\int \frac{g(x)-f'(x)}{f(x)}dx} \end{align*} Let $A(x)=f(x) - \frac{g(x)-f'(x)}{f(x)}$ and $B(x)=C e^{-\int \frac{g(x)-f'(x)}{f(x)}dx}$ so that, \begin{align*} y' + A(x)y &= B(x) \:\:\:,\:\:\: \text{IF}=e^{\int A(x)dx}\\ \frac{d}{dx}\left[ ye^{\int A(x)dx} \right]&= B(x)e^{\int A(x)dx}\\ y &=e^{-\int A(x)dx}\int B(x)e^{\int A(x)dx}dx\\ y &=e^{-\int \left[ f(x) - \frac{g(x)-f'(x)}{f(x)} \right] dx}\int C e^{-\int \frac{g(x)-f'(x)}{f(x)}dx} e^{\int \left[ f(x) - \frac{g(x)-f'(x)}{f(x)} \right] dx}dx \end{align*}

\begin{align} \boxed{y =Ce^{-\int \left[ f(x) - \frac{g(x)-f'(x)}{f(x)} \right] dx}\int e^{\int \left[ f(x) - 2\frac{g(x)-f'(x)}{f(x)} \right] dx}dx } \end{align} This is a solution to Equation (5) if and only if there exists a $p(x)$ such that $p'(x)=\textit{constant}$. Concisely, we can inspect (5) using Equation (6) and check whether such a $p(x)$ could be found. \begin{align}\label{A7} \begin{split} \frac{d}{dx}\left[e^{\int \frac{g(x)-f'(x)}{f(x)}dx} \right] &=c_1\\ e^{\int \frac{g(x)-f'(x)}{f(x)}dx} &=c_1x\\ \int \frac{g(x)-f'(x)}{f(x)}dx &=\ln (c_1x)\\ \frac{g(x)-f'(x)}{f(x)} &=\frac{1}{x}\\ f(x)+xf'(x)&= xg(x) \\ \frac{d}{dx}\left[ xf(x) \right]&= xg(x) \\ f(x)&=\frac{1}{x}\left[\int xg(x) dx +c_1 \right] \:\:\:,\:\:\: c_1\in \mathbb{R} \end{split} \end{align}

Thus, Equation (5) is equivalent to \begin{align*} y'' + f(x)y' + \frac{1}{x}\frac{d}{dx}\left[ xf(x)\right] y &=0\\ xy'' + \frac{d}{dx}\left[ \left( xf(x)\right)y \right] &=0 \:\:\:,\:\:\: q(x)=xf(x) \\ \end{align*}