内容发布更新时间 : 2024/12/23 20:37:17星期一 下面是文章的全部内容请认真阅读。
常微分方程的边值问题和本征值问题
一、问题描述
利用搜索法和弦割法，得到该常微分方程的本征值，再利用打靶法计算多个本征值。
二、解决方法 （一）搜索法
1.先随便猜测k的一个试验值,程序中令k=1 2.由Numerov算法
根据本题的条件，kn+1=kn=kn-1=k,s=0,得到yn+2,yn+1,yn间的迭代公式 令con=(k*h)^2/12
yn+2=2*(1-5*con)*yn+1/(1+con)-yn
3自己给定φ的初始条件，然后利用公式得到边界值φ(1)
4.然后以小的步长dk增加k值，这里令dk=1,每当φ(1)改变符号时，就将步长减半后倒退回来重复
5.当步长小于所要求的容许误差时终止程序，此时的k值即为所求。
（二）弦割法
1. 随便猜测两个k值，这里令k0=1,k1=2
2. 自己给定φ的初始条件，对两个k值分别利用上述公式进行迭代，得到边界值y1(1)和y2(1)。
3. 比较y1(1)和y2(1)的绝对值大小。若绝对值大，说明对应的k值距离本征值距离较远。 4. 将（k0+k1）/2赋给k2,边界值绝对值小的对应的k值保持不变，边界值绝对值大的对应k值重新定位k2的值。
5. 重复进行实验，当y1(1)和y(2)的差的绝对值小于容许误差时终止程序。此时k1的值即为所求。
当搜索法和弦割法大致求出了一个本征值后，利用打靶法，调整k值再度进行搜索，得到多个本征值，绘出其中一个本征值对应的函数图像，观察其性质。
三、程序实现 1. 搜索法
subroutine add(t,y0,y1) !利用子程序表示函数值的迭代 implicit none
real(8)::t,h,con,y0,y1,y2
integer::i,n n=10000 h=1.0/n
con=(t*h)**2/12 do i=1,n-1
y2=2*(1-5*con)*y1/(1+con)-y0 !利用Numerov算法，得到迭代公式 y0=y1 !向前迭代 y1=y2 end do return
end subroutine add program zy3 implicit none
real(8)::diffk,dk,yold,k,b0,b1 integer::s
b0=0.01 b1=0.01 s=1
k=s dk=1 diffk=0.0000001 call add(k,b0,b1)
yold=b1 do while(abs(dk)>diffk) k=k+dk b0=0.01 b1=0.01
call add(k,b0,b1)
if(yold*b1<0)then k=k-dk dk=dk/2.0 end if end do
write(*,*)k end 2. 弦割法
subroutine add(t,b) implicit none
real(8)::t,h,con,y0,y1,y2,b integer::i,n b=0 y0=0.01
取初始值，根据题目条件，令y0=y1,来保证x=0的位置
导数为0
!给定一个猜测的k值,此为搜索的初值 !给定步长
!给定步长最后达到的误差范围 !通过运行子程序，得到由初始值积到x=1时的不为0
的函数值
!开始搜索 !在k中走一步 !若果y1变号 !后退 !步长减半 !写出求得的本征值
!利用子程序表示函数值的迭代 ! y1=0.01 n=10000 h=1.0/n
con=(t*h)**2/12 do i=1,n-1
y2=2*(1-5*con)*y1/(1+con)-y0 !利用Numerov算法，得到迭代公式 y0=y1 !向前迭代 y1=y2 end do
b=abs(y1) !得到x=1处函数值的绝对值，为确定k2点的
位置做准备
return
end subroutine add program zy3
real(8)::a,k0,k1,k2,dk,m1,m2,dm integer::i,n
k0=1 !给两个启动值 k1=2 k2=0
dm=0.00000001 !表示k0和k1对应的函数值相等时允许的误差 m1=0 !此值表示k值取k0时，x=1处函数值的绝对值 m2=0 !此值表示k值取k1时，x=1处函数值的绝对值 do while (.true.) call add(k0,m1)
call add(k1,m2) !运行子程序，分别得到两个k值对应的x=1处的
函数值的绝对值
if(abs(m1-m2)
k0=k2 !当k0对应的函数值绝对值较大时，表示其离本
征值较远，而将其舍弃不用，赋k2值
k1=k1 !此时k1距离本征值较近，不用变 else
k2=(k1+k0)/2.0 !反之，舍去k1，k0不变 k0=k0 k1=k2 end if end do
write(*,*)k2 !得到本征值k end program zy3 3. 打靶法得多个本征值
subroutine add(t,y0,y1) !利用子程序表示函数值的迭代 implicit none