Return to computing page for the first course APMA0330
Return to computing page for the second course APMA0340
Return to Mathematica tutorial for the first course APMA0330
Return to Mathematica tutorial for the second course APMA0340
Return to the main page for the course APMA0330
Return to the main page for the course APMA0340
Return to Part V of the course APMA0330


This secton is devoted to one of the most important differential equations---Bessel equation. Its solutions were named for Friedrich Wilhelm Bessel (1784--1846); however, Daniel Bernoulli is generally credited with being the first to introduce the concept of Bessels functions in 1732. He used the function of zero order as a solution to the problem of an oscillating chain suspended at one end. In 1764 Leonhard Euler employed Bessel functions of both zero and integral orders in an analysis of vibrations of a stretched membrane, an investigation which was further developed by Lord Rayleigh in 1878, where he demonstrated that Bessels functions are particular cases of Laplaces functions.

Bessel, while receiving named credit for these functions, did not incorporate them into his work as an astronomer until 1817. The Bessel function was the result of Bessels study of a problem of Kepler for determining the motion of three bodies moving under mutual gravita- tion. In 1824, he incorporated Bessel functions in a study of planetary perturbations where the Bessel functions appear as coefficients in a series expansion of the indirect perturbation of a planet, that is the motion of the Sun caused by the perturbing body. It was likely Lagrange’s work on elliptical orbits that first suggested to Bessel to work on the Bessel functions.

The notation Jz,n was first used by Hansen (1843) and subsequently by Schlomilch (1857) and later modified to Jn(2z) by Watson (1922). Subsequent studies of Bessel functions included the works of Mathews in 1895, “A treatise on Bessel functions and their applications to physics” written in collaboration with Andrew Gray. It was the first major treatise on Bessel functions in English and covered topics such as applications of Bessel functions to electricity, hydrodynamics and diffraction. In 1922, Watson first published his comprehensive examination of Bessel functions “A Treatise on the Theory of Bessel Functions”.

Bessel equation

The Bessel differential equation is the linear second-order ordinary differential equation given by
\[ x^2 \frac{{\text d}^2 y}{{\text d} x^2} + x\,\frac{{\text d} y}{{\text d} x} + \left( x^2 - \nu^2 \right) y(x) = 0 \qquad \mbox{or} \qquad \frac{\text d}{{\text d}x} \left( x \,\frac{{\text d} y}{{\text d} x} \right) + \left( x - \frac{\nu^2}{x} \right) y(x) = 0 , \]
where ν is a real constant, called the order of the Bessel equation. Upon substitution \( y = u\, x^{-1/2} , \) it is reduced to the equation without first derivative:
\[ u'' + \left( 1 + \frac{1 - 4\nu^2}{4x^2} \right) u(x) = 0 . \]
A Bessel equation is a special case of a confluent hypergeometric equation.
Out[1]= x/2-x3/16+x5/384
psi5=FullSimplify[phi5*Integrate[phi5^-2 *(1/x),x]]
(1/4608)x (192 - 24 x^2 + x^4) ((576 (-8 + x^2))/( x^2 (192 - 24 x^2 + x^4)) +
12 Log[x] + (-3 - 3 I Sqrt[3]) Log[-12 - 4 I Sqrt[3] + x^2] +
3 I (I + Sqrt[3]) Log[-12 + 4 I Sqrt[3] + x^2])

Plot[{psi5, BesselY[1, x]}, {x, 0.1, 6},
PlotLegends -> {psi, BesselY}, PlotStyle -> {Blue, Black}]
Plot[{phi5, BesselJ[1, x]}, {x, 0.1, 6},
PlotLegends -> {phi, BesselJ}, PlotStyle -> {Blue, Black}]

The Bessel function of the first kind of order ν:
\[ J_{\nu} (x) = \sum_{k\ge 0} \frac{(-1)^k}{k!\,\Gamma (k+\nu +1)} \left( \frac{x}{2} \right)^{2k+\nu} , \]
where \( \Gamma (z) = \int_0^{\infty} x^{z-1} e^{-x} {\text d} x \) is the gamma function.


There are two Bessel functions of the second kind of order ν: one is called the Weber function:
\[ Y_{\nu} (x) = \frac{\cos \nu \pi \,J_{\nu} (x) - J_{-\nu} (x)}{\sin \nu\pi} . \]
Its product with π is called the Neumann functions: Nν(x) = πYν(x).


The Bessel function of the third kind of order ν or Hankel functions:
\[ \begin{split} H_{\nu}^{(1)} (x) &= J_{\nu} (x) + {\bf j}\,Y_{\nu} (x) , \\ H_{\nu}^{(2)} (x) &= J_{\nu} (x) - {\bf j}\,Y_{\nu} (x) , \end{split} \]
where j is the unit vector in positive vertical direction on the complex plane ℂ. Because of the linear independence of the Bessel function of the first and second kind, the Hankel functions provide an alternative pair of solutions to the Bessel differential equation.
================================ to be modified ==============


Power series solution

The coefficients of the Frobenius solution

\[ y(x) = \sum_{k\ge 0} a_k x^{k+\sigma} \]
\[ \sum_{k\ge 0} a_k \left[ (k+\sigma )^2 - \nu^2 \right] x^{k+\sigma} + \sum_{k\ge 2} a_{k-2} x^{k+\sigma} =0. \]
The indicial equation \( \sigma^2 = \nu^2 \) has two roots \( \sigma = \pm \nu \) to which correspond two linearly independent solutions.

The solution that is finite at x = 0 can be represented by

\[ J_{\nu} (x) = \left( \frac{x}{2} \right)^{\nu} \sum_{k\ge 0} \frac{(-1)^k x^{2k}}{4^k k! \Gamma (k+\nu +1)} , \]
where Γ is the Gamma function of Euler. The above series (which converges for all real or complex x) defines the Bessel function of the first kind. If \( \nu \) is not an integer, then the function \( J_{-\nu} (x) \) is linearly independent of \( J_{\nu} (x) . \) Therefore these two functions constitute the fundamental set of solutions to Bessel's equation. When \( \nu =n \) is an integer,
\[ J_{-n} (x) = (-1)^n \,J_{n} (x) , \]
and we need to find another linearly independent solution to Bessel's equation. This can be done using reduction of order, which leads to another linearly independent solution:
\[ J_{n} (x) \,\int \frac{{\text d} x}{x \left[ J_n (x) \right]^2} . \]
Upon multiplication by a suitable constant, this expression defines either \( Y_n (x) = \pi N_n (x), \) the Weber function, or \( N_n (x) , \) the Neumann function. Heinrich Martin Weber (1842--1913) and Carl (also Karl) Gottfried Neumann (1832--1925) were German mathematicians. The Neumann function can be defined by a single expression:
\[ N_{\nu} (x) = \frac{\cos \nu\pi \,J_{\nu} (x) - J_{-\nu} (x)}{\sin \nu\pi} , \]
which is valid for all values of \( \nu . \) If \( \nu =n , \) an integer, then l'Hopital's rule yields
\begin{align*} Y_n (x) &= \pi\,N_n (x) = 2\,J_n (x) \left( \ln \frac{x}{2} + \gamma \right) - \sum_{k=0}^{n-1} \frac{(n-k-1)!}{k!} \left( \frac{x}{2} \right)^{2k-n} \\ &\quad - \sum_{k\ge 0} \frac{(-1)^k x^{2k+n}}{2^{2k+n} k! (n+k)!} \left[ 2 \left( 1 + \frac{1}{2} + \frac{1}{3} + \cdots + \frac{1}{n} \right) + \frac{1}{k+1} + \frac{1}{k+2} + \cdots + \frac{1}{k+n} \right] , \end{align*}
where \( \gamma = \Gamma '(1) \approx 0.5772156 \ldots \) is the Euler constant. Bessel functions can be expressed through the single Bessel function of order zero:
\[ J_n (x) = {\bf j}^n T_n \left( {\bf j}\,\frac{\text d}{{\text d}x} \right) J_0 (x) , \qquad n=1,2,\ldots , \]
where Tn is the Chebychev polynomial of the first kind and j is the unit vector in positive vertical direction on the complex plane. The parametric Bessel differential equation
\[ x^2 y'' + x\,y' + \left( \alpha^2 x^2 - \nu^2 \right) y =0 \]
has the general solution on \( (0, \infty ) : \)
\[ y(x) = A\,J_{\nu} (\alpha x ) + B\, Y_{\nu} (\alpha x ) \qquad \mbox{or} \qquad y(x) = A\,J_{\nu} (\alpha x ) + B\, N_{\nu} (\alpha x ) , \]
with arbitrary constants A and B. When \( \alpha = {\bf j} , \) the imaginary unit (\( {\bf j}^2 =-1 \) ), we get the modified Bessel equation
\[ x^2 y'' + x\,y' - \left( x^2 + \nu^2 \right) y =0 \qquad\mbox{or in self-adjoint form}\qquad \frac{\text d}{{\text d}x} \left( x\,\frac{{\text d} y}{{\text d} x} \right) - \left( x + \frac{\nu^2}{x}\right) =0 \]
that has two linearly independent solutions:
\begin{align*} I_{\nu} (x) &= {\bf j}^{-\nu} \,J_{\nu} \left( {\bf j} x \right) = \sum_{k\ge 0} \frac{1}{k! \,\Gamma (k + \nu +1)} \left( \frac{x}{2} \right)^{2k+\nu} , \\ K_{\nu} (x) &= \frac{\pi}{2} \,\frac{I_{-\nu}(x) - I_{\nu} (x)}{\sin \pi\nu} , \end{align*}
called the modified Bessel functions of the first and second kind. The function \( K_{\nu} (x) \) is also referred to as Macdonald's function.


  1. Dettman, J.W., Power Series Solutions of Ordinary Differential Equations, The American Mathematical Monthly, 1967, Vol. 74, No. 3, pp. 428--430.
  2. Grigorieva, E., Methods of Solving Sequence and Series Problems, Birkhäuser; 1st ed. 2016.


Return to Mathematica page
Return to the main page (APMA0330)
Return to the Part 1 (Plotting)
Return to the Part 2 (First Order ODEs)
Return to the Part 3 (Numerical Methods)
Return to the Part 4 (Second and Higher Order ODEs)
Return to the Part 5 (Series and Recurrences)
Return to the Part 6 (Laplace Transform)
Return to the Part 7 (Boundary Value Problems)