SciPy中的Bessel函数计算
中篇中将处理的基本问题:
1. SciPy和MATLAB文档中Bessel微分方程的不同
2. SciPy中Bessel函数的功能框架及其与MATLAB的对应关系
3. 在物理问题计算中的数据维度问题
4. 微积分计算
下篇中将处理的高阶问题:
5. 球坐标Bessel函数
6. 如何使用Cython加速计算
7. 如何使用numpy避免循环
8. 如何使用并行和GIL技术
1. SciPy和MATLAB文档中Bessel微分方程的不同
非常有必要首先补充一些基本的数学物理知识,Bessel方程是描述各类和波动有关的物理现象的,在数学物理中极其重要。数值计算工具大大降低了使用这个数理工具完成各类本质上是数学求解问题的开发和研究任务,因此,我们有必要了解这些工具如何使用。怎么算了解?我们设立四个指标:
计算准确性
计算效率
适用性
可解释性
所有这些标准的订立和提升,都要求我们在使用计算工具之前,对数学本身和问题本身有深刻和正确的理解与表述。为此,我们首先从两个文档对Bessel方程的表述的不同,讲一讲这个方程的推导过程。
Bessel函数的符号在SciPy的文档里与MATLAB不同,输入/自变量为x,系数 ν \nu ν写成了 α \alpha α,第三项导致符号也不一样:
x 2 d 2 y d x 2 + x d y d x + ( x 2 − α 2 ) y = 0 x^{2} \frac{d^{2} y}{d x^{2}}+x \frac{d y}{d x}+\left(x^{2}-\alpha^{2}\right) y=0 x2dx2d2y+xdxdy+(x2−α2)y=0
MATLAB文档中的Bessel函数表达式:
z 2 d 2 y d z 2 + z d y d z − ( z 2 + ν 2 ) y = 0 z^{2} \frac{d^{2} y}{d z^{2}}+z \frac{d y}{d z}-\left(z^{2}+\nu^{2}\right) y=0 z2dz2d2y+zdzdy−(z2+ν2)y=0
Bessel函数在数学物理中,是Bessel方程的解,复变函数在Bessel微分方程上的应用,使得 α \alpha α可以是复数。Bessel微分方程在历史上有两个推导方法,其一是从亥姆霍兹方程导出:
∇ 2 u + k 2 u = 0 \nabla^{2} u+k^{2} u=0 ∇2u+k2u=0
经过曲线坐标变换得到:
1 ρ ∂ ∂ ρ ( ρ ∂ u ∂ ρ ) + 1 ρ 2 ∂ 2 u ∂ ϕ 2 + ∂ 2 u ∂ z 2 + k 2 u = 0 \frac{1}{\rho} \frac{\partial}{\partial \rho}\left(\rho \frac{\partial u}{\partial \rho}\right)+\frac{1}{\rho^{2}} \frac{\partial^{2} u}{\partial \phi^{2}}+\frac{\partial^{2} u}{\partial z^{2}}+k^{2} u=0 ρ1∂ρ∂(ρ∂ρ∂u)+ρ21∂ϕ2∂2u+∂z2∂2u+k2u=0
取 u ( ρ , ϕ , z ) = R ( ρ ) Φ ( ϕ ) Z ( z ) u(\rho, \phi, z)=R(\rho) \Phi(\phi) Z(z) u(ρ,ϕ,z)=R(ρ)Φ(ϕ)Z(z):
{ Φ ′ ′ ( ϕ ) + m 2 Φ ( ϕ ) = 0 Z ′ ′ ( z ) − ω 2 Z ( z ) = 0 ρ 2 d 2 R d ρ 2 + ρ d R d ρ + [ ( k 2 + ω 2 ) ρ 2 − m 2 ] R = 0 \left\{\begin{array}{l}{\Phi^{\prime \prime}(\phi)