上海电机学院文理学院 上海 201306

摘 要：利用Matlab软件的griddata函数对平均场方法获取的原子核的势能面进行了插值平滑，通过绘制等高线，得到等势能面，并以核素74Ge等势能面的绘制为例，分析了74Ge的形变。

关键词：原子核;等势能面;Matlab

中图分类号：O571.21

绘制原子核的等势能面是分析核素形变的重要的手段。从系统的哈密顿量出发，利用受限的Hartree-Fock方法，[1]可以计算β自由度（描述原子核的轴对称形变）和γ自由度（描述原子核的三轴形变）下的势能。由于计算的核势能是散点图，还需通过数值方法对散点进行插值计算，得到光滑的势能面，并在β-γ极坐标平面上通过绘制等高线得到等势能面，从而进一步判断原子核的形变情况。

Matlab作为一个具有较强功能的科学计算软件，能够直接调用数值方法的计算程序，可以快速高效地实现各种数值计算过程。本文将利用Matlab软件中的griddata函数对离散的原子核势能进行插值平滑，具体说明Matlab软件在原子核等势能面绘制中的应用。

1 griddata函数简介

Matlab的griddata函数可以对散点数据进行两维和三维插值。以两维插值为例，其命令语法如下：

zq= griddata（x，y，z，xq，yq，method）

上述命令可以实现将（x，y，z）的散点数据平滑成z=f（x，y）形式的曲面，在（xq，yq）指定的查询点进行插值，并返回插入的值zq，其中xq和yq为两维查询点网格，zq则为两维数组。method为插值方法，可以选择linear、nearest、natural、cubic、linear 或 v4，其中前四个选项是基于三角剖分的插值方法，而 v4则为双调和样条插值方法，默认方法为linear。本文中，我们采用cubic选项，即基于三角剖分的三次插值方法。griddata函数的详细语法说明可见Matlab软件的帮助文件。

下面我们利用griddata函数来绘制核素74Ge的等势能面。

2 核素74Ge等势能面绘制

在原子核的壳模型中，以核素56Ni为不活跃的壳芯，74Ge的价核子可以在pf5/2g9/2的模型空间中被激发。利用受限的Hartree-Fock方法，通过JUN45有效相互作用，我们计算了β和γ自由度下，体系总角动量为0时，74Ge的势能。相关的计算方法参见文献[2]的讨论，这里不再赘述。计算得到的势能分布可见图（a）中的散点。利用griddata函数可对散点进行插值平滑，再进一步地将平滑后的势能面投影到β-γ极坐标平面上，就可以做出等势能面图。具体的计算过程见下文的Matlab程序代码和对应的注释：

r=load（'fort.30'）;%读取计算的离散的核势能

subplot（2，1，1）;

plot3（r（：，1），r（：，2），r（：，3），'r.'）;%绘制势能的三维散点图

xlabel（'＼＼it x'）; % x=β

ylabel（'＼＼it y'）; % y=β *sin γ

zlabel（'{＼＼it E}（MeV）'）; %E为势能

title（'（a）'）; hold on;

xlin=linspace（min（r（：，1）），max（r（：，1）），800）;

ylin=linspace（min（r（：，2）），max（r（：，2）），800）;

[X，Y]=meshgrid（xlin，ylin）; %将插值计算的区域网格化

Z=griddata（r（：，1），r（：，2），r（：，3），X，Y，'cubic'）; %利用griddata函数插值平滑离散点，并在网格化的查询点[X，Y]上，返回插入的值Z

mesh（X，Y，Z）; %绘制平滑后的势能面

grid on;box on;

subplot（2，1，2）;

contourf（X，Y，Z，12）; %将势能面投影到β-γ极坐标平面上，绘制等高线，得到等势能面

title（'（b）'）;h=colorbar;

set（get（h，'Title'），'string'，'{＼＼it E}（MeV）'）;hold on;

sector（0，0，0，60，0.25）; % 利用自定義的函数sector 绘制β-γ平面极坐标

function sector（x0，y0，angle1，angle2，R） %定义sector函数

t=linspace（angle1/180*pi，angle2/180*pi，1000）;

x=R*cos（t）;y=R*sin（t）;x（1001）=0;y（1001）=0;

plot（x，y，'k-'）;axis（[x0，x0+R，y0，y0+R]）;

set（gca，'YColor'，'white'）;box off;

x1=R*cos（30/180*pi）;y1=0.01+R*sin（30/180*pi）;

x2=R*cos（60/180*pi）;y2=0.01+R*sin（60/180*pi）;

xlabel（{'＼＼it ＼＼beta'}）;text（x1，y1，'{＼＼it ＼＼gamma}=30^{＼＼circ}'）;

text（x2，y2，'60^{＼＼circ}'）;text（R，0，'0^{＼＼circ}'）;

end

（a）插值平滑后的势能面

（b）对应的等势能面

执行上述代码后，得到的计算结果如图所示。在图（a）中可以看到，通过griddata函数，利用基于三角剖分的三次插值方法能够很好地平滑离散点，得到光滑连续的势能面。图（b）中，在β-γ极坐标平面上所绘制的等势能面也能够准确地描述核素74Ge的形变。当γ处于450到550的范围内时，存在势能的极小值，表明74Ge具有一定的三轴形变，这与实验上的分析是相一致的。[3]

3 结语

利用Matlab软件的griddata函数对受限的Hartree-Fock方法计算得到的核素74Ge的离散的势能面进行了基于三角剖分的三次插值平滑，绘制的等势能面能够很好地反映出74Ge的结构特征。由于Matlab语言编写简单，函数功能丰富，本文提供的方法不失为绘制原子核等势能面的一种有效途径。

参考文献：

[1]RING P，SCHUCK P.The Nuclear Many-Body Problem[M].New York：Spring-Verlag，1980：266.

[2]金华.76Ge 形变的壳模型分析.上海电机学院学报[J].2017，20（5）：301-305.

[3]SUN J J，SHI Z，LI X Q，et al.Spectroscopy of 74Ge：From soft to rigid triaxiality[J].Physics Letters B，2014，734：308-313.

作者简介：金华（1979-），男，上海电机学院教师，主要从事物理学的教学与研究工作。