基于神经网络波函数的粲偶素势模型

史焕辰; 崔二亮; 周丹

doi:10.7498/aps.75.20260493

基于神经网络波函数的粲偶素势模型

西北农林科技大学理学院, 杨凌　712100

河北师范大学物理学院, 石家庄　050024

通讯作者: E-mail: erliang.cui@nwafu.edu.cn; E-mail: danzhou@hebtu.edu.cn

中图分类号: 14.40.Lb, 12.39.Pn, 02.70.Ss, 84.35.+i

Charmonium potential models based on neural network wave functions

College of Science, Northwest A&F University, Yangling 712100, China

Department of Physics, Hebei Normal University, Shijiazhuang 050024, China

Corresponding authors: E-mail: erliang.cui@nwafu.edu.cn; E-mail: danzhou@hebtu.edu.cn

MSC: 14.40.Lb, 12.39.Pn, 02.70.Ss, 84.35.+i

摘要: 基于神经网络波函数发展了一套完整的粲偶素势模型参数拟合框架, 旨在解决传统唯象势模型参数难以快速响应实验进展的问题. 该方法将神经网络波函数作为径向薛定谔方程的变分解, 通过显式构造节点因子与残差网络结构, 精确求解粲偶素体系的基态与激发态; 在此基础上, 采用最大似然估计以及拉普拉斯近似对势模型参数进行推断, 以获得最优参数及其不确定度. 应用该框架对粲偶素势模型的五个势参数进行联合优化, 优化后参数为: 正反粲夸克质量$m_{{\mathrm{c}}} = m_{\bar{{\mathrm{c}}}} = (1.4801 \pm 0.0055)\ \text{GeV}$, 强相互作用耦合常数$\alpha_{\mathrm{s}} = $$ 0.5459 \pm 0.0120$, 弦张力$\sigma_{\text{lin}} = (0.1442 \pm 0.0028)\ \text{GeV}^2$, 涂抹半径$\sigma_{\text{smear}} = (1.1047 \pm 0.0568)\ \text{GeV}$. 相比于原始参数, 优化后模型对16个已知粲偶素态的均方根偏差(RMSD)从$26.26\ \text{MeV}$降至$24.62\ \text{MeV}$. 基于优化参数预测了包括S, P, D, F波在内的粲偶素质量谱, 结果符合势模型理论预期, 预测的高角动量态为未来实验搜寻提供了理论参考.

关键词:

Abstract: This study presents a comprehensive framework for fitting the phenomenological potential model parameters of charmonium by integrating neural network wavefunctions. This approach overcomes the limitations of traditional fitting methods, which adapt slowly to new experimental data. In this method, the radial Schrödinger equation is solved variationally using a neural network wavefunction ansatz. The wavefunction incorporates an explicitly constructed node factor that accurately represents the nodal structure of excited states. It also employs a residual network architecture to enhance expressiveness and training stability. The number of nodes is specified according to the target radial excitation, while the node positions are treated as trainable parameters, enabling the neural network to optimize the excited-state wavefunctions without preassigning their nodal locations. The energy variance is used as the training objective, and automatic differentiation is employed to calculate the derivatives required in the variational optimization and subsequent parameter inference. After determining the eigenwavefunctions and eigenenergies, the potential parameters are inferred via maximum likelihood estimation with the Laplace approximation, yielding optimal parameter values and their uncertainties. In the likelihood construction, a fixed theoretical uncertainty is introduced together with the experimental errors so that the fitting procedure reflects both experimental precision and the intrinsic approximation of the potential model. The method is applied to a phenomenological charmonium potential that includes a Coulomb term, a linear confinement term, a smeared spin-spin interaction, and spin-orbit and tensor forces treated perturbatively. The fitted parameters are the charm quark mass $m_c = m_{\bar{c}} = (1.4801 \pm 0.0055)\ \text{GeV}$, the strong coupling constant $\alpha_s = 0.5459 \pm 0.0120$, the string tension $\sigma_{\text{lin}} = (0.1442 \pm 0.0028)\ \text{GeV}^2$, and the smearing radius $\sigma_{\text{smear}} = (1.1047 \pm 0.0568)\ \text{GeV}$. Compared with the original parameters, the optimized model improves the agreement with the 16 known charmonium states from the Particle Data Group. The root-mean-square deviation is reduced from 26.26 MeV to 24.62 MeV. Based on the optimized parameters, the mass spectrum of charmonium is predicted up to high excitation levels, including S, P, D, and F waves. The predicted energy level spacings and fine structures are consistent with the expected behavior of the confining potential, indicating that the optimized parameters retain the main physical features of the conventional potential model. Notably, the framework provides predictions for unobserved high-angular-momentum states such as the $n{\,}^3{\mathrm{F}}_J$ and $n{\,}^1{\mathrm{F}}_3$ series, which may serve as references for future experimental searches at facilities such as BESIII, Belle II and LHCb. The main innovation of this work is to combine neural network wavefunctions with statistical parameter inference, so that the solution of the Schrödinger equation, the fitting of potential parameters, and the quantification of parameter uncertainties are incorporated into an automated and extensible workflow. The proposed framework can be readily extended to other heavy quark systems such as bottomonium and $B_c$ mesons, as well as to more complex spectroscopic problems involving coupled-channel effects.

Key words:

粒子态

实验值

计算值 (原始参数)

计算值 (新参数)

相对偏差 (原始参数)

相对偏差 (新参数)

$ \eta_{\mathrm{c}}(1{\,}^1 {\mathrm{S}}_0) $

$ 2984.1\pm0.4 $

2982.35

2985.08

0.0586%

0.0328% $ \downarrow $

$ J/\psi(1{\,}^3 {\mathrm{S}}_1) $

$ 3096.900\pm0.006 $

3090.46

3094.57

0.2079%

0.0751% $ \downarrow $

$ \eta_{\mathrm{c}}(2{\,}^1 {\mathrm{S}}_0) $

$ 3637.8\pm0.6 $

3630.73

3637.54

0.1943%

0.0072% $ \downarrow $

$ \psi(2{\,}^3 {\mathrm{S}}_1) $

$ 3686.097\pm0.011 $

3672.28

3679.73

0.3750%

0.1727% $ \downarrow $

$ \psi(3{\,}^3 {\mathrm{S}}_1) $

$ 4039.6\pm4.3 $

4071.65

4082.67

0.7933%

1.0663% $ \uparrow $

$ \psi(4{\,}^3 {\mathrm{S}}_1) $

$ 4415\pm5 $

4406.52

4419.31

0.1921%

0.0975% $ \downarrow $

$ h_{\mathrm{c}}(1{\,}^1 {\mathrm{P}}_1) $

$ 3525.37\pm0.14 $

3515.18

3521.81

0.2889%

0.1009% $ \downarrow $

$ \chi_{{\mathrm{c}}0}(1{\,}^3 {\mathrm{P}}_0) $

$ 3414.71\pm0.30 $

3427.04

3418.48

0.3612%

0.1104% $ \downarrow $

$ \chi_{{\mathrm{c}}1}(1{\,}^3 {\mathrm{P}}_1) $

$ 3510.67\pm0.05 $

3505.71

3509.20

0.1414%

0.0419% $ \downarrow $

$ \chi_{{\mathrm{c}}2}(1{\,}^3 {\mathrm{P}}_2) $

$ 3556.17\pm0.07 $

3554.59

3567.64

0.0443%

0.3225% $ \uparrow $

$ \chi_{{\mathrm{c}}0}(2{\,}^3 {\mathrm{P}}_0) $

$ 3922.1\pm1.8 $

3867.16

3881.39

1.4009%

1.0380% $ \downarrow $

$ \chi_{{\mathrm{c}}2}(2{\,}^3 {\mathrm{P}}_2) $

$ 3922.5\pm1.0 $

3967.51

3974.43

1.1474%

1.3238% $ \uparrow $

$ \psi(1{\,}^3 {\mathrm{D}}_1) $

$ 3773.7\pm0.7 $

3785.88

3794.86

0.3228%

0.5608% $ \uparrow $

$ \psi_2(1{\,}^3 {\mathrm{D}}_2) $

$ 3823.51\pm0.34 $

3800.36

3808.69

0.6054%

0.3875% $ \downarrow $

$ \psi_3(1{\,}^3 {\mathrm{D}}_3) $

$ 3842.71\pm0.20 $

3806.08

3813.89

0.9533%

0.7499% $ \downarrow $

$ \psi(2{\,}^3 {\mathrm{D}}_1) $

$ 4191\pm5 $

4143.05

4149.04

1.1441%

1.0012% $ \downarrow $

注: 所有质量的单位均为MeV, 不同参数下的计算值均来自本文所使用的计算框架. $ \downarrow $为相对偏差下降的态, $ \uparrow $为相对偏差上升的态.

粒子态

质量/MeV

粒子态

质量/MeV

粒子态

质量/MeV

$ \boldsymbol{\eta_c}(\mathbf{1{\,}^1S_0}) $

2985.08

$ h_{\mathrm{c}}(4{\,}^1 {\mathrm{P}}_1) $

4599.04

$ \eta_{{\mathrm{c}}3}(1{\,}^1 {\mathrm{F}}_3) $

4035.37

$ \boldsymbol{J/\psi}(\mathbf{1{\,}^3S_1}) $

3094.57

$ \chi_{{\mathrm{c}}0}(4{\,}^3 {\mathrm{P}}_0) $

4548.98

$ \chi_{{\mathrm{c}}2}(1{\,}^3 {\mathrm{F}}_2) $

4040.40

$ \boldsymbol{\eta_c}(\mathbf{2{\,}^1S_0}) $

3637.54

$ \chi_{{\mathrm{c}}1}(4{\,}^3 {\mathrm{P}}_1) $

4597.01

$ \chi_{{\mathrm{c}}3}(1{\,}^3 {\mathrm{F}}_3) $

4038.76

$ \boldsymbol{\psi}(\mathbf{2{\,}^3S_1}) $

3679.73

$ \chi_{{\mathrm{c}}2}(4{\,}^3 {\mathrm{P}}_2) $

4627.69

$ \chi_{{\mathrm{c}}4}(1{\,}^3 {\mathrm{F}}_4) $

4030.69

$ \eta_{\mathrm{c}}(3{\,}^1 {\mathrm{S}}_0) $

4053.60

$ \eta_{{\mathrm{c}}2}(1{\,}^1 {\mathrm{D}}_2) $

3807.92

$ \eta_{{\mathrm{c}}3}(2{\,}^1 {\mathrm{F}}_3) $

4362.51

$ \boldsymbol{\psi}(\mathbf{3{\,}^3S_1}) $

4082.67

$ \boldsymbol{\psi}(\mathbf{1{\,}^3D_1}) $

3794.86

$ \chi_{{\mathrm{c}}2}(2{\,}^3 {\mathrm{F}}_2) $

4362.90

$ \eta_{\mathrm{c}}(4{\,}^1{\mathrm{ S}}_0) $

4396.46

$ \boldsymbol{\psi_2}(\mathbf{1{\,}^3D_2}) $

3808.69

$ \chi_{{\mathrm{c}}3}(2{\,}^3 {\mathrm{F}}_3) $

4364.60

$ \boldsymbol{\psi}(\mathbf{4{\,}^3S_1}) $

4419.31

$ \boldsymbol{\psi_3}(\mathbf{1{\,}^3D_3}) $

3813.89

$ \chi_{{\mathrm{c}}4}(2{\,}^3 {\mathrm{F}}_4) $

4360.61

$ \boldsymbol{h_c}(\mathbf{1{\,}^1P_1}) $

3521.81

$ \eta_{{\mathrm{c}}2}(2{\,}^1 {\mathrm{D}}_2) $

4168.68

$ \eta_{{\mathrm{c}}3}(3{\,}^1 {\mathrm{F}}_3) $

4657.63

$ \boldsymbol{\chi_{c0}}(\mathbf{1{\,}^3P_0}) $

3418.48

$ \boldsymbol{\psi}(\mathbf{2{\,}^3D_1}) $

4149.04

$ \chi_{{\mathrm{c}}2}(3{\,}^3 {\mathrm{F}}_2) $

4660.06

$ \boldsymbol{\chi_{c1}}(\mathbf{1{\,}^3P_1}) $

3509.20

$ \psi_2(2{\,}^3 {\mathrm{D}}_2) $

4169.32

$ \chi_{{\mathrm{c}}3}(3{\,}^3 {\mathrm{F}}_3) $

4660.33

$ \boldsymbol{\chi_{c2}}(\mathbf{1{\,}^3P_2}) $

3567.64

$ \psi_3(2{\,}^3 {\mathrm{D}}_3) $

4179.08

$ \chi_{{\mathrm{c}}4}(3{\,}^3 {\mathrm{F}}_4) $

4655.67

$ h_{\mathrm{c}}(2{\,}^1 {\mathrm{P}}_1) $

3942.71

$ \eta_{{\mathrm{c}}2}(3{\,}^1 {\mathrm{D}}_2) $

4485.30

$ \eta_{{\mathrm{c}}3}(4{\,}^1 {\mathrm{F}}_3) $

4930.43

$ \boldsymbol{\chi_{c0}}(\mathbf{2{\,}^3P_0}) $

3881.39

$ \psi(3{\,}^3 {\mathrm{D}}_1) $

4463.47

$ \chi_{{\mathrm{c}}2}(4{\,}^3 {\mathrm{F}}_2) $

4929.10

$ \chi_{{\mathrm{c}}1}(2{\,}^3 {\mathrm{P}}_1) $

3932.72

$ \psi_2(3{\,}^3 {\mathrm{D}}_2) $

4484.51

$ \chi_{{\mathrm{c}}3}(4{\,}^3 {\mathrm{F}}_3) $

4932.13

$ \boldsymbol{\chi_{c2}}(\mathbf{2{\,}^3P_2}) $

3974.43

$ \psi_3(3{\,}^3 {\mathrm{D}}_3) $

4498.95

$ \chi_{{\mathrm{c}}4}(4{\,}^3 {\mathrm{F}}_4) $

4930.71

$ h_{\mathrm{c}}(3{\,}^1 {\mathrm{P}}_1) $

4290.97

$ \eta_{{\mathrm{c}}2}(4{\,}^1 {\mathrm{D}}_2) $

4772.61

$ \chi_{{\mathrm{c}}0}(3{\,}^3 {\mathrm{P}}_0) $

4190.25

$ \psi(4{\,}^3 {\mathrm{D}}_1) $

4747.05

$ \chi_{{\mathrm{c}}1}(3{\,}^3 {\mathrm{P}}_1) $

4276.53

$ \psi_2(4{\,}^3 {\mathrm{D}}_2) $

4771.57

$ \chi_{{\mathrm{c}}2}(3{\,}^3 {\mathrm{P}}_2) $

4336.74

$ \psi_3(4{\,}^3 {\mathrm{D}}_3) $

4789.99

基于神经网络波函数的粲偶素势模型

通讯作者: E-mail: erliang.cui@nwafu.edu.cn;

通讯作者: E-mail: danzhou@hebtu.edu.cn

1. 西北农林科技大学理学院, 杨凌　712100

2. 河北师范大学物理学院, 石家庄　050024

收稿日期: 2026-04-09

关键词:

Charmonium potential models based on neural network wave functions

Corresponding author: E-mail: erliang.cui@nwafu.edu.cn;

Corresponding authors: E-mail: danzhou@hebtu.edu.cn

1. 西北农林科技大学理学院, 杨凌　712100

2. 河北师范大学物理学院, 石家庄　050024

Received Date: 2026-04-09

Keywords:

全文HTML

1. 引　言

粲偶素作为由正反粲夸克构成的束缚态体系, 是研究量子色动力学(quantum chromodynamics, QCD)在非微扰能区的重要系统. 实验上, 粲偶素谱的研究经历了从早期$ {\mathrm{e}}^+{\mathrm{e}}^- $对撞实验到近年来北京谱仪(BESIII)、大型强子对撞机(LHCb)等高精度实验的长期积累, 获得了大量关于粲偶素态的质量、宽度、衰变分支比等精确数据^[1-7]. 与此同时, 理论上的夸克势模型通过引入唯象的势函数, 在描述粲偶素能谱方面取得了显著成功^[8-12]. 然而, 传统的势参数拟合通常依赖于手动调整或$ \chi^2 $最小化^[13], 不同势函数形式的引入使得参数空间更为复杂, 传统拟合方法的计算流程相对烦琐, 难以快速响应实验进展的实时更新.

近年来, 机器学习(machine learning, ML)技术, 尤其是深度神经网络在求解偏微分方程和量子力学方程中展现出强大的能力^[14-20]. 其中, 神经网络波函数(neural network wavefunction, NNWF)作为一种变分拟设形式, 能够以高精度逼近量子系统的基态和激发态, 并被成功应用于强子物理的质量计算^[21-23]. 相比传统的变分方法, NNWF结合自动微分(automatic differentiation, AD)技术, 可以精确计算能量期望值关于势参数的梯度, 消除了有限差分带来的截断误差, 从而使得基于梯度的全局优化算法能够高效地在高维参数空间中搜索最优解, 这一特性为构建自动化势参数拟合框架提供了数学基础.

本文在一种基于神经网络波函数的量子求解器中^[24], 发展了一套完整势参数拟合工具, 用于粲偶素质量谱及其他重介子的唯象学研究. 该工具将机器学习与传统的联合拟合方法结合, 能够自动完成神经网络训练与势参数优化, 为强子质量谱的高精度计算提供了新思路和方法. 当实验数据更新时, 可快速重新拟合, 使理论计算始终与实验保持同步. 这一特性对于理论上解释当前高能物理不断涌现的强子谱实验结果和理论上预言新强子态尤为重要.

2. 理论与方法

2.1. 唯象势模型

在非相对论性夸克模型中, 粲偶素体系由一对正反粲夸克通过唯象势相互作用描述. 势能项$ V_{{\mathrm{c }}\bar{{\mathrm{c}}}}(r) $包含了从量子色动力学在非微扰能区唯象提取的相互作用, 其形式综合考虑了短程的单胶子交换、长程的禁闭效应以及自旋相关的精细结构^[10]. 具体表达式如下:

其中r为正反夸克之间的相对距离, $ m_{\mathrm{c}} $和$ m_{\bar{{\mathrm{c}}}} $分别为粲夸克和反粲夸克的质量, $ C_{\mathrm{f}} = -4/3 $为色因子, $ \alpha_s{\mathrm{}} $为强相互作用耦合常数, $ \sigma_{\text{lin}} $为弦张力常数, $ \sigma_{\text{smear}} $为涂抹半径参数, $ {\boldsymbol{s}}_i \cdot {\boldsymbol{s}}_j $为自旋-自旋耦合算符, $ {\boldsymbol{l}} \cdot {\boldsymbol{s}} $为自旋-轨道耦合算符, $ T_{\text{force}} $为张量相互作用算符.

精细作用项中, 自旋-轨道耦合项以及张量力项在原点附近具有$ 1/r^3 $的发散行为, 直接作为势能项纳入薛定谔方程求解时会使基态能量发散^[25]. 由于波函数在原点附近受离心势抑制, 微扰修正有限且收敛, 因此, 本文采用微扰论方法. 在仅包含$ V_{\text{cornell}} = C_{\mathrm{f}}\alpha _{\mathrm{s}}/r+\sigma _{\text{lin}}r $和$ V_{\mathrm{ss}}=\dfrac{32\text{π} \alpha _{\mathrm{s}}}{9 m_{\mathrm{c}}m_{\bar{{\mathrm{c}}}}} \dfrac{\sigma _{\text{smear}}^{3}}{\text{π} ^{3/2}}\times {\mathrm{e}}^{-\sigma _{\text{smear}}^{2}r^2} ({\boldsymbol{s}}_i\cdot {\boldsymbol{s}}_j) $的哈密顿量下求解体系的本征波函数, 然后将$ {\boldsymbol{l}} \cdot {\boldsymbol{s}} $耦合项与张量力项在该波函数下的期望值作为能级的微扰修正.

上述势模型共包含5个可调参数: 正反粲夸克质量$ m_{\mathrm{c}} $与$ m_{\bar{{\mathrm{c}}}} $、强相互作用耦合常数$ \alpha_{\mathrm{s }}$、弦张力$ \sigma_{\text{lin}} $、高斯涂抹半径$ \sigma_{\text{smear}} $. 在传统方法中, 这些参数通过拟合实验能谱得到, 通常采用手动调整或$ \chi^2 $最小化流程, 本文引入基于神经网络波函数的变分方法, 将势参数与波函数进行联合优化, 实现从实验能谱到唯象势模型自动化拟合.

2.2. 神经网络波函数优化方法

为了精确求解径向薛定谔方程并获得高精度本征波函数, 引入神经网络波函数作为变分解, 并特别引入了可显式控制节点数目的构造方法, 以精确描述激发态波函数的节点结构. 设约化径向波函数$ u(r;{\boldsymbol{w}}) $由神经网络$ {\cal{N}}(r; {\boldsymbol{w}}) $表示:

其中$ r = \|r\| $为径向距离, $ {\boldsymbol{w}} $为神经网络中的待优化参数, $ {\cal{C}}(r) $为约束函数, $ {\cal{N}}(r; {\boldsymbol{w}}) $为可训练的神经网络, 选用全连接网络的多层感知机结构, 网络输入为径向坐标r, 输出为一个标量值, $ {\cal{P}}(r) $为节点因子. 为保证波函数在无穷远处的渐近行为符合物理要求, 引入约束函数

其中α为可调衰减系数. 该高斯包络形式确保了波函数在远离原点时指数衰减, 符合束缚态波函数的边界条件.

波函数的节点(即零点)是激发态的关键特征, 直接影响能量本征值. 本文通过可训练的节点参数实现节点结构的灵活构造:

其中M为节点数量, $ \{r_i\}_{i=1}^{M} $为可训练节点参数, γ为缩放因子, 用于控制节点处波函数的过渡宽度. 可根据目标激发态的节点结构自由设置节点数目, 通常第n个激发态波函数在径向具有$ n-1 $个节点, 节点位置$ r_i $作为网络参数参与训练, 使得网络能够自动优化节点分布使能量最小化, 无需预先指定节点位置. 使用$ \tanh $函数替代简单的线性因子$ r - r_i $, 使得节点处波函数光滑过零, 避免线性节点带来的梯度突变, 提高训练稳定性.

为了增强网络表达能力和训练稳定性, 采用残差连接结构. 网络由多个隐藏层构成, 每层神经元数量可自定义. 第一层为线性变换层, 将输入映射到隐藏空间; 后续层使用残差块, 其结构为

其中σ为激活函数. 本文采用高斯激活函数, 定义为

该激活函数具有局部响应特性, 适合量子态波函数的光滑性要求.

采用能量方差最小化优化方案, 损失函数定义为

其中N为采样点数, $ E_{{\mathrm{L}}, i}({\boldsymbol{w}}) $为第i个采样点处的局部能量, 定义为$E_{{\mathrm{L}}, i}({\boldsymbol{w}})=Hu_i(r;{\boldsymbol{w}})/u_i(r;{\boldsymbol{w}})$, H为径向哈密顿算符, $ \bar{E}_{\mathrm{L }}$为平均局部能量. 通过规定节点数M及轨道角动量l, 本方法可自动计算对应激发态.

2.3. 势模型参数拟合

得到本征波函数及本征能谱后, 采用基于高斯似然的最大似然估计进行势模型参数拟合. 为量化参数的不确定度, 进一步引入拉普拉斯近似: 在最大似然估计点处对负对数似然进行二阶泰勒展开, 将似然函数近似为高斯分布, 从而获得参数的最优估计与协方差矩阵.

由于实验误差跨度极大, 直接采用$ \chi^2 $拟合会导致高精度数据在拟合中占主导地位, 低精度数据贡献被严重压制. 考虑到势模型的近似性以及计算本身的不确定性, 选取5 MeV作为固定理论误差, 与实验误差的平方和作为总体误差项, 避免高精度数据在拟合中过度主导, 为势模型下把一系列粲偶素态作为整体进行联合分析提供了计算依据.

在高斯似然假设下, 定义负对数似然函数为

其中$ {\boldsymbol{\theta}}=(m_{\mathrm{c}}, m_{\bar{{\mathrm{c}}}}, \alpha_{\mathrm{s}}, \sigma_{\text{lin}}, \sigma_{\text{smear}}) $为待优化参数, N为实验粒子总数, $ E_i^{\text{th}}({\boldsymbol{\theta}}) $为质量计算值, $ \sigma_i^{\text{exp}} $为实验误差, $ \sigma_{\text{th}} $为固定理论误差. 通过最小化$ {\cal{L}}_{\text{NLL}}({\boldsymbol{\theta}}) $得到最大似然估计点$ {\boldsymbol{\theta}}_{\text{MLE}} $作为势模型的最优参数组合.

在$ {\boldsymbol{\theta}}_{\text{MLE}} $点处计算黑塞矩阵(Hessian matrix)$ {\boldsymbol{H}} = \nabla_{{\boldsymbol{\theta}}}^2 {\cal{L}}_{\text{NLL}}({\boldsymbol{\theta}})|_{{\boldsymbol{\theta}}_{\text{MLE}}} $, 则参数协方差矩阵$ {\boldsymbol{\varSigma}} = {\boldsymbol{H}}^{-1} $, 参数不确定度由$ \sqrt{\text{diag}({\boldsymbol{\varSigma}})} $给出. 该方法利用自动微分完成优化与不确定度评估, 输出势模型参数的最优值及其协方差.

3. 结果与讨论

3.1. 训练设置

所有计算均基于统一的计算平台: CPU采用Intel Core i7-13700K, GPU使用NVIDIA GeForce RTX 4060 Ti. 神经网络模型与训练流程完全基于PyTorch框架实现, 并利用CUDA加速张量运算. 网络输入为标量相对距离r, 随后经过6个隐藏层、每层128个神经元, 其中第1层为标准线性层, 后5层为残差块, 最终输出1个波函数标量值.

约束函数((3)式)中的衰减系数取固定值$ \alpha=0.008 $, 节点因子((4)式)中的缩放因子同样取固定值$ \gamma=1.0 $; 高斯激活函数((6)式)中的参数β为可训练参数, 初始值设为1.0. 蒙特卡罗采样点数$ N=1000 $, 训练采用AdamW优化器, 初始学习率为$ 10^{-3} $, 权重衰减系数为$ 10^{-4} $. 学习率调度器采用自定义的$ \mathrm{LambdaLR}$, 前2000步保持不变, 之后每2000步衰减为原来的$ 0.1 $倍. 为兼顾计算精度与计算成本, 在100步预热后, 当能量方差降至$ 3\times10^{-5} $时停止训练, 典型训练步数为1000—3000步. 为避免由特定随机种子选择引入人为偏好, 本文主计算未预先固定随机种子; 因此, 神经网络参数初始化和蒙特卡罗采样过程均保持其自然随机性.

3.2. 势参数优化

势模型的原始参数组合来自文献[13], 参数分别取为$ m_{\mathrm{c}} = m_{\bar{{\mathrm{c}}}} = 1.4794\ \text{GeV} $, $ \alpha_{\mathrm{s}} = 0.5461 $, $ \sigma_{\text{lin}} = 0.1425\ \text{GeV}^2 $, $ \sigma_{\text{smear}} = 1.09467\ \text{GeV} $. 优化后的参数组合为$ m_{\mathrm{c}} = m_{\bar{{\mathrm{c}}}} = 1.4800811\ \text{GeV} $, $ \alpha_{\mathrm{s }}= 0.5458532 $,$ \sigma_{\text{lin}} = 0.1441820\ \text{GeV}^2 $, $ \sigma_{\text{smear}} = 1.1046847\ \text{GeV} $. 不确定度分别为$ \Delta m_{\mathrm{c }}= \Delta m_{\bar{{\mathrm{c}}}} = 0.0054706\ \text{GeV} $, $ \Delta \alpha_{\mathrm{s}} = 0.0120337 $, $ \Delta \sigma_{\text{lin}} = 0.0028496\ \text{GeV}^2 $, $ \Delta \sigma_{\text{smear}} = 0.0567969\ \text{GeV} $. 其中$ m_{\mathrm{c}} $以及$ \sigma_{\text{lin}} $的不确定度较低, 偏差较小, 而$ \sigma_{\text{smear}} $的不确定度较高, 约$ 56.8\ \text{MeV} $.

选取Particle Data Group发布的最新粲偶素实验值^[1] (截止至2025年1月), 使用原始参数组合以及优化后新参数组合分别计算出能谱并给出相对偏差$ (E^{\text{th}} - E^{\text{exp}})/E^{\text{exp}} \times 100{\text{%}}$, 其中$ E^{\text{th}} $为计算值, $ E^{\text{exp}} $为实验值.

表1对比了新参数组合与原始参数组合下的能谱相对偏差, 优化后多数已知态的相对偏差有所降低. 为了量化势参数优化对能谱描述的改进, 计算了16个已知粲偶素态的计算值与实验值的均方根偏差(root mean square deviation, RMSD), 其定义为

其中$ E_i^{\text{th}} $为计算值, $ E_i^{\text{exp}} $为实验值. 优化前RMSD为$ 26.26\ \text{MeV} $, 优化后降至$ 24.62\ \text{MeV} $, 这表明模型参数经过优化后, 整体偏差水平些许下降. 此外, 部分高激发态的相对偏差有所增加, 分别为$ 1{\,}^3 {\mathrm{P}}_2 $, $ 1{\,}^3 {\mathrm{D}}_1 $, $ 2{\,}^3 {\mathrm{P}}_2 $以及$ 3{\,}^3 {\mathrm{S}}_1 $, 该偏差的增加可能是由于提升低精度数据的贡献而引入的. 另外, 由于理论中通常将实验值$ \psi(4040) $归类为$ 3{\,}^3 {\mathrm{S}}_1 $, 然而实验中明确指出$ 4.0\ \text{GeV} $附近存在显著的阈值效应(例如$ {\mathrm{D}}\bar{{\mathrm{D}}} $, $ {\mathrm{D}}^*\bar{{\mathrm{D}}} $等开粲道的阈值), 导致该能区可能存在多个共振峰的混合或复杂干涉, 因此, 该态在拟合中可能不适宜作为纯$ 3{\,}^3 {\mathrm{S}}_1 $处理, 未来可考虑引入耦合道效应^[1,26-28].

考虑到参数不确定性以及计算精度, 本文取小数点后4位作为最终结果: $ m_{\mathrm{c}} = m_{\bar{{\mathrm{c}}}} = (1.4801 \pm 0.0055)\ \text{GeV} $, $ \alpha_{\mathrm{s}} = 0.5459 \pm 0.0120 $, $ \sigma_{\text{lin}} = (0.1442 \pm 0.0028)\ \text{GeV}^2 $, $ \sigma_{\text{smear}} = (1.1047 \pm 0.0568)\ \text{GeV} $.

3.3. 粲偶素质量谱预测

基于上述优化后的势参数, 进一步计算了更多可能的粲偶素质量谱, 涵盖S, P, D, F波等多个轨道角动量激发态, 结果列于表2, 通过神经网络波函数方法获得的各轨道角动量通道的径向概率密度分布$ \left|u\left(r\right)\right|^2 $见图1. 为便于在统一谱表中比较已知态与待观测态的位置关系, 表2保留了实验已确认的粒子态并以粗体标出, 同时列出了优化参数下预测的其余高激发态和高角动量态.

从能级结构来看, $ n{\,}^3 {\mathrm{S}}_1 $系列($ J/\psi $, $ \psi(2{\mathrm{ S}}) $, $ \psi(3 {\mathrm{S}}) $, $ \psi(4 {\mathrm{S}}) $)的相邻能级间隔分别为$ 585.16$,$ 402.94 $和$ 336.63\; \text{MeV} $, 随主量子数增大逐渐减小, 符合势模型中夸克禁闭势的预期行为. $ n{\,}^1 {\mathrm{S}}_0 $系列($ \eta_{\mathrm{c}} $粒子)亦呈现相同趋势. 对于P波态, 精细结构能级顺序为$ \chi_{{\mathrm{c}}0} \lt \chi_{{\mathrm{c}}1} \lt h_{\mathrm{c}} \lt \chi_{{\mathrm{c}}2} $, 与实验观测的能级顺序一致, 表明模型对自旋相关相互作用的描述是合理的.

在粲偶素态预言方面, 本模型预测了$ \chi_{{\mathrm{c}}J}(n{\,}^3 {\mathrm{F}}_J) $和$ \eta_{{\mathrm{c}}3}(n{\,}^1 {\mathrm{F}}_3) $系列的高角动量态($ l=3 $), 其质量分布在$ 4030 $—$ 4930\ \text{MeV} $区间, 与其他模型的理论预言基本一致^[29-33]. 这些态尚未被实验确认, 未来在BESIII, Belle II, LHCb等高亮度对撞实验中, 通过合适的产生过程(如双光子产生、辐射跃迁等)有望被观测到.

[1]	Navas S, Amsler C, Gutsche T, et al 2024 Phys. Rev. D 110 030001 doi: 10.1103/PhysRevD.110.030001
[2]	Aubert J J, Becker U, Biggs P J, et al 1974 Phys. Rev. Lett. 33 1404 doi: 10.1103/PhysRevLett.33.1404
[3]	Augustin J E, Boyarski A M, Breidenbach M, et al 1974 Phys. Rev. Lett. 33 1406 doi: 10.1103/PhysRevLett.33.1406
[4]	Appelquist T, Politzer H D 1975 Phys. Rev. Lett. 34 43 doi: 10.1103/PhysRevLett.34.43
[5]	Bai J Z, et al 1998 Phys. Rev. Lett. 81 3091 doi: 10.1103/PhysRevLett.81.3091
[6]	Yuan C Z 2021 Nat. Sci. Rev. 8 nwab182 doi: 10.1093/nsr/nwab182
[7]	Wang X, Liu X, Gao Y 2026 Rev. Mod. Phys. 98 021001 doi: 10.1103/2mrp-chly
[8]	Eichten E, Gottfried K, Kinoshita T, Lane K D, Yan T M 1978 Phys. Rev. D 17 3090 doi: 10.1103/PhysRevD.17.3090
[9]	Bykov A A, Dremin I M, Leonidov A V 1984 Phys.-Usp. 27 321 doi: 10.1070/pu1984v027n05abeh004291
[10]	Godfrey S, Isgur N 1985 Phys. Rev. D 32 189 doi: 10.1103/PhysRevD.32.189
[11]	Lucha W, Schöberl F F, Gromes D 1991 Phys. Rep. 200 127 doi: 10.1016/0370-1573(91)90001-3
[12]	Bokade C A, Bhaghyesh 2025 Phys. Rev. D 111 014030 doi: 10.1103/PhysRevD.111.014030
[13]	Barnes T, Godfrey S, Swanson E S 2005 Phys. Rev. D 72 054026 doi: 10.1103/PhysRevD.72.054026
[14]	Carleo G, Troyer M 2017 Science 355 602 doi: 10.1126/science.aag2302
[15]	Pfau D, Spencer J S, Matthews A G D G, Foulkes W M C 2020 Phys. Rev. Res. 2 033429 doi: 10.1103/PhysRevResearch.2.033429
[16]	Raissi M, Perdikaris P, Karniadakis G E 2019 J. Comput. Phys. 378 686 doi: 10.1016/j.jcp.2018.10.045
[17]	Han J, Zhang L, E W 2019 J. Comput. Phys. 399 108929 doi: 10.1016/j.jcp.2019.108929
[18]	Vicentini F, Hofmann D, Szabó A, Wu D, Roth C, Giuliani C, Pescia G, Nys J, Vargas-Calderon V, Astrakhantsev N, Carleo G 2022 SciPost Phys. Codebases 7 doi: 10.21468/SciPostPhysCodeb.7
[19]	Manzhos S 2020 Mach. Learn.-Sci. Techn. 1 013002 doi: 10.1088/2632-2153/ab7d30
[20]	Han J, Lu J, Zhou M 2020 J. Comput. Phys. 423 109792 doi: 10.1016/j.jcp.2020.109792
[21]	Mutuk H 2019 Adv. High Energy Phys. 2019 3105373 doi: 10.1155/2019/3105373
[22]	Mutuk H 2019 Chin. Phys. C 43 093103 doi: 10.1088/1674-1137/43/9/093103
[23]	Wu W L, Meng L, Zhu S L 2026 Phys. Rev. Lett. 136 071901 doi: 10.1103/ckpr-s876
[24]	Shi H C, Cui E L, Zhou D 2025 arXiv: 2510.17490 v2 [quant-ph]
[25]	Ding Y B, Li X Q, Shen P N 2000 Commun. Theor. Phys. 33 613 doi: 10.1088/0253-6102/33/4/613
[26]	Ablikim M, et al 2008 Phys. Lett. B 660 315 doi: 10.1016/j.physletb.2007.11.100
[27]	Seth K K 2005 Phys. Rev. D 72 017501 doi: 10.1103/PhysRevD.72.017501
[28]	Mo X H, Yuan C Z, Wang P 2010 Phys. Rev. D 82 077501 doi: 10.1103/PhysRevD.82.077501
[29]	Radford S F, Repko W W 2007 Phys. Rev. D 75 074031 doi: 10.1103/PhysRevD.75.074031
[30]	Mahmood T, Darwish J A, Hussain T, Ahmed M, Sherwani R A K 2024 Adv. High Energy Phys. 2024 5195790 doi: 10.1155/2024/5195790
[31]	Cao L, Yang Y C, Chen H 2012 Few-Body Systems 53 327 doi: 10.1007/s00601-012-0478-z
[32]	Liao X, Manke T 2002 arXiv: hep-lat/0210030 [hep-lat]
[33]	Liu L, Moir G, Peardon M, Ryan S M, Thomas C E, Vilaseca P, Dudek J J, Edwards R G, Joó B, Richards D G 2012 J. High Energy Phys. 2012 126 doi: 10.1007/JHEP07(2012)126