In [1]:
using ApproxFun, Plots, ComplexPhasePortrait, ApproxFun,  SingularIntegralEquations,
        SpecialFunctions, DualNumbers, SO

using SingularIntegralEquations.HypergeometricFunctions
gr();

M3M6: Methods of Mathematical Physics

$$ \def\dashint{{\int\!\!\!\!\!\!-\,}} \def\infdashint{\dashint_{\!\!\!-\infty}^{\,\infty}} \def\D{\,{\rm d}} \def\E{{\rm e}} \def\dx{\D x} \def\dt{\D t} \def\dz{\D z} \def\ds{\D s} \def\C{{\mathbb C}} \def\R{{\mathbb R}} \def\H{{\mathbb H}} \def\CC{{\cal C}} \def\HH{{\cal H}} \def\FF{{\cal F}} \def\I{{\rm i}} \def\Ei{{\rm Ei}\,} \def\qqqquad{\qquad\qquad} \def\qqand{\qquad\hbox{and}\qquad} \def\qqfor{\qquad\hbox{for}\qquad} \def\qqwhere{\qquad\hbox{where}\qquad} \def\Res_#1{\underset{#1}{\rm Res}}\, \def\sech{{\rm sech}\,} \def\acos{\,{\rm acos}\,} \def\erfc{\,{\rm erfc}\,} \def\vc#1{{\mathbf #1}} \def\ip<#1,#2>{\left\langle#1,#2\right\rangle} \def\br[#1]{\left[#1\right]} \def\norm#1{\left\|#1\right\|} \def\half{{1 \over 2}} \def\fL{f_{\rm L}} \def\fR{f_{\rm R}} \def\HF(#1,#2,#3,#4){{}_2F_1\!\left({#1,#2 \atop #3}; #4 \right)} \def\questionequals{= \!\!\!\!\!\!{\scriptstyle ? \atop }\,\,\,} $$

Dr Sheehan Olver
s.olver@imperial.ac.uk

Office Hours: 3-4pm Mondays, 11-12am Thursdays, Huxley 6M40
Website: https://github.com/dlfivefifty/M3M6LectureNotes

Lecture 27: Singular points of ordinary differential equations

  1. Singular points of ODEs
    • regular singular point
    • irregular singular point
    • singular points at ∞

Singular points

We will consider the first order ODE $$ u' + a(z) u = 0 $$ and the second order ODE $$ u'' + a(z) u' + b(z) u = 0 $$ where $a$ and $b$ are rational.

The simplest point is an ordinary point, where we know that $u$ is analytic:

Definition (ordinary point) $z_0$ is a ordinary point of the ODE if $a(z)$ and $b(z)$ are analytic at $z_0$.

Regular singular point

Suppose $a$ has a simple pole, and w.l.o.g. take $z_0=0$, that is: $$ a(z) = {a_{-1} \over z} + a_0 + a_1 z + a_2 z^2 + \cdots $$ where $a_0 + a_1 z + a_2 z^2 + \cdots $ is analytic. Take as an ansatz to the first order ODE $$u(z) = z^\alpha ( 1 + z v(z)) = z^\alpha + z^{\alpha+1} v$$ where $\alpha = a_{-1}$. Note that $c u(z)$ is a solution of $u'(z) = a(z) u(z)$: the choice of $1$ above makes the solution unique.

Note that $$ u'(z) = \alpha z^{\alpha-1} + z^\alpha ((\alpha+1) v + z v') $$ The ODE $u'(z) = a(z) u(z)$ becomes an inhomogenous equation: \begin{align*} z^\alpha ((\alpha+1) v + z v' - z a(z) v) &= z^{\alpha} (a(z) - {\alpha \over z} ) & \Leftrightarrow\\ z v'(z) + (\alpha+1 - z a(z)) v(z) &= a_0 + a_1 z + a_2 z^2 + \cdots \end{align*} In operator form, expressing $$ v(z) = v_0 + v_1z + v_2 z^2 = (1,z,z^2,\ldots)\begin{pmatrix}v_0\\v_1\\\vdots \end{pmatrix} $$ we have \begin{align*} z v' &=(1,z,z^2,\ldots) \begin{pmatrix} 0 \cr & 1 \cr && 2 \cr &&& 3 \cr &&&&\ddots \end{pmatrix} \begin{pmatrix}v_0\\v_1\\\vdots \end{pmatrix} \\ (\alpha+1 - z a(z)) v(z) &= (1,z,z^2,\ldots) \begin{pmatrix} 1 \cr -a_0 & 1 \cr -a_1 & -a_0& 1 \cr -a_2 & -a_1&-a_0& 1 \cr \vdots&\ddots&\ddots&\ddots&\ddots \end{pmatrix} \begin{pmatrix}v_0\\v_1\\\vdots \end{pmatrix} \end{align*} Thus the equation for our unknown coefficients becomes: $$ \begin{pmatrix} 1 \cr -a_0 & 2 \cr -a_1 & -a_0& 3 \cr -a_2 & -a_1&-a_0& 4 \cr \vdots&\ddots&\ddots&\ddots&\ddots \end{pmatrix} \begin{pmatrix}v_0\\v_1\\\vdots \end{pmatrix} = \begin{pmatrix}a_0\\a_1\\\vdots \end{pmatrix} $$ This equation induces decay in $|v_k|$ just like last lecture, hence we know its analytic with the same radius of convergence as $a$.

Here is an example: as we can see, v(z) is analytic in a disk of radius r where r is dictated by the function a.


In [2]:
r = 0.9  # radius we are solving in, can be up to but not including 1
S = Taylor(Circle(r))
α = 0.3
a = z -> α/z + 1/(1-z)
ã = Fun( z -> a(z) - α/z, S)
za = Fun( z -> z*a(z), S)
z = Fun(S)
D = Derivative() : S  S
L = z*D + (α+1-za)
v = L \ ã

# v is analytic in disk of radius r
phaseplot(v, (-1,1), (-1,1))


Out[2]:
-1.0 -0.5 0.0 0.5 1.0 -1.0 -0.5 0.0 0.5 1.0

Here is a plot of u, which shows that it has a branch cut:


In [3]:
# u has a branch 
u = z -> z^α + z^(α+1)*v(z)
phaseplot(u, (-1,1), (-1,1))


Out[3]:
-1.0 -0.5 0.0 0.5 1.0 -1.0 -0.5 0.0 0.5 1.0

We can confirm it solves the ODE:


In [4]:
up = dualpart(u(dual(0.1,1.0))) # uses DualNumbers.jl to calculate derivative

up - a(0.1)u(0.1)


Out[4]:
4.884981308350689e-15 - 7.521914234669674e-17im

Definition (first order regular singular point) $z_0$ is a regular singular point of the first order ODE if $(z-z_0)a(z)$ is analytic at $z_0$.

For second order equations, we will arrive at the following analogue:

Definition (second order regular singular point) $z_0$ is a regular singular point of the second order ODE if $(z-z_0) a(z)$ and $(z-z_0)^2b(z)$ are analytic at $z_0$.

Assume $$ a(z) = {a_{-1} \over z} + a_0 + a_1 z + \cdots $$ and $$ b(z) = {b_{-2} \over z^2} + {b_{-1} \over z} + b_0 + b_1 z + \cdots $$ To determine the right ansatz, consider $$ u'' + {a_{-1} \over z} u' + {b_{-2} \over z^2} u'' = 0 $$ which has solutions $z^\alpha$ solving the indicial equation:

Definition (indicial equation) $$ \alpha(\alpha-1) + a_{-1} \alpha + b_{-2} = 0 $$

Now take as an ansatz: $$ u(z) = z^\alpha (1 + z v(z)) = z^\alpha + z^{\alpha+1} v(z) $$ Differentiating we have \begin{align*} u'(z) &= \alpha z^{\alpha-1} + z^\alpha((\alpha+1) v + z v') \\ u''(z) &= \alpha(\alpha-1) z^{\alpha-2} + z^{\alpha-1}( \alpha(\alpha+1) v + \alpha z v' + (\alpha+1) z v' + z v' + z^2 v'' ) \\ &= \alpha(\alpha-1) z^{\alpha-2} + z^{\alpha-1}( \alpha(\alpha+1) v + 2(\alpha+1) z v' + z^2 v'' ) \end{align*} Thus $v$ solves $$ z^2 v''+ z (2(\alpha+1) + z a(z)) v' + ( \alpha(\alpha+1) + z a(z) (\alpha+1) + z^2 b(z)) v = - {\alpha(\alpha-1) \over z} - a(z) \alpha -z b(z) = \underbrace{- \alpha a_0 - b_{-1} - (\alpha a_1 + b_0) z - (\alpha a_2 + b_1) z^2 + \cdots}_{f(z) = f_0 + f_1 z + \cdots} $$

The question is when can forward substitution proceed. Let's write out the operators, first we have: $$ z^2 v'' \equiv \begin{pmatrix} 0 \cr & 0 \cr &&2 \cr &&& 6 \cr &&&&12 \cr &&&&&\ddots \end{pmatrix} $$ Then we have $$ 2 (\alpha+1) z v' \equiv \begin{pmatrix} 0\cr & 2 (\alpha+1) \cr &&&4 (\alpha+1) \cr &&&& 6(\alpha+1) \cr &&&&&8 (\alpha+1) \cr &&&&&&\ddots \end{pmatrix} $$ and $$ z^2 a(z) v' = (z a(z))(z v') \equiv \begin{pmatrix} a_{-1} \\ a_0 & a_{-1} \\ a_1 & a_0 & a_{-1} \\ a_2 & a_1 & a_0 & a_{-1} \\ \vdots &\ddots&\ddots&\ddots&\ddots \end{pmatrix} \begin{pmatrix} 0 \cr & 1 \cr &&2 \cr &&& 3 \cr &&&&4 \cr &&&&&\ddots \end{pmatrix} $$ Finally since \begin{align*} \alpha(\alpha+1) + z (\alpha+1) a(z) + z^2 b(z) &= \alpha(\alpha-1) + \alpha a_{-1} + b_{-2} + 2 \alpha + a_{-1} + ((\alpha+1) a_0 + b_{-1}) z + ((\alpha+2) a_1 + b_0) z^2 + \cdots \\ & = 2 \alpha + a_{-1} + ((\alpha+1) a_0 + b_{-1}) z + ((\alpha+2) a_1 + b_0) z^2 + \cdots \end{align*} we have $$ (\alpha(\alpha+1) + z (\alpha+1) a(z) + z^2 b(z) ) v \equiv \alpha \begin{pmatrix} 2 \alpha + a_{-1} \\ (\alpha+1) a_0 + b_{-1} & 2 \alpha + a_{-1} \\ (\alpha+1) a_1 + b_0 & (\alpha+1) a_0 + b_{-1} & 2 \alpha + a_{-1} \\ (\alpha+1) a_2 + b_1 & (\alpha+1) a_1 + b_0 & (\alpha+1) a_0 + b_{-1} & 2 \alpha + a_{-1} \\ \vdots &\ddots&\ddots&\ddots&\ddots \end{pmatrix} $$

We need to only worry about the diagonal to determine if forward substitution proceeds, in this case we have $$ \begin{pmatrix} 2 \alpha + a_{-1} \\ \times & 4 \alpha + 2 + 2a_{-1} \\ \times & \times &6 \alpha + 3a_{-1} + 6 \\ \times & \times & \times & 8 \alpha + 4a_{-1} + 12 \\ \times & \times & \times & \times & 10 \alpha + 5a_{-1} + 20\\ \vdots & \ddots & \ddots & \ddots & \ddots & \ddots \end{pmatrix} $$ Thus the first condition for this to give us a solution is that $2 (k+1)\alpha + (k+1)a_{-1} + k(k+1) \neq 0$ for all $k=1,2,\ldots$. This can be thought of another way: note that plugging $\alpha+k+1$ in to the indicial equation gives us: $$ (\alpha+k+1) (\alpha+k) + a_{-1} (\alpha+k+1) + b_{-2} = (\alpha) (\alpha-1) + a_{-1} \alpha + b_{-2}+ 2 (k+1)\alpha + (k+1)a_{-1} + k(k+1) = 2 (k+1)\alpha + (k+1)a_{-1} + k(k+1) $$ Thus forward substitution fails if and only $\alpha+k+1$ is the second root to the indicial equation for some $k=0,1,2,\ldots$.

Remark If this condition is satisfied, then $v$ is analytic near the singular point, which follows by bounding the forward elimination as before, but we'll skip this for now.

Remark If $\alpha+k+1$ is an integer, then using the ansatz $z^{\alpha+k+1}(1 + z v(z))$ proceeds without problem. To get the second solution, we need to introduce logarithms, but we'll skip this for now.

We now check the results of solving $$ z^2 v''+ z (2(\alpha+1) + z a(z)) v' + ( \alpha(\alpha+1) + z a(z) (\alpha+1) + z^2 b(z)) v = - {\alpha(\alpha-1) \over z} - a(z) \alpha -z b(z) $$ numerically.

We first choose a solutino to the an indicial equation:


In [5]:
S = Taylor()
z = Fun(S)
D = Derivative() : S

a₋₁ = -0.6
b₋₂ = 0.3

a = z -> a₋₁/z + 1/(2-z)
ã = Fun( z -> a(z) -  a₋₁/z, S)
za = Fun( z -> z*a(z), S)

b = z -> b₋₂/z^2 + 0.4/z + 1/(3-z)
zb̃ = Fun(z -> z*b(z) - b₋₂/z, S)
z²b = Fun( z -> z^2*b(z), S)

α = ((1-a₋₁) +  sqrt((a₋₁-1)^2-4b₋₂))/2
@show α
α*(α-1) + α*a₋₁ + b₋₂


α = 1.3830951894845303
Out[5]:
2.7755575615628914e-16

We can solve for v, which is analytic in a neighbourhood of 0:


In [6]:
L = z^2*D^2 + z*(2(α+1) + za)*D + (α*(α+1) + za*(α+1) + z²b)
v = L \ (-α*ã - zb̃)

phaseplot(v, (-3,3), (-3,3))


Out[6]:
-3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3

Thus one solution is:


In [7]:
u = z -> z^α + z^(α+1)*v(z)

phaseplot(u, (-3,3), (-3,3))


Out[7]:
-3 -2 -1 0 1 2 3 -3 -2 -1 0 1 2 3