由柱函数的定义(递推关系)可得:{()"2()+x2(0)=x2z(0) 有(1),(2)消去z得:J(x)=-J(x)-J-1(x) 递推关系中最大与最小下标差2,J应该满足二阶微分方程。为此 (2)式对x求导得:-z(x)+z"(x)+xZ"(x)=-2+(x)-xZ+1(x) (1)式v→+1得:(v+1)Z+(x)+xZ+1(x)=x2(x) 上两式消去Z1得:-Zx)+z(x)+xZ(x)=vZ(x)-xZ(x) 再利用(2)消去Z+1:-xz(x)+xz(x)+x2zr(x)=vvZ,(x)-xZ(x)-x2z(x) 整理:x2z(x)+xz2(x)+(x2-y2)Zx)=0此即Bese方程 反之, Bessel函数J(x)满足 Bessel方程,C;(x)=vJ(x)也必满足 Bessel方程 但:[x"Cx)=[x"vJ(x)=vx"J(x)=vx"Jy-1(x) x"C-1(x)≠xCr-1(x) 尽管Cx)=vJ(x)满足 Bessel方程,但它并非柱函数 目例2.计算积分:1=0(1-)xdx,其中4a)=0 解:利用柱函数递推关系:[x"Jx)=xyJ1-1(x)一 x"J(ax)=x"J-1(ax):x的幂次比J的阶数多1,必能写成全微分形式 =[(-x2)4ax)xdx=[(1-x)-[ x,(ax)ldx分部积分 (1-2)xJ(axl x,(ax)(2x)dx x2x)dx幂次比阶数多1 r22(axdx=J2(a) 接下去利用已知条件J(a)=0和递推关系降阶 由递推关系:J(x)+小2(x)=-J1(x)→J2(a 故:l=-J1(a) 目例3.试求:=Jx) cos xdx 解:l=xJ(x)cos x o(r)cos rdx =xJo(x)cos x- [x Jo(r)cos x-xJo(r)sin x]dx FIH: Jo(x)=1(r) x Jo(x)cos x+IxJ(x)]cosx+rJo(x)sinx/dx xJo(x)cosx+ I[xJ(x)]sin.xy' dx xJo(x)cosx +x,(r)sin x
由柱函数的定义 （递推关系） 可得： (1) v Zv(x) + x Zv ′(x) = x Zv-1(x) (2) - v Zv(x) + x Zv ′(x) = -x Zv+1(x) 有(1), (2) 消去 Zv ′ 得：Jv+1(x) = 2 v x Jv(x) - Jv-1(x) 递推关系中最大与最小下标差 2，Jv 应该满足二阶微分方程 。为此 (2) 式对 x 求导得： -v Zv ′(x) + Zv ′(x) + x Zv ″(x) = -Zv+1(x) -x Zv+1 ′ (x) (1) 式 v ⟶ v + 1 得：(v + 1) Zv+1(x) +x Zv+1 ′ (x) = x Zv(x) 上两式消去 Zv+1 ′ 得：-v Zv ′(x) + Zv ′(x) + x Zv ″(x) = v Zv+1(x) - x Zv(x) 再利用 (2) 消去 Zv+1：-v x Zv ′(x) + x Zv ′(x) + x2 Zv ″(x) = v [v Zv (x) - x Zv ′(x)] - x2 Zv(x) 整理：x2 Zv ″(x) + x Zv ′(x) + x2 - v2 Zv(x) = 0 此即 Bessel 方程。 反之，Bessel函数 Jv(x) 满足Bessel方程 ，Cv(x) = v Jv(x) 也必满足Bessel方程 ， 但：[xv Cv(x)]′ = [xv v Jv(x)]′ = v[xv Jv(x)]′ = v xv Jv-1(x) = v v - 1 xv Cv-1(x) ≠ xv Cv-1(x) 尽管 Cv(x) = v Jv(x) 满足Bessel方程 ，但它并非柱函数 。 ☺ 例 2. 计算积分： I = ∫0 1 (1 -x2 J0(a x) x x，其中 J0(a) = 0 解：利用柱函数递推关系 ： [xv Jv(x)]′ = xv Jv-1(x) ⟶ 1 a [xv Jv(a x)] ′ = xv Jv-1(a x)：x 的幂次比 J 的阶数多 1，必能写成全微分形式 I = 0 1 1 - x2 J0(a x) x x = 0 1 1 - x2 1 a [x J1(a x)] ′ x 分部积分 = 1 a 1 - x2 x J1(a x) 0 1 - 1 a 0 1 x J1(a x) (-2 x) x = 2 a 0 1 x2 J1(x) x = 2 a 0 1 1 a x2 J2(a x) ′ x = 2 a2 J2(a) 幂次比阶数多 1 接下去利用已知条件 J0(a) = 0 和递推关系降阶 ， 由递推关系 ：J0(x) + J2(x) = 2 x J1(x) ⟶ J2(a) = 2 a J1(a) 故：I = 4 a3 J1(a) ☺ 例 3. 试求：I = 0 x J0(x) cos x x 解：I = x J0(x) cos x 0 x - 0 x x [J0(x) cos x]′ x = x J0(x) cos x - 0 x [x J0 ′ (x) cos x - x J0(x) sin x] x 利用：J0 ′ (x) = -J1(x) = x J0(x) cos x + 0 x [x J1(x)] cos x + x J0 (x) [x J1(x)]′ sin x x = x J0(x) cos x + 0 x {[x J1(x)] sin x}′ x = x J0(x) cos x + x J1(x)sin x 6 z14a.nb
14a. nb Clear [x]i Integrate[ BesselJ[0,x]c。s【x],x] FunctionExpand [% x HypergeometricPEQ BesselJ [0, x] Cos [x]+x BesselJ [l, x sin[x] 目例4.试证明:6(x)+2∑(x)=1 解:由生成函数 (x)a 类似地:c(r1-川|=Jm()厂m级数绝对收敛 两式相乘:左边=1=右边= Jm(x)Jn(x)A k ∑ JA()Jk+dx)/ 展开的唯一性:t的零次幂系数为1 6(x)+2)(x)=1 展开的唯一性:t的1次幂系数为0 J(x)Jk+1(x)=0 其实:A(x)4x1(x)=Sk(x)4t(x)+A(x)1(x),第二项令k=-k 并利用Jk=(-14:=)Jx)k+(x)+)(-)(x)(-1-1k-1(x) =SJ(x)k(x)-S(x)k1()=0 Sum[BesselJ [n, x] BesselJ[n+l, x], [n,-10001, 100001] Q数值计算 目例5.试设计一种算法计算 Bessel函数J(x),n为整数。 2 解:由递推关系:J-1(x)+J+(x)=-J(x)→Jm+1(x)=-J(x)-Jm1(x)
Clear[x]; Integrate[BesselJ[0, x] Cos[x], x] FunctionExpand[%] x HypergeometricPFQ 1 4 , 3 4 , 1 2 , 1, 3 2 , -x2 x BesselJ[0, x] Cos[x] + x BesselJ[1, x] Sin[x] ☺ 例 4. 试证明：J0 2(x) + 2 n=1 ∞ Jn 2( x) = 1 解：由生成函数 ：exp x 2 t - t -1 = n=-∞ ∞ Jn(x) t n 类似地： exp x 2 t -1 - t = m=-∞ ∞ Jm(x) t -m 级数绝对收敛 两式相乘 ：左边 = 1 = 右边 = n=-∞ ∞ m=-∞ ∞ Jm(x) Jn(x) t n-m 令 m = k n - m = l ⟹ 1 = l=-∞ ∞ k=-∞ ∞ Jk(x) Jk+l(x) t l 展开的唯一性 ：t 的零次幂系数为 1 ⟹ k=-∞ ∞ Jk(x) Jk(x) = 1 J-n(x) = (-1)n Jn(x) J0 2(x) + 2 n=1 ∞ Jn 2( x) = 1 展开的唯一性 ：t 的 1 次幂系数为 0 ⟹ k=-∞ ∞ Jk(x) Jk+1(x) = 0 其实： k=-∞ ∞ Jk(x) Jk+1(x) = k=0 ∞ Jk(x) Jk+1(x) + k=-∞ -1 Jk(x) Jk+1(x), 第二项令 k = -k′ 并利用 J-k = (-1)k Jk： = k=0 ∞ Jk(x) Jk+1(x) + k=∞ 1 (-1)k′ Jk′(x) (-1)k′-1 Jk′-1(x) = k=0 ∞ Jk(x) Jk+1(x) - k=0 ∞ Jk(x) Jk+1(x) = 0 Sum[BesselJ[n, x] BesselJ[n + 1, x], {n, -10 001, 10 000}] 0 数值计算 ☺ 例 5. 试设计一种算法计算Bessel函数 Jn(x)，n 为整数。 解：由递推关系 ：Jv-1(x) + Jv+1(x) = 2 v x Jv(x) ⟶ Jn+1(x) = 2 n x Jn(x) - Jn-1(x) z14a.nb 7
可否像 Legendre多项式那样迭代:Pn1(x)= xp(x) Pr-I(r) 即由Jx),J1(x)求出2(x),以此类推?No!以x=1为例 正向迭代: 先利用番 Mathematica代码先算出Jo()及J(1) nn(x) 再由迭代关系:Jn+1(x)=-Jm(x)-Jl-1(x)一直计算到/2s(1) 记由此算法计算得到的Bese函数为:J(x) 令误差mn(x)=Logn (x)-J(x) 右图画出x=1时n(x)与n的关系 由图可见,这种算法在n=20时n20~30,意味着相对误差达10°,结果完全没有意义。 失败原因?看A1)=Lg/a)--( B,(r) Jm-(r) 分子是利用递推关系计算J+(x)时的两项之差 分母是分子所含两项中的一项 °oo 右图画出x=1时B(x)与n的关系 由图可见,这种算法在n=10时B10~2,意味着在从Jx)与o(x)计算1(x)丢失了两位有效数字 注:假设有一个有7位有效数字的数a=1.234567且 相当于b≈1.234567 101≈1.22234(这里保留7位有效数字) a,b相减的运算a-b=0.012223仅有5位有效数字,丢失了2位有效数字 即:当(a-b)/b=10-2时,a-b的值相比于a或b丢失了2位有效数字 从图还可以看到,βn(x)随n单调下降,意味着这种丢失有效数字的现象随着迭代愈发严重 因此,我们不可能期望在做了十几次迭代后,还能得到精确的值。——一正向迭代是不稳定的 到这种情况,通常采用逆向迭代法,即由n较大时的函数值,逆向求解n较小时的函数值 先利用番 Mathematica代码先算出J3o(1)及小31(1) 再由选代关系:J(m)25从(1…1(x,一直计算到b(1) 由此算法计算得到的Bee函数记为:Jn(x) -14.0 令相对误差nx)=Log 3(x)-J(x) Jn(x) 右图画出x=1时nx)与n的关系 由图可见,这种算法在保证结果的相对误差都在10-14逆向迭代是稳定的 如何取起始值 给定一个x值,由J(x)的性质可知,logx n(x)随n增大而下降, lim J(x)=0 右图为x=1时M(x随n的变化 可看出随n增大,n(x迅速下降
可否像Legendre多项式那样迭代 ：Pn+1(x) = (2 n + 1) (n + 1) x Pn(x) - n (n + 1) Pn-1(x) 即由 J0(x), J1(x) 求出 J2(x), 以此类推 ？No！ 以 x = 1 为例。 正向迭代 ： 先利用 Mathematica 代码先算出 J0(1) 及 J1(1) 再由迭代关系 ：Jn+1(x) = 2 n x Jn(x) -Jn-1(x) 一直计算到 J25(1) 记由此算法计算得到的Bessel函数为 ：J n(x) 令误差 ηn(x) = Log10 J n(x) - Jn(x) Jn(x) 右图画出 x = 1 时 ηn(x) 与 n 的关系 ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ 5 10 15 20 n -10 10 20 30 40 ηn(x) 由图可见，这种算法在 n = 20 时 η20 ∼ 30, 意味着相对误差达 1030，结果完全没有意义 。 失败原因 ？看 βn(x) = Log 2 n x Jn(x) - Jn-1(x) Jn-1(x) 分子是利用递推关系计算 Jn+1(x) 时的两项之差 分母是分子所含两项中的一项 右图画出 x = 1 时 βn(x) 与 n 的关系 ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ ◦ 5 10 15 20 n -3 -2 -1 0 1 βn(x) 由图可见，这种算法在 n = 10 时 β10 ∼ 2, 意味着在从 J 9(x) 与 J 10(x) 计算 J 11(x) 丢失了两位有效数字 。 注：假设有一个有7位有效数字的数 a = 1.234567 且： (a - b) b = 10-2; 相当于 b = 1.234567 1.01 = 1.22234 （这里保留7位有效数字 ）， a, b 相减的运算 a - b = 0.012223 仅有5位有效数字 ，丢失了 2 位有效数字 即：当 (a - b)/b = 10-2 时，a - b 的值相比于 a 或 b 丢失了 2 位有效数字 。 从图还可以看到 ， βn(x) 随 n 单调下降，意味着这种丢失有效数字的现象随着迭代愈发严重 。 因此，我们不可能期望在做了十几次迭代后 ，还能得到精确的值 。—— 正向迭代是不稳定的 逆向迭代 ： 遇到这种情况 ，通常采用逆向迭代法 ，即由 n 较大时的函数值 ，逆向求解 n 较小时的函数值 。 先利用 Mathematica 代码先算出 J30(1) 及 J31(1) 再由迭代关系 ：Jv-1(x) = 2 v x Jv(x) - Jv+1(x)，一直计算到 J0(1) 由此算法计算得到的Bessel函数记为 ：J n(x) 令相对误差 ηn(x) = Log J n(x) - Jn(x) Jn(x) 右图画出 x = 1 时 ηn(x) 与 n 的关系 ◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦ ◦◦ ◦ ◦ 0 5 10 15 20 25 n -14.5 -14.0 -13.5 -13.0 ηn(x) 由图可见，这种算法在保证结果的相对误差都在 10-14 —— 逆向迭代是稳定的 如何取起始值 ？ 给定一个 x 值，由 Jn(x) 的性质可知 ， Jn(x) 随 n 增大而下降 ，lim n∞Jn(x) = 0 右图为 x = 1 时 Jn(x) 随 n 的变化 可看出随 n 增大，Jn(x) 迅速下降 ◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦◦ 5 10 15 20 n -20 -15 -10 -5 0 log(Jn(x)) 8 z14a.nb
z14a. nb 因此,若要计算某个m值的Jm(xo) ,--+1(xo)=0 再由 Bessell函数的逆向迭代关系 求出所有阶数的值:Jn(xo),n=0,1,2,…,m 尽管2(0)≠m(x0,但m(0)女 =a与n无关 因而从J(x0)到Jxo)只需除以一个与n无关的常数,也就是重整( renormalize) 如何重整? 利用Bes函数的生成函数 cr=∑4(Pm 令 J(x)+2)J2k(x)=1 而:S=30x)+2)2(x) O,(x) 给出 Bessel函数 由图可见,这种算法在保证结果的相对误差<10-13(图中取nmax=150) 文件z14-02包含这些迭代计算的备 Mathematica代码 也可利用:56(x)+2(x)=1重整 当然,也可利用其它算法计算Jo(x),利用、Jo(x)的值来重整 143 Bessel函数的其他性质 Q微商表示 类似于 Legendre多项式的 Rodrigues表示,整数阶Bess!函数,也可以表示为微商形式。 出发点依旧是柱函数的两个基本递推关系 对整数阶 Bessel函数:利用柱函数的递推关系式 令:U(x)≡ d Un(x) d[x"Jn(x) =-x"n+(x)=-x Un(r) 即:Um+1(x) [Un(x)l r dx 进一步化为:Jn(x)=(-1)x 1 d y 对 Neumman函数也有类似公式,以为 Neumman函数也是柱函数 即:Na()=(-yx2 对半奇数阶Besl函数,先化为球 Bessel函数,再推导
因此，若要计算某个 m 值的 Jm(x0) 可取 nmax ≫ m, 令： J nmax+1(x0) = 0 J nmax (x0) = ε ∼ 10-100 再由Bessel函数的逆向迭代关系 ： J n-1(x0) = 2 n x J n(x0) - J n+1(x0) 求出所有阶数的值 ：J n(x0), n = 0, 1, 2, … m 0 10 20 30 40 50 60 70 n 220 240 260 280 300 rn(x) 尽管 J n(x0) ≠ Jn(x0)，但 rn(x0) = J n(x0) Jn(x0) = α 与 n 无关， 因而从 J n(x0) 到 Jn(x0) 只需除以一个与 n 无关的常数 ，也就是重整 （renormalize） 如何重整 ？ 利用Bessel函数的生成函数 (x/2) (t-t-1) = n=-∞ ∞ Jn(x) t n 令 t = 1 ⟶ J0 (x) + 2 k=1 ∞ J2 k(x) = 1 而：S = J 0(x) + 2 k=1 ∞ J 2 k(x) 则： J n(x) S 给出Bessel函数 Jn(x) 的值。 0 10 20 30 40 50 60 70 n -15.5 -15.0 -14.5 -14.0 -13.5 -13.0 ηn(x) 由图可见 ，这种算法在保证结果的相对误差 < 10-13 （图中取 nmax = 150） 文件 z14-02 包含这些迭代计算的 Mathematica代码 也可利用 ：J0 2(x) + 2 n=1 ∞ Jn 2( x) = 1 重整 当然，也可利用其它算法计算 J0(x)，利用 J0(x) 的值来重整 。 14.3 Bessel 函数的其他性质 微商表示 类似于Legendre多项式的 Rodrigues 表示，整数阶Bessel函数，也可以表示为微商形式。 出发点依旧是柱函数的两个基本递推关系。 对整数阶Bessel函数：利用柱函数的递推关系式 令：Un(x) ≡ Jn(x) xn ⟶ Un(x) x = [x-n Jn(x)] x = -x-n Jn+1(x) = -x Un+1(x) 即：Un+1(x) = - 1 x x [Un(x)] = - 1 x x [Un(x)] = (-1)n+1 1 x x n+1 [U0(x)] 进一步化为 ：Jn(x) = (-1)n xn 1 x x n J0(x) 对Neumman函数也有类似公式 ，以为Neumman函数也是柱函数 即：Nn(x) = (-1) n xn 1 x x n N0(x) 对半奇数阶Bessel函数，先化为球Bessel函数，再推导。 z14a.nb 9
10 z14anb x'J(x)=x"J-1(x) I**l n(x)=r"+lim-1( xJ(x)]=-x-"J+1( Ix" jn(x)=-x"n+I() 令:U(x)≡ Jn(x) d Un(x) drain(x) -x"n+I(x)=-x Un+l(x) 即:Un+1(x)= Un(x)]=(-1)m Uo(x) 进一步化为:(x)=(-1)x Jo(x)比较:4(x)=(-/d 对第二类球Bese函数也有类似公式。 J(x)无法表示为初等函数,(x)是否可以表为初等函数? 台k++)(2)中取r=1 在级数表示:J()=(-1) (-1) 2 =(yx=如x,其中利用了:(=V,2k(2k+1)=(2k+D x台d(2k+1)! J1(x)-J1(x) -N1(x) 利用了M( COS vT)J-J smnv丌 k!rk+1 (-1)x2k k!22k+1 xk!24(2k-1) =-—,其中利用了 Vr,2k!(2k-1)!1=(k)! sinx n(x)=(-1y J(x)=(-1)x2 Rayleigh公式 nn(x)=(-1)2x2 n(x)=-(-1yx 目例1.在§13.1中,我们讨论了 Jacobl-Anger公式,即:把平面波展开为柱面波之叠加。 试证明 Rayleigh方程,即:将平面波展开为球面波之叠加 exp[ i y]=)i(21+ 1)in(kr)P (cos y) Rayleigh方程 比较: exp[i k pcos=)pJkp)cn6 Jacobi- Anger公式
jn(x) = π 2 x J n+1 2 (x) 由： x [xv Jv(x)] = xv Jv-1(x) x [x-v Jv (x)] = -x-v Jv+1(x) 令：v=n+1/2 x xn+1 jn(x) = xn+1 jn-1(x) x [x-n jn(x)] = -x-n jn+1(x) 令：Un(x) ≡ jn(x) xn ⟶ Un(x) x = [x-n jn(x)] x = -x-n jn+1(x) = -x Un+1(x) 即：Un+1(x) = - 1 x x [Un(x)] = (-1)n+1 1 x x n+1 [U0(x)] 进一步化为 ：jn(x) = (-1) n xn 1 x x n j0 (x) 比较： Jn(x) = (-1)n xn 1 x x n J0(x) 对第二类球Bessel函数也有类似公式 。 J0(x) 无法表示为初等函数 ， j0(x) 是否可以表为初等函数 ？ 在级数表示 ： Jv(x) = k=0 ∞ (-1)k k ! Γ(k + v + 1) x 2 2 k+v 中取 v = 1 2 ，则 j0(x) = π 2 x J 1 2 (x) = π 2 k=0 ∞ (-1)k k ! Γ k + 1 + 1 2 x 2 2 k = π 2 k=0 ∞ (-1)k x2 k k ! 22 k Γ k + 1 + 1 2 = π 2 k=0 ∞ (-1)k x2 k k ! 22 k 2-k-1 (2 k + 1)!! Γ 1 2 = 1 x k=0 ∞ (-1) k x2 k+1 (2 k + 1)! = sin x x , 其中利用了 ：Γ 1 2 = π , 2k k ! (2 k + 1)!! = (2 k + 1)! n0(x) = π 2 x N1 2 (x) = cos π 2 J 1 2 (x) - J -1 2 (x) sin π 2 , 利用了 Nv = (cos v π) Jv - J-v sin v π = - π 2 x J -1 2 (x) = - π 2 x k=0 ∞ (-1)k k ! Γ k + 1 - 1 2 x 2 2 k-1/2 = - π x k=0 ∞ (-1)k x2 k k ! 22 k Γ k + 1 - 1 2 = - π x k=0 ∞ (-1)k x2 k k ! 2k (2 k - 1)!! Γ 1 2 = -1 x k=0 ∞ (-1)k x2 k (2 k)! = - cos x x , 其中利用了 ：Γ 1 2 = π , 2k k ! (2 k - 1)!! = (2 k)! 从而： jn(x) = (-1)n xn 1 x x n j0 (x) = (-1)n xn 1 x x n sin x x nn(x) = (-1)n xn 1 x x n n0(x) = -(-1)n xn 1 x x n cos x x —— Rayleigh 公式 ☺ 例 1. 在§ 13.1中，我们讨论了 Jacobi-Anger 公式，即：把平面波展开为柱面波之叠加。 试证明 Rayleigh方程，即：将平面波展开为球面波之叠加。 exp[ k r cos γ] = l=0 ∞ l (2 l + 1) jl(k r) Pl(cos γ) —— Rayleigh 方程 (1.5) 比较：exp[ k ρ cos θ] = n=-∞ ∞ n Jn(k ρ) n θ —— Jacobi - Anger 公式 10 z14a.nb