氢原子径向波函数数值求解 | 寒山烟雨
现在的位置: 首页 > 小试身手 > 正文

氢原子径向波函数数值求解

2012年06月03日 小试身手 ⁄ 共 727字 ⁄ 字号 评论 8 条 ⁄ 阅读 11,059 views 次

氢原子镜像波函数方程(化为原子单位制)为:

\left( -\frac{d^2}{2dr^2}+\frac{l(l+1)}{2r^2}-\frac{1}{r} \right) P(r)=EP(r)

其边界条件为:P(0)=0,P(∞)=0,l为轨道量子数,l=0,1,2,3…分别对应s轨道,p轨道,d轨道,f轨道。。。

本意使用Nomerov方法,数值求解径向波函数,但是对于r=0处有奇点,此处采用差分法求解这个二阶线性非齐次微分方程的本征值问题。Nomerov有五阶精度,而差分只有二阶精度,出现的问题稍后再说。

由于有无穷这个奇点,此处去取无穷等于40作为截短点(l=0时,r=40个玻尔半径时波函数已趋近于零)。利用打靶法确定本征能量E。
程序:(此处取l=1,求p轨道基态波函数,程序未做归一化处理)

 Fortran | 
 
 copy code |
?

01
Module Number
02
INTEGER,PARAMETER :: N=10000
03
INTEGER,PARAMETER :: K=5000
04
INTEGER,PARAMETER :: L=1
05
REAL(8),PARAMETER :: h=400.0d0/N
06
INTEGER :: i,j
07
END
08
 
09
PROGRAM main
10
use Number
11
REAL(8) :: JF,Sum,Y(0:N),E1,E2,Yk,Yk1,EPS
12
OPEN(10,FILE='Data.txt')
13
Y(0)=0.d0
14
Y(1)=0.005d0
15
EPS=5.d0
16
Yk=0.d0
17
E1=-1.d0
18
E2=0.d0 
19
 
20
j=2
21
CALL CF(Y,E1)
22
Yk1=Y(k)
23
E=-1.d0/(2*J*j)
24
 
25
DO WHILE (ABS(E1-E2)>1.d-6)
26
    E=(E1+E2)/2.d0
27
    PRINT*,E
28
    CALL CF(Y,E)
29
    IF (ABS(Y(k)-Yk)<EPS)THEN
30
        EXIT
31
    ENDIF
32
    IF ((Y(k)-Yk)*(Yk1-Yk)<0) THEN
33
        E1=E
34
    ELSE
35
        E2=E
36
    ENDIF
37
ENDDO
38
CALL CF(Y,e)
39
 
40
Sum=JF(Y(0:K))
41
PRINT*,E !,Y(k)/SQRT(Sum)
42
DO i=0,K
43
    WRITE(10,*)i*h,Y(i)
44
ENDDO
45
END
46
 
47
SUBROUTINE CF(Y,e)
48
USE Number
49
real(8) :: r,Y(0:N)
50
DO i=1,N-1
51
    r=i*h
52
    Y(i+1)=2*Y(i)+h*h*(L*(L+1)/(r*r)-2.d0/r-2*e)*Y(i)-Y(i-1)
53
ENDDO
54
END
55
 
56
FUNCTION JF(Y)
57
USE Number
58
REAL(8) :: JF,Y(0:k)
59
JF=0.d0
60
DO i=0,k-1
61
    JF=JF+(Y(i)**2+Y(i+1)**2)*h/2.d0
62
ENDDO
63
END

运行结果如图:

E=-0.1250

明显,在25个玻尔半径之后,波函数开始上翘,主要原因在于截短点的选取有问题。至于为什么没能达到40个玻尔半径,主要原因则是差分法精度不够。而能量则基本接近解析值。

 

0
【上篇】
【下篇】

目前有 8 条留言    访客:5 条, 博主:3 条

  1. Una 2013年05月07日 上午8:12  @回复  Δ-49楼 回复
    Unknown Unknown Unknown Unknown

    I value the data on your site. Thnx!

  2. zhengheng-li 2012年06月06日 下午10:44  @回复  Δ-48楼 回复
    Firefox Firefox GNU/Linux GNU/Linux

    公式很好看。


    • 管理员
      hsyyf 2012年06月06日 下午10:56  @回复  ∇地下1层 回复
      Firefox Firefox Windows Windows

      公式是Late生成的,你的Ctex怎么样了?

      • zhengheng-li 2012年06月08日 上午12:57  @回复  ∇地下2层 回复
        Firefox Firefox Ubuntu Ubuntu

        “公式很好看。”这一条不是我发的,很奇怪


        • 管理员
          hsyyf 2012年06月08日 上午8:19  @回复  ∇地下3层 回复
          UC Browser UC Browser GNU/Linux GNU/Linux

          cdn搞得鬼,中午回去关了去。。。

  3. zhengheng-li 2012年06月06日 下午10:38  @回复  Δ-47楼 回复
    Firefox Firefox Ubuntu Ubuntu

    图是用什么画的?Fortran吗?


    • 管理员
      hsyyf 2012年06月06日 下午10:55  @回复  ∇地下1层 回复
      Firefox Firefox Windows Windows

      fortran得到数据,origin画图,matlab或者mathematica也行

  4. 冷轩信 2012年06月04日 上午9:39  @回复  Δ-46楼 回复
    Google Chrome Google Chrome GNU/Linux GNU/Linux

    南瓜说的啥……………… :shock:

给我留言

留言无头像?


×