的特解进行讨论并用Matlab程序实现绘制特解的变化图像

1、、前言： 到目前为止已经有不少人对Bessel方程进行了讨论，但对Bessel方程的特定边界条件下的特解情况的讨论却很少。本文主要是针对一类零阶Bessel方程在如下边界条件的特解进行讨论并用Matlab程序实现绘制特解的变化图像。Matlab是类似于Fortran和C的一种语言，虽然很难简短描述其特点，但对于科学计算而言，其主要特色有：1) 易于编程2) 整数、实数、复数之间的统一性3) 高精度及扩充的数值范围4) 综合性的数学工具库5) 包含图形用户界面在内的功能完备的图形工具6) 与传统编程语言的接口7) Matlab程序的可移植性本文便是利用的强大图形工具绘制出方程特解的图像，以便形象

2、地描述Bessel方程特解的振荡性、衰减性、周期性等。、数学模型的求解及分析： 数学模型：1.1 方程通解的求解过程 首先对方程 (1)进行化简。引入，则有 = (2)=将方程(1)两边同时乘以后，经整理后得: (3)将方程(2)代入到方程(3)有: (4)对比Bessel标准方程（这里是非负常数，不一定是正整数）知，方程（4）其实就是一个的Bessel方程；即零阶的Bessel方程。根据参考文献158页定理11:若方程 中系数具有这样的性质,即和均能展成的幂级数,且收敛区间为,则方程有形如 即 的特解,这里,是一个特定的常数。级数也以为收敛区间。将方程（4）改写成 易见，它满足定理11的条件

3、，且, , 按展开成幂级数,它的收敛区间为，由定理11，方程有形如 (5)的解.这里。而和是待定常数。将（5）代如（4）有 按t的同幂次项整理有：令各项系数等于零，得一系列的代数方程： (6)因为，故从(6)的第一个方程解得：考虑时方程(4)的一个特解.这时总可以从(6)中逐个地确定所有的系数。把代入（6）得到 k=2,3,对讨论奇偶两种情况，分别有 k=1,2,3,. 从而求得 k=1,2,3,. k=1,2,3,将各代入方程(5)得到方程(4)的一个解 (7)不妨令 <2>则(7)变为 根据参考文献P404 (7.9): (n=0时,右边第二项不存在)可以得到 (8)当为实数时

4、, 是一个衰减振荡函数, 在有奇性。 、都是Bessel方程的解，且其中任意两个都是线性无关的。所以Bessel方程（4）的通解可写成： (9)1.2 在左右边界条件下方程特解的求解过程1.2.1左边界条件根据Bessel函数性质有 即有 (12)现在将(8)和(12)代入(11)有 (13) 当右边界条件为时有 =0 (14)将(13)和(14)组成方程组联合求解和有 (15)代入（13）得： 即将(16)代入(15)有 (17)将A和B代入（9）有当右边界条件是Y(R)=0时有令有 并且 (18)将方程（13）和（18）组成方程组求解A和B由（18）得： (19)代入(13)有即将(19)

5、代入(20)将A和B代入（9）有、Matlab程序实现2.1 当右边界条件为时 用Matlab对方程在边界条件下根据a的变化来分析方程的解程序代码a=-10,-6,-4,0,6,10,u=40,R=13; for i=1:6 title('Bessel函数在a变化时的情形') xlabel('自变量X') ylabel('函数值Y') b=1-a(i);t=sqrt(u); Q=60;B=besselj(1,t*R)*Q/(besselj(1,t*R)*(a(i)*bessely(0,t)-b*t*bessely(1,t)-bessely(1,t

6、*R)*(a(i)*besselj(0,t)-b*t*besselj(1,t);A=-bessely(1,t*R)*Q/(besselj(1,t*R)*(a(i)*bessely(0,t)-b*t*bessely(1,t)-bessely(1,t*R)*(a(i)*besselj(0,t)-b*t*besselj(1,t);x=0:0.05:50;y=A*besselj(0,x)+B*bessely(0,x);switch icase 1,plot(x,y,'r-');case 2,plot(x,y,'g:');case 3,plot(x,y,'b-.&

7、#39;); case 4,plot(x,y,'y-');case 5,plot(x,y,'m:');case 6,plot(x,y,'c-.'); end;hold on;end;hold off;legend ('a=-10','a=-6','a=-4','a=0','a=6','a=10',0);形成方程解的变化图像:图像分析1、 从方程解的变化图像中很容易看出Bessel方程（1）在a变化时y随x的增加具有衰减振荡性，即y的振幅随x值增大而减小

8、。2、 并且方程解y随x的变化是成周期性的，变化的周期T的大小跟a的大小无关。3、 方程（1）在a变化时，方程解y随x的变化图像有固定交点，并且交点位置与a的大小无关。用Matlab对方程在边界条件下根据R的变化来分析方程的解程序代码a=6,u=50,R=30,40,100,130,170,200; for i=1:6 title('bessel函数在R变化时的情形') xlabel('自变量X') ylabel('函数值Y') b=1-a;t=sqrt(u); Q=60;B=besselj(1,t*R(i)*Q/(besselj(1,t*R(i

9、)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(1,t*R(i)*(a*besselj(0,t)-b*t*besselj(1,t);A=-bessely(1,t*R(i)*Q/(besselj(1,t*R(i)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(1,t*R(i)*(a*besselj(0,t)-b*t*besselj(1,t);x=0:0.05:50;y=A*besselj(0,x)+B*bessely(0,x);switch icase 1,plot(x,y,'r-');case 2,plo

10、t(x,y,'g:');case 3,plot(x,y,'b-.'); case 4,plot(x,y,'y-');case 5,plot(x,y,'m:');case 6,plot(x,y,'c-.'); end;hold on;end;hold off;legend ('R=30','R=40','R=100','R=130','R=170','R=200',0)形成方程解的变化图像:图像分析1、从方程解的变化图像中

11、很容易看出Bessel方程（1）在R变化时y随x的增加具有衰减振荡性，即y的振幅随x值增大而减小。2、并且方程解y随x变化是成周期性的，变化的周期T的大小跟R的大小无关。3、方程（1）在R变化时，方程解y随x的变化图像有固定交点，并且交点位置与R的取值无关。用Matlab对方程在边界条件下根据的变化来分析方程的解程序代码a=33,R=135,u=30,40,50,70,80,100; for i=1:6 title('bessel函数在u变化时的情形') xlabel('自变量X') ylabel('函数值Y') b=1-a; t=sqrt(u(

12、i); Q=60;B=besselj(1,t*R)*Q/(besselj(1,t*R)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(1,t*R)*(a*besselj(0,t)-b*t*besselj(1,t);A=-bessely(1,t*R)*Q/(besselj(1,t*R)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(1,t*R)*(a*besselj(0,t)-b*t*besselj(1,t);x=0:0.05:50;y=A*besselj(0,x)+B*bessely(0,x);switch icase 1

13、,plot(x,y,'r-');case 2,plot(x,y,'g:');case 3,plot(x,y,'b-.'); case 4,plot(x,y,'y-');case 5,plot(x,y,'m:');case 6,plot(x,y,'c-.'); end;hold on;end;hold off;legend ('u=30','u=40','u=50','u=70','u=80','u=100'

14、;,0)形成方程解的变化图像:图像分析1、 从方程解的变化图像中很容易看出Bessel方程（1）在变化时y随x的增加具有衰减振荡性，即y的振幅随x值增大而减小。2、 并且方程解y随x的变化是成周期性，变化的周期T的大小跟的大小无关。3、 方程（1）在变化时，方程解y随x的变化图像无固定交点，可知交点位置与值有关。2.2 当右边界条件是Y(R)=0时 用Matlab对方程在边界条件下根据a的变化来分析方程的解程序代码a=-10,-6,-4,0,6,10,u=40,R=13; for i=1:6 title('bessel函数在a变化时的情形') xlabel('自变量X&

15、#39;) ylabel('函数值Y') b=1-a(i);t=sqrt(u); Q=60;B=besselj(0,t*R)*Q/(besselj(0,t*R)*(a(i)*bessely(0,t)-b*t*bessely(1,t)-bessely(0,t*R)*(a(i)*besselj(0,t)-b*t*besselj(1,t);A=-bessely(0,t*R)*Q/(besselj(0,t*R)*(a(i)*bessely(0,t)-b*t*bessely(1,t)-bessely(0,t*R)*(a(i)*besselj(0,t)-b*t*besselj(1,t);x

16、=0:0.05:50;y=A*besselj(0,x)+B*bessely(0,x);switch icase 1,plot(x,y,'r-');case 2,plot(x,y,'g:');case 3,plot(x,y,'b-.'); case 4,plot(x,y,'y-');case 5,plot(x,y,'m:');case 6,plot(x,y,'c-.'); end;hold on;end;hold off;legend ('a=-10','a=-6',&

17、#39;a=-4','a=0','a=6','a=10',0);形成方程解的变化图像:分析图像1、 从方程解的变化图像中很容易看出Bessel方程（1）在a变化时y随x的增加具有衰减振荡性，即y的振幅随x值增大而减小。2、 并且方程解y随x的变化是成周期性，变化的周期T的大小跟a的大小无关。3、 方程（1）在a变化时，方程解y随x的变化图像有固定交点，并且交点位置与a的大小无关。2.2.2 用Matlab对方程在边界条件下根据R的变化来分析方程的解程序代码a=6,u=50,R=30,40,100,130,170,200; for i=1:

18、6 title('bessel函数在R变化时的情形') xlabel('自变量X') ylabel('函数值Y') b=1-a;t=sqrt(u); Q=60;B=besselj(0,t*R(i)*Q/(besselj(0,t*R(i)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(0,t*R(i)*(a*besselj(0,t)-b*t*besselj(1,t);A=-bessely(0,t*R(i)*Q/(besselj(0,t*R(i)*(a*bessely(0,t)-b*t*bessely(1,t)-b

19、essely(0,t*R(i)*(a*besselj(0,t)-b*t*besselj(1,t);x=0:0.05:50;y=A*besselj(0,x)+B*bessely(0,x);switch icase 1,plot(x,y,'r-');case 2,plot(x,y,'g:');case 3,plot(x,y,'b-.'); case 4,plot(x,y,'y-');case 5,plot(x,y,'m:');case 6,plot(x,y,'c-.'); end;hold on;end

20、;hold off;legend ('R=30','R=40','R=100','R=130','R=170','R=200',0)形成方程解的变化图像:分析图像1、 从方程解的变化图像中很容易看出Bessel方程（1）在R变化时y随x的增加具有衰减振荡性，即y的振幅随x值增大而减小。2、 并且方程解y随x的变化是成周期性的，变化的周期T的大小跟R的大小无关。3、 方程（1）在R变化时，方程解y随x的变化图像有固定交点，可知交点位置与R值无关。2.2.3 用Matlab对方程在边界条件下根据的变化来

21、分析方程的解程序代码a=33,R=135,u=30,40,50,70,80,100; for i=1:6 title('bessel函数在u变化时的情形') xlabel('自变量X') ylabel('函数值Y') b=1-a; t=sqrt(u(i); Q=60;B=besselj(0,t*R)*Q/(besselj(0,t*R)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(0,t*R)*(a*besselj(0,t)-b*t*besselj(1,t);A=-bessely(0,t*R)*Q/(besse

22、lj(0,t*R)*(a*bessely(0,t)-b*t*bessely(1,t)-bessely(0,t*R)*(a*besselj(0,t)-b*t*besselj(1,t);x=0:0.05:50;y=A*besselj(0,x)+B*bessely(0,x);switch icase 1,plot(x,y,'r-');case 2,plot(x,y,'g:');case 3,plot(x,y,'b-.'); case 4,plot(x,y,'y-');case 5,plot(x,y,'m:');case 6,plot(x,y,'c-.'); end;hold on;end;hold off;legend ('u=30','u=40','u=50',