氘核中子与质子的相互作用

氘核中子与质子的相互作用摘 要： 核子间相互作用力简称核力。

研究中子和质子的相互作用时，氘核是一个最简单而且有用的例子。

本文从量子力学里最常见的方势井和有心力模型出发，通过求解氘核基态波函数来分析核子之间的相互作用。

通过和现有实验数据的对比，说明模型的正确性和不足，进而采用高斯作用势，利用数值方法求解波函数。

关键词：有心力；方势井；薛定谔方程；数值解.1. 球形方势井假定核力是有心力。

在质心坐标系中，氘核基态波函数ψ满足[ℏ22μ∇2+V (r )]ψ(r )=Eψ(r ) 其中，μ=M n M p M n +M p，因为中子和质子的质量近似相等，所以μ≈M p 2。

V (r )表示质子与中子的间位能。

束缚态能量E =−B =−2.226MeV ，所以B 表示核子间的结合能。

由于基态的s 波部分只包括径向波函数，所以可以设ψ(r )=R (r )=u (r )r。

若只考虑径向部分，拉普拉斯运算满足∇2ψ(r )=1r 2ddr (r 2dψ(r )dr)。

所以，∇2[u (r )r]=1rd 2u (r )dr 2。

薛定谔方程(1)可以化简为d 2u (r )dr 2−2μV (r )ℏ2u (r )+2μEℏ2u (r )=0令α2=−2μE ℏ2，v (r )=−2μV (r )ℏ2，得到方程{d 2dr2+[−α2+v (r )]}u (r )=0 u (r )满足边界条件u (0)=0，u (∞)=0。

简单起见，假定V (r )为球形方势井，即V (r )={−V 0,r <b0,r >b图一 势能分布函数(4)(1)(2)(3)方程(3)变为{ (d 2dr 2+K 2)u (r )=0,r <b (d 2dr2−α2)u (r )=0,r >b其中，K 2=2μV 0ℏ2−α2=2μℏ2(V 0−B )。

由此，可以得到解u (r )={C 1sin (Kr ),r <bC 2e −αr ,r >b考虑到波函数的物理意义，概率密度函数只与|u (r )r|有关，所以C 1和C 2可以取正实数。

图二 约化径向波函数为了确定待定常数C 1和C 2，需运用波函数的连续性和归一化条件。

u (r )在r =b 处连续，所以C 1sin (Kb )=C 2e −αb基态波函数ψ(r )=u (r )r={C 1sin (Kr )/r,r <bC 2e −αr /r,r >b对|ψ(r )|2在全空间进行积分，即∫dφ∫sin θdθ∫|ψ(r )|2r 2dr =1∞π2π得到4π(C 12∫sin 2(Kr )dr+C 22∫e −2αr dr ∞bb)=1解出C 122[b −sin (2Kb )2K ]+C 222αe −2αb =14π结合之前得到的关系式，可以得到关于C 12和C 22的二元方程组，即{sin 2(Kb )C 12−e −2αb C 22=012[b −sin (2Kb )2K ]C 12+e −2αb 2αC 22=14π(5)(6)(7)(8)(9)(10)(11)(12)通过二阶行列式可以很容易解出C 12和C 22的值，进而求得C 1和C 2的值。

另外，由于势井的深度有限，所以 u (r )的一阶导数在r =b 处也连续。

所以有kC 1cos (Kb )=−αC 2e −αb结合式和式，不但可以确定待定系数C 1和C 2的值，还可以解出束缚态能量。

式(7)除以式(13)，得到tan (Kb )=−K/α令Kb =x ，可以得到关于x 的超越方程tan x =−xαb 。

利用图解法，在坐标系中画出函数y =tan x 和y =−x αb左右两边两个函数的图像能够求到数值解，也可以采用计算机进行数值求解。

（见附录）值得注意的是，外部波函数的梯度∇u (r )=−αC 2e −αr <0。

假设氘核只有一个束缚态，则内部波函数C 1sin (Kr )在r =b 处必定是刚刚开始下降，即Kb =π2+δ，δ为一小量。

如果B =0，则Kb =π2，此时V 0与b 之间满足V 0b 2=(π2)2ℏ22μ≈π2ℏ24M p。

对于实际的B 来说，V 0b 2比π2ℏ24M p稍大，也就是说，势井的深度与宽度之间满足一定的约束条件。

以上所讨论的模型都十分粗略和简单，我们可以结合已知的实验数据对建立的模型加以修正。

1) 氘核是由一个质子和一个中子组成的稳定核，它只有一个束缚态（这点与之前的假设一致），核自旋J =1；2) 结合能B =2.226±0.002MeV ； 3) 磁矩μD =0.857411±0.000019核磁子； 4) 电四极矩Q D =(2.735±0.014)×10−27cm 2。

由此，能够看出实验结果和模型假设之间的一些矛盾。

因为既然认为氘核处于基态，其s 波函数是球对称的，就不会有电四极矩。

如果核子间相互作用是有心力场，则轨道角动量是运动常数，氘核基态必然只能有一定的轨道角动量。

S 态不可能和其他态混起来，所以Q D ≠0表示核子间位能不可能全是有心力场。

但由于Q D 和μD −μp −μn 的数值比较小，为了估计有心力的大小，略去这些是合理的。

2. 高斯型作用势以上的讨论比较简单，当质子与中子的相互作用势变得相对复杂时，想要求得波函数的解析解是几乎不可能的。

所以以下给出氘核基态波函数的一种数值解法。

选取高斯形式的相互作用势，即V (r )=−V 0e −(r/r 0)2结合相关的文献[3]，V 0=72.194MeV ，r 0=1.484fm 。

根据式(2)得到方程u ′′(r )=Au (r )(13)(14)(15)(16)其中，A =2μℏ2[V (r )−E ]，基态能量E =−2.226MeV 。

由于氘核的半径很小，对上式的求解取r 的范围为：r ～(0,25fm)。

把r 所属的区间进行均匀的n 等分，每个区间的长度为h ，第i 个区间的起点为r i ，终点为r i+1。

将r i+1和r i−1处的约化径向波函数在r i 点进行泰勒展开，u (r i+1)=u (r i )+ℎu ′(r i )+ℎ22u ′′(r i )+ℎ36u ′′′(r i )+ℎ424u (4)(r i )+ℎ5120u (5)(r i )+ℎ6720u (6)(r i )+⋯u (r i−1)=u (r i )−ℎu ′(r i )+ℎ22u ′′(r i )−ℎ36u ′′′(r i )+ℎ424u (4)(r i )−ℎ5120u (5)(r i )+ℎ6720u (6)(r i )+⋯两式相加，得u (r i+1)−2u (r i )+ u (r i−1)=ℎ2u ′′(r i )+ℎ412u (4)(r i )+ℎ6360u (6)(r i )+⋯经过简单的推导，可以将式的微分方程化成如下：u (r i+1)≈[2+56ℎ2A (r i )]u (r i )−[1−ℎ212A (r i−1)]u (r i−1)[1−ℎ212A (r i+1)]即约化径向波函数满足的三点中央差分递推公式。

要确定这个波函数的具体数值，需要知道u (r 1)和u (r 2)两点的数值，u (r 1)=0是已知的，但u (r 2)的数值不好确定，简单起见，取一个小的整数来近似表示u (r 2)。

通过Matlab 软件，可以通过循环迭代方式求出每一点的波函数，然后画出u (r )的图像（见附录）。

3. 结论本文介绍了两种作用势模型下质子与中子的相互作用，并给出了两种模型的计算方法，分析了两种方法的合理性和优缺点，着重介绍了常用的科学计算方法。

但由于缺乏数据和相关的研究，并不能和真实的情况进行对比。

参考文献：[1] 杨立铭，于敏．原子核物理讲义[M]. 北京：北京大学出版社，2014：6-8. [2] 杨福家. 原子物理学[M]. 北京：高等教育出版社，2008：300-307.[3] 李小华，张贵清等. 氘核基态波函数的数值解[J]. 广西物理，2011，32(4)：11-13.(17)(18)附录图三函数y=tanx与y=x/2的图象数值方法求解超越方程，可以使用简单的二分法进行计算。

但利用二分法的前提是待求解的方程必须在给定的区间内单调连续，否则会得不到正确的解，所以可以先做出相应的图形，通过简单的计算大致判断根所在的区间。

本文利用C语言来编写此算法，并以一个具体的函数作为例子进行求解。

源代码如下：#include<stdio.h>#include<math.h>double fun(double x);double findroot(double a,double b);int main(void){double a=1.8,b=2.5;//选定的初值区间为[1.8,2.5]double value=findroot(a,b);printf("%g\n",value);return 0;}//待求解的方程为tanx=-x/2double fun(double x){double y;y=tan(x)+x/2.0;return y;}//二分法求解double findroot(double a,double b){double x1,x2,x0;double root=0;x1=a;x2=b;if(fun(x1)*fun(x2)>0){printf("重新选择计算区间！\n");}else{while(fabs(x1-x2)>1e-6)//计算结果的精度{x0=(x1+x2)/2.0;if(fun(x1)*fun(x0)<0){x2=x0;}else{x1=x0;}}root=(x1+x2)/2.0;}return root;}运行结果为：2.28893%---------------数值求解高斯形式相互作用势下的波函数-----------%clear;clc;n=500;r=linspace(0,25,n+1);%等间隔插入个n+1个点，将区间均分成n段h=25/n;%步长%质子质量为938MeV/c2，ℏc=197fm*MeVA=938/197^2*(-72.194*exp(-r.^2/1.484^2)+2.226);u=zeros(1,n+1);u(2)=0.0025;for i=2:nu(i+1)=((2+5/6*h^2*A(i))*u(i)-(1-h^2./12*A(i-1))*u(i-1))/(1-h^2./12*A(i+1)); endu=u/max(u);plot(r,u,'.-r'),xlabel('r/fm'),ylabel('u(r)'),set(gca,'ytick',[]);grid on,box on,title('约化径向波函数');图像如下：图四波函数图像值得说明的是，由于带入数据的精度和初始条件的设定合理性，会使得求解出的函数存在很大的偏差，原因是波函数的边界条件应该是u(0)=0，u(∞)=0。