Extending pairing energy density functional using pairing rotational moments of inertia

TL;DR

This paper enhances pairing energy density functionals by incorporating derivative terms, improving the modeling of pairing rotational moments of inertia in isotopes, with adjustments based on experimental data.

Contribution

It introduces pair-density derivative terms into the pairing EDF and demonstrates their systematic improvement of rotational moments of inertia calculations.

Findings

Pair-density derivative terms improve the accuracy of rotational moments of inertia.

Adjusted coupling constants reproduce experimental data.

Conventional pairing gaps overestimate odd-even mass staggering.

Abstract

The pairing energy density functionals (EDFs) that include the spatial derivative and kinetic terms of the pair densities are discussed. The coupling constants of the pairing EDF are adjusted to reproduce the experimental pairing rotational moment of inertia, and the pair-density derivative terms are shown to systematically improve the values of the pairing rotational moments of inertia in Sn and Pb isotopes. It is pointed out that the conventional average pairing gaps overestimate the experimental odd-even mass staggering in the presence of the pair-density derivative terms.

Tables

Table 1. Table 1: Pairing coupling constants V n subscript 𝑉 𝑛 V_{n} fitted to the neutron pairing rotational moment of inertia in 120 Sn for each value of η 𝜂 \eta and C ~ n Δ ​ ρ superscript subscript ~ 𝐶 𝑛 Δ 𝜌 \tilde{C}_{n}^{\Delta\rho} .

${\tilde{C}}_n^{\Delta \rho }$ (MeV fm5)	${\tilde{C}}_n^{\tau }$ (MeV fm5)	$V_n$ (MeV fm3)
		$\eta =0$	$\eta =0.5$	$\eta =1$
0	0	$-$148.6	$-$225	$-$389
$-$20	80	$-$593	$-$820	$-$1093
$-$40	160	$-$774	$-$1072	$-$1450
$-$60	240	$-$720	$-$993	$-$1499
$-$80	320	$-$354	$-$521	$-$924

Equations

$${\Delta }_n^{(3)}(N)=\frac{(-1{)}^N}{2}[E(N+1)-2E(N)+E(N-1)],$$

$${\tilde{\Delta }}_n^{(3)}(N)=\frac{1}{2}[{\Delta }_n^{(3)}(N-1)+{\Delta }_n^{(3)}(N+1)].$$

$${\Delta }_n(\rho )\equiv \frac{\mathrm{Tr}{\tilde{h}}_n{\rho }_n}{\mathrm{Tr}{\rho }_n}=-\frac{2}{N}\int \mathrm{d}\mathbf{r}{\tilde{C}}_n^{\rho }[{\rho }_0]{\rho }_n(\mathbf{r}){\tilde{\rho }}_n(\mathbf{r}),$$

$${\Delta }_n(\tilde{\rho })\equiv \frac{\mathrm{Tr}{\tilde{h}}_n{\tilde{\rho }}_n}{\mathrm{Tr}{\tilde{\rho }}_n}=-\frac{2}{\tilde{N}}\int \mathrm{d}\mathbf{r}{\tilde{C}}_n^{\rho }[{\rho }_0]{\tilde{\rho }}_n^2(\mathbf{r}),$$

$$E(N)=E(N_0)+{\lambda }_n(N_0)(N-N_0)+\frac{1}{2{\mathcal{J}}_n(N_0)}(N-N_0{)}^2,$$

$${\lambda }_n(N)=-\frac{1}{4}[S_{2n}(N+2)+S_{2n}(N)],$$

$${\mathcal{J}}_n(N)=\frac{4}{E(N+2)-2E(N)+E(N-2)}=\frac{4}{S_{2n}(N)-S_{2n}(N+2)}.$$

$$\left(\begin{array}{cc}A& B\\ B^{\ast }& A^{\ast }\end{array}\right)\left(\begin{array}{c}{P}_i\\ -P_i^{\ast }\end{array}\right)=$$

$$\left(\begin{array}{cc}A& B\\ B^{\ast }& A^{\ast }\end{array}\right)\left(\begin{array}{c}{Q}_i\\ -Q_i^{\ast }\end{array}\right)=$$

$${\mathcal{J}}_n=2N_n^{20}(A+B{)}^{-1}{N}_n^{20},$$

$$\mathrm{i}\frac{\partial }{\partial t}{a}_{\mu }(t)=[\hat{H}(t)+\hat{F}(t),a_{\mu }(t)]$$

$$\hat{F}(t)=$$

$$\hat{F}=$$

$$\delta a_{\mu }(t)=$$

$$\delta \hat{H}(t)=$$

$$\delta H(\omega )=$$

$$(E_{\mu }+E_{\nu }-\omega )X_{\mu \nu }(\omega )+\delta H_{\mu \nu }^{20}(\omega )=-F_{\mu \nu }^{20},$$

$$(E_{\mu }+E_{\nu }+\omega )Y_{\mu \nu }(\omega )+\delta H_{\mu \nu }^{02}(\omega )=-F_{\mu \nu }^{02},$$

$$\left(\begin{array}{c}X(\omega )\\ Y(\omega )\end{array}\right)=-{\left[\left(\begin{array}{cc}A& B\\ B^{\ast }& A^{\ast }\end{array}\right)-\omega \left(\begin{array}{cc}1& 0\\ 0& -1\end{array}\right)\right]}^{-1}\left(\begin{array}{c}{F}^{20}\\ F^{02}\end{array}\right).$$

$$\left(\begin{array}{cc}\delta H^{11}& \delta H^{20}\\ -\delta H^{02}& -(\delta H^{11}{)}^T\end{array}\right)=\left(\begin{array}{cc}{U}^{†}& V^{†}\\ V^T& U^T\end{array}\right)\left(\begin{array}{cc}\delta h& \delta {\tilde{h}}^{(+)}\\ -\delta {\tilde{h}}^{(-)\ast }& -\delta h^T\end{array}\right)\left(\begin{array}{cc}U& V^{\ast }\\ V& U^{\ast }\end{array}\right),$$

$$\delta h=$$

$$\delta {\tilde{h}}^{(+)}=$$

$$\delta {\tilde{h}}^{(-)}=$$

$$U_a=$$

$$U_b=$$

$$S(\hat{F},\omega )=$$

$$=$$

$$+\mathrm{i}\omega \left(⟨Q_i|\hat{F}|0{⟩}^{\ast }⟨P_i|\hat{F}|0⟩-⟨P_i|\hat{F}|0{⟩}^{\ast }⟨Q_i|\hat{F}|0⟩\right)\},$$

$$S({\hat{N}}_n,\omega =0)=$$

$$=$$

$${\tilde{\chi }}_t(\mathbf{r})=\frac{V_t}{4}\left[1-\eta {\left(\frac{{\rho }_0(\mathbf{r})}{{\rho }_c}\right)}^{\beta }\right]|{\tilde{\rho }}_t(\mathbf{r}){|}^2,$$

$${\tilde{\chi }}_t(\mathbf{r})=$$

$$+{\tilde{C}}_t^{J0}|{\tilde{J}}_t{|}^2+{\tilde{C}}_t^{J1}|{\tilde{\mathbf{J}}}_t{|}^2+{\tilde{C}}_t^{J2}|{\widetilde{\underline{\mathsf{J}}}}_t{|}^2+{\tilde{C}}_t^{\nabla J}\mathrm{Re}({\tilde{\rho }}_t^{\ast }\nabla \cdot {\tilde{\mathbf{J}}}_t),$$

$${\tilde{\rho }}_t(\mathbf{r})=$$

$${\tilde{\tau }}_t(\mathbf{r})=$$

$${\tilde{\mathsf{J}}}_{tab}(\mathbf{r})=$$

$${\tilde{J}}_t(\mathbf{r})=$$

$${\tilde{\mathbf{J}}}_{ta}(\mathbf{r})=$$

Full text

Keywords: spontaneous symmetry breaking, pairing rotation, QRPA, energy density functional

1 Introduction

Pairing is one of the important many-body correlations in finite nuclei and other quantum systems [2, 3, 4]. The existence of the pairing-type correlation in finite nuclei has been known by a number of experimental data. The odd-even staggering of the binding energy (OES) is the most commonly used observable relevant to the nuclear pairing. In systems with even number of neutrons (or protons), all the neutrons can form pairs and thus obtain maximal pairing correlations, while in systems with odd number of neutrons, one last neutron cannot form a pair, and the systems are less bound than the neighboring systems with even number of neutrons.

Pairing correlations can be taken into account within the mean-field approximations such as Bardeen-Cooper-Schrieffer (BCS) and Hartree-Fock-Bogoliubov (HFB) theories [5]. Quasiparticles which are linear combination of the particle creation and annihilation allow us to describe a superconducting ground state as a one-body field. For the qualitative and quantitative description of the pairing properties of the atomic nuclei, it is essential to identify appropriate observable that is sensitive to the detail of the nuclear pairing. Conventionally the theoretical pairing gap is compared with the experimental OES, since the direct computation of the odd-mass systems in the mean-field approximation is not as accurate as that of the even-mass systems. As the relation between the pairing gap and OES is indirect, it is not easy to discuss the detailed information on the pairing correlation from the pairing gap without ambiguity.

In the mean-field approximation based on a nuclear effective interaction, the particle-hole and pairing potentials are derived from a single nuclear effective Hamiltonian. In actual applications based on the zero-range Skyrme effective interactions (with a few exceptions such as SkP [6] and application to the nuclear matter [7, 8]), however, Skyrme effective interactions are used only for the particle-hole channel, and simplified forms are assumed for the pairing interaction. It is because that the pairing properties derived from the Skyrme interactions are found to be unrealistic [9]. Such a current standard prescription for the mean-field theory with the Skyrme interactions is closely connected to the nuclear density functional theory (DFT). In nuclear DFT, the correspondence between the energy density functional (EDF) and the effective interaction is not required, and the particle-hole and pairing parts of the EDF can be independent.

The coupling constants in the nuclear EDF is phenomenologically optimized using a selected experimental data set. Recently new EDF parametrizations such as SVmin [10] and UNEDF [11, 12, 13, 14], based on the optimizations using a large number of experimental data set have been proposed. Although many coupling constants are considered in the particle-hole channel, the pairing functionals used there are still limited to a simple form. This is because the optimization of the pairing part relies mainly on the theoretical pairing gap and experimental OES. As commented in [11], a proper treatment of odd-mass systems including the evaluation of the time-odd fields is unavoidable in the optimization of the pairing functional in a more general form using the OES.

The pairing correlation introduces spontaneous breaking of the gauge symmetry. A zero-energy Nambu-Goldstone (NG) mode [15, 16, 17], called pairing rotation appears associated with the broken gauge symmetry [18]. In our previous work [19], we have proposed that the moments of inertia for the pairing rotation that are related to the two-nucleon shell gap indicators $\delta_{2n}$ and $\delta_{2p}$ [20] and the proton-neutron interaction energy $\delta V_{pn}$ [21], are excellent pairing indicators. By employing the pairing rotational moments of inertia as pairing observables, the issues related to the evaluation of the odd-mass systems can be avoided, and the details of the pairing part of the nuclear EDF is expected to be constrained.

The aim of this article is to show the usefulness of the pairing rotational moments of inertia as nuclear pairing observables, and to demonstrate its sensitivity to the form of the pairing EDF. We extend the pairing EDF by including the spacial derivative terms of the pair density. These terms are originated from the momentum-dependent terms in the Skyrme effective interaction. In section 2, we describe the problem of using the conventional pairing gaps as approximations to the experimental OES. Then we introduce the pairing rotational moment of inertia as a pairing observable. In section 3, an efficient microscopic method for computing the pairing rotational moment of inertia within the linear response theory for nuclear DFT is recapitulated [22, 19]. The extended form of the pairing EDF including the pair-density derivative terms is presented in section 4, and the results of the numerical calculations for neutron pairing in Sn and Pb isotopes are discussed in section 5. Conclusions are given in section 6.

2 Pairing observable

2.1 Pairing gap and OES

Conventionally the OES is used for the observable related to the nuclear pairing property [9]. The simplest three-point formula for the OES is defined by the ground-state energy differences of the three neighboring isotopes [23, 24]

[TABLE]

where $E(N)$ is the ground-state energy of an isotope with $N$ neutrons. In order to reduce the local fluctuation in the OES, the average of the two OES values of even-odd and odd-even nuclei are used [13]

[TABLE]

In this paper we use definition (2) for the experimental OES.

In the mean-field calculation, instead of evaluating the OES by either equation  (1) or (2), the theoretical pairing gap is usually evaluated and compared with the OES. For the conventional pairing EDF $\tilde{\chi}_{n}[\rho,\tilde{\rho},\tilde{\rho}^{\ast}]=\tilde{C}_{n}^{\rho}[\rho_{0}]|\tilde{\rho}_{n}(\bm{r})|^{2}$ with an isoscalar-density-dependent coupling constant $\tilde{C}_{n}^{\rho}[\rho_{0}]$, the pairing gaps are defined by averaging the pairing potential $\tilde{h}_{n}(\bm{r})=\partial\tilde{\chi}_{n}/\partial\tilde{\rho}_{n}^{\ast}$ with the particle-hole or pair density. One definition of the neutron pairing gap is

[TABLE]

where $N$ is the neutron number. This is the pairing potential averaged with the particle-hole density, and is often used as the theoretical pairing gap [6, 25]. Another definition of the pairing gap is the pairing potential averaged with the pair density and thus the focus is put more on the properties around the Fermi surface [26, 27]

[TABLE]

where $\tilde{N}=\int{\rm d}\bm{r}\tilde{\rho}_{n}(\bm{r})$.

In figure 1, we show typical examples of the neutron pairing gap and OES for Sn and Pb isotopes. The experimental data are evaluated with equation (2) using the values in [28], while the average pairing gap (3) is computed using the UNEDF1-HFB EDF [14] with three kinds of density-dependent pairing functionals. The pairing strengths are fitted to reproduce the neutron OES at 120Sn. The general trends agree quite well, but the pairing gap does not have sensitivity to determine the density dependence of the pairing functional through the OES.

It has been pointed out many times that the pairing gap is not an appropriate quantity to estimate the OES. Existence of the multiple definitions of the pairing gaps such as equations (3) and (4), and also in the experimental OES (such as four-point, five-point formulae) makes the unique comparison impossible. In order to avoid the above issue, one has to evaluate the OES directly by computing the odd-mass systems. As the time-reversal symmetry is broken in the odd-mass systems, the precise computation of the odd-mass systems is still a challenging problem in the nuclear DFT [29]. In addition to the blocking prescription, the time-odd functionals are activated, and their coupling constants are less known than those of the time-even functionals [30].

2.2 Pairing rotational moment of inertia

Another important aspect of the pairing correlation is the spontaneous breaking of the gauge symmetry. The mean-field solution with the pairing correlations becomes superconducting, and is no longer the eigenstate of the particle number operator. There appears a new collective degree of freedom in the symmetry-broken state, called zero-energy NG mode [15, 16, 17]. The pairing NG mode, generated from the broken U(1) gauge symmetry is called pairing rotation that corresponds to the rotation of a complex phase [3, 18]. The pairing rotation is a zero-energy excitation in the intrinsic superconducting frame, while in the laboratory frame where the particle number is conserved, the ground states of the even-even nuclei can be regarded as a pairing rotational band.

For a system in which the proton shell is closed and the neutron shell is open and superconducting, the ground-state energy for a system with $N$ neutrons can be described from a reference system with $N_{0}$ neutrons by expanding the energy up to the second order with respect to the particle number fluctuation $N-N_{0}$,

[TABLE]

where $\lambda_{n}(N_{0})={\rm d}E/{\rm d}N|_{N=N_{0}}$ is the neutron chemical potential, and ${\cal J}_{n}(N_{0})={\rm d}^{2}E/{\rm d}N^{2}|_{N=N_{0}}^{-1}$ is the neutron pairing rotational moment of inertia. The second order term in $N-N_{0}$ is the pairing rotational energy.

Experimentally, the neutron chemical potential $\lambda_{n}(N)$ can be extracted from the two-neutron separation energy

[TABLE]

while the neutron pairing rotational moment of inertia ${\cal J}_{n}(N)$ can be evaluated by the inverse of double binding-energy differences.

[TABLE]

We have proposed that the pairing rotational moment of inertia is an excellent pairing observable [19]. The moment of inertia of the pairing rotation is a quantity related to the gauge symmetry breaking due to the pairing correlations. The inertia is evaluated using the ground-state energies of even-$N$ systems only, and thus there is no ambiguity from the unknown time-odd functionals or the approximations when computing the odd-mass systems. We note the neutron pairing rotational moment of inertia is proportional to the inverse of the two-neutron shell gap indicator $\delta_{2n}$.

3 Finite-amplitude method (FAM) for NG modes

The pairing rotational moment of inertia is originally defined as a second derivative of the energy with respect to the particle number. Theoretically it can be evaluated from three HFB calculations using equation (7), but this expression is regarded as an approximation to the second derivative. The quasiparticle random-phase approximation (QRPA) allows us to derive the pairing rotational moment of inertia as an exact second derivative (Thouless-Valatin moment of inertia [31]). The zero-energy NG modes are the solutions of the QRPA equations. The QRPA equations in the $PQ$ (momentum-coordinate) representation are written as [5]

[TABLE]

where $A$ and $B$ are QRPA matrices that include the residual interactions, $Q_{i}$ and $P_{i}$ are the two-quasiparticle matrix elements of the coordinate and momentum operators $\hat{Q}_{i}$ and $\hat{P}_{i}$ (normalized with $\langle[\hat{Q}_{i},\hat{P}_{j}]\rangle={\rm i}\delta_{ij}$) for the eigen mode with the frequency $\Omega_{i}$ and the inertia ${\cal J}_{i}$. For the zero-energy neutron pairing rotational mode in proton-shell closed nuclei, $\Omega_{i}=0$, the momentum and coordinate operators become the neutron particle number operator $\hat{N}_{n}$ and its canonically conjugate neutron gauge angle operator $\hat{\Theta}_{n}$, respectively. The Thouless-Valatin moment of inertia for the neutron pairing rotational mode is given by

[TABLE]

where $N^{20}_{n}$ is the two-quasiparticle matrix elements of the $\hat{N}_{n}$ operator (assumed to be real).

Because the dimensions of the two-quasiparticles (the dimensions of the $A$ and $B$ matrices) are very large, and additional model space truncation is necessary for calculating the Thouless-Valatin inertia directly from equation (10). Instead, we evaluate the Thouless-Valatin inertia in an equivalent alternative approach, the FAM for the linear response theory of nuclear DFT [32, 33].

The FAM is derived from the small-amplitude limit of the time-dependent HFB (TDHFB) equation, thus it is formally equivalent to the QRPA. The TDHFB equation is written as

[TABLE]

with a weak external field that depends on time

[TABLE]

This introduces the oscillations with the same frequencies for the quasiparticle wave functions, and the mean-field Hamiltonian. The oscillating part of the quasiparticle $\delta a_{\mu}(t)=a_{\mu}(t)e^{-{\rm i}E_{\mu}t}-a_{\mu}$ and the Hamiltonian $\delta\hat{H}(t)=\hat{H}(t)-\hat{H}_{\rm HFB}$ are expressed as

[TABLE]

From the first order terms in the small parameter $\eta$ in the TDHFB equation (11), we have the linear response equations

[TABLE]

or equivalently

[TABLE]

The FAM enables us to evaluate the fluctuations of the Hamiltonian, $\delta H^{20}$ and $\delta H^{02}$, very efficiently. These quantities are written in terms of $A$ and $B$ matrices and QRPA wave functions $X$ and $Y$ (see equation (9) in [34]). In the FAM, they are evaluated as follows. The fluctuations, $\delta H^{20}$ and $\delta H^{02}$, together with other one-quasiparticle-one-quasihole part $\delta H^{11}$ are connected to the particle-hole and particle-particle fluctuations by a Bogoliubov transformation

[TABLE]

where $U$ and $V$ are the HFB wave functions. Then the fluctuations of the particle-hole and pairing parts of the Hamiltonian matrix are written in terms of the wave functions

[TABLE]

with

[TABLE]

where $h[U^{\ast},V^{\ast},U,V]$ and $\tilde{h}[U^{\ast},V^{\ast},U,V]$ are the particle-hole and pairing part of the Hamiltonian as functions of wave functions. From equations (21)-(23) the fluctuations of the particle-hole and pairing parts can be evaluated essentially using the routine that computes $h$ and $\tilde{h}$ from the HFB wave function, and the functional derivative is replaced by a numerical derivative with a finite-amplitude small parameter $\eta$. This is the reason the method is called finite-amplitude method. Using the technique above, we can solve equations (17) and (18) iteratively.

The output of the FAM is the strength functions. The expression of the strength function in the $PQ$ representation is given by [35, 22]

[TABLE]

where $\langle P_{i}|\hat{F}|0\rangle=\langle[\hat{P}_{i},\hat{F}]\rangle$ and $\langle Q_{i}|\hat{F}|0\rangle=\langle[\hat{Q}_{i},\hat{F}]\rangle$. In the standard calculations using the FAM, the imaginary part of the strength function $-{\rm Im}S(\hat{F},\omega_{\gamma})/\pi$ for a complex frequency $\omega_{\gamma}=\omega+{\rm i}\gamma$ provides the Lorentzian-smeared strength distribution.

The expression above is not valid for $\omega=0$ when the NG mode is present ($\Omega_{\rm NG}=0$). The strength function for the momentum operator of the NG mode is computed directly from equation (19)

[TABLE]

Therefore in the FAM, the Thouless-Valatin moment of inertia is given from the strength function for the zero-energy response with the broken-symmetry operator. This expression has been checked numerically for the center-of-mass motion and the pairing rotations [22].

We have compared the neutron pairing rotational moment of inertia evaluated from the Thouless-Valatin inertia using equation (27) and from the inverse of the double binding-energy differences (7) using three HFB states in figure 2 in order to see the consistency between these two evaluations. The calculations are performed for Sn and Pb isotopes, and the pairing functional of the volume type is employed with the parameters described in the next section. The agreement between the double binding energy differences and the Thouless-Valatin inertia is generally good, but in some nuclei there is a discrepancy. We also evaluate the pairing rotational moments of inertia using the same HFB states but with the relation between the particle number and chemical potential ${\cal J}_{n}(N)=[{\rm d}\lambda_{n}(N)/{\rm d}N]^{-1}\approx(2\Delta N)/[\lambda_{n}(N+\Delta N)-\lambda_{n}(N-\Delta N)]$. The pairing rotational moment of inertia derived from the chemical potential agrees with the Thouless-Valatin inertia much better than that from the double binding-energy differences of the three HFB states does, because in finite difference the numerical error is generally larger in evaluating the second derivative than in evaluating the first derivative. The total energy is affected more by the fluctuation around the mean field than the chemical potential is. However next to the closed shell, the evaluation from the chemical potential deviates completely because the chemical potential at the closed shell is not uniquely defined (in this calculation the energy of the last occupied orbit is used as the chemical potential when the pairing collapses).

4 Pairing EDF

Skyrme EDF is a zero-range energy density functional, widely used for the nuclear EDF calculations. It is originally based on the Skyrme effective interaction [36, 37]. Therefore the coupling constants of the Skyrme EDF is derived from the parameters of the Skyrme effective interactions. In the concept of the generalized EDF, however, the EDF can be parametrized independent of the effective interaction. The bilinear density form together with the density-dependent term is assumed in the particle-hole Skyrme EDF, and the coupling constants are optimized using a large number of experimental data sets. Recent effort in this direction produces a series of the Skyrme-EDF parameter sets: SVmin [10], UNEDF0 [11], UNEDF1 [12], UNEDF1-HFB [14], and UNEDF2 [13].

In the Skyrme EDF, a simple zero-range pairing force with density dependence is commonly used [9]

[TABLE]

where $\tilde{\rho}_{t}(\bm{r})$ is a pair density of neutron or proton, $\rho_{0}(\bm{r})=\rho_{n}(\bm{r})+\rho_{p}(\bm{r})$ is the isoscalar particle-hole density, $\rho_{c}=0.16$ fm*-3*, and $V_{t}$ is the strength parameter. The parameter $\beta$ is usually taken to be one. The parameter $\eta$ controls the density dependence of the pairing functional; $\eta=0$ for volume pairing without density dependence, $\eta=1$ for surface pairing, and $\eta=0.5$ for mixed pairing. The density dependence of the pairing interaction is not fully determined yet, and is an open problem. In the recent optimized Skyrme functionals, values around the mixed pairing is adopted [10, 13].

There are many attempts to extend the pairing functional from the standard one (28). The isoscalar and isovector density dependence are proposed and investigated to improve the behavior of the pairing in the nuclear matter [38] and neutron-rich nuclei [39, 40, 41, 42]. Pairing functionals of the Fayans model [43, 44] with density gradient terms are recently discussed [45].

In this paper, we extend the pairing EDF in a straightforward way. The form of the pairing EDF in equation (28) corresponds to the $t_{0}$ and $t_{3}$ terms in the original Skyrme effective interaction. More terms are derived to the pairing EDF from the original Skyrme effective interaction. In the Skyrme SkP parametrization [6], both the particle-hole and pairing parts are derived from a single Skyrme effective interaction. Although the pairing interaction derived from the Skyrme force is known to be unrealistic, the concept of the generalized DFT allows us to fit the pairing coupling constants independently.

Within the bilinear form of the pair density (with isoscalar particle-hole density dependence), the most general form of the isovector (nn and pp) pairing functional is given by [46]

[TABLE]

where independent coupling constants are assumed for neutrons and protons, and the local pair densities are given from the non-local pair densities $\tilde{\rho}_{t}(\bm{r},\bm{r}^{\prime})$ and $\tilde{\bm{s}}_{t}(\bm{r},\bm{r}^{\prime})$ as

[TABLE]

In this paper, we reduce the number of degrees of freedom to make the discussion simpler. We first assume the local gauge invariance [46, 47, 48], which is the symmetry associated with the local gauge transformation of the wave function

[TABLE]

The local gauge invariance of the pairing EDF provides two conditions to the coupling constants

[TABLE]

We do not consider the tensor pairing functional in this paper for simplicity ($\tilde{C}^{J0}_{t}=\tilde{C}^{J1}_{t}=\tilde{C}^{J2}_{t}=0$). Finally we have a pairing EDF of the following form

[TABLE]

The second term originally corresponds to the momentum-dependent terms with $t_{1}$ parameter in the Skyrme effective interaction.

The pairing potential for the pairing EDF (38) is given by

[TABLE]

5 Pairing rotational moments of inertia in Sn and Pb

We discuss the effect of the pair-density derivative terms introduced in the previous section on the pairing rotational moments of inertia. As typical examples of the pairing rotations, we focus on the neutron pairing rotations in Sn and Pb isotopes.

The calculation for the DFT and FAM for the Thouless-Valatin moments of inertia is performed using the HFBTHO code [49, 50] and its extension to the FAM with the module developed in [51, 22]. Twenty harmonic oscillator major shells together with 60 MeV pairing cutoff energy is used. The parameters for the Gauss quadratures are $N_{\rm GH}=N_{\rm GL}=40$, and $N_{\rm Leg}=80$. The UNEDF1-HFB functional [14] is used for the particle-hole part of the EDF. In the pairing functional (38), there are two pairing coupling constants. In order to see the dependence on the pair-density derivative terms, we use several sets of the coupling constant $\tilde{C}^{\Delta\rho}_{n}=0,-20,-40,-60$, and $-80$ MeV fm5. For each value of $\tilde{C}^{\Delta\rho}_{n}$, we fit the pairing strength $V_{n}$ using the experimental value of the neutron pairing rotational moment of inertia at 120Sn (6.64 MeV*-1*). Three types of density dependence, volume, mixed, and surface are considered. The coupling constants used in this analysis are summarized in table 1.

We analyze the effect of the pair-density derivative coupling constants on the basic quantities, rms radius and HFB energies in 120Sn and 198Pb. For the pairing functionals listed in table 1, the rms radius agrees within 0.03 fm and 0.02 fm for 120Sn and 198Pb, that is about 0.6% and 0.4%, respectively, and the HFB energy agrees within 1.7 MeV and 5.0 MeV that is 0.2% and 0.3%, respectively. The effect of the pair-density derivative coupling constants on these quantities are similar to the effect from different density dependence of the pairing functional. Therefore we focus on the pairing properties in the following.

We begin with the effect of the pair-density derivative terms in 120Sn. Figure 3 shows the neutron pair-density distribution. In the comparison of the isoscalar density dependence of the pairing EDF, the pair density is rather constant inside the nucleus when the volume pairing is employed, while it has a peak in the nuclear surface region with the surface pairing. The pair-density derivative terms in the EDF generally provide a peak in the surface region of the pair density. This effect is similar for volume, mixed and surface pairing, and with stronger pair-density derivative terms, almost similar pair-density distribution are obtained for three kinds of isoscalar density dependence. The same trend is also seen in the pair kinetic density distribution. With the pair-density derivative terms, the sign of the pair kinetic density can change.

The pair-density distribution can be compared with figure 2 in [52], where the neutron pair density for 120Sn is calculated with full SkP interaction (pairing part is from SkP), and the particle-hole SkP interaction with the volume pairing (SkPδ). The pair density obtained by UNEDF1-HFB with volume pairing without pair-density derivative terms is close to the one obtained by SkPδ. The pair density obtained with the full SkP interaction has a strong peak in the surface area. Such a pair density is obtained either with the surface pairing, or with the pair-density derivative terms. We note that in the full SkP results, isoscalar density dependence with $\beta=1/6$ and pair-density derivative terms are included, and moreover, pair tensor density terms are included.

In figure 4 the neutron pairing rotational moment of inertia for Sn and Pb are shown. The calculations without pair-density derivative terms ($\tilde{C}_{n}^{\Delta\rho}=0$) are comparable with the previous results in [19], where the same pairing functionals are used, but the coupling constants are fitted to the neutron pairing gap in 120Sn. The results of the volume and mixed pairing are almost the same as the previous one, because the coupling constants listed in Table 1 are similar to the previous values, while the surface pairing strength fitted to the neutron pairing rotational moment of inertia is smaller than the previous value fitted to the neutron pairing gap ($V_{n}=-474.32$ MeV fm3), and the neutron pairing rotational moment of inertia is systematically larger than the values obtained in  [19]. In all cases, the results without the pair-density derivative terms cannot reproduce the trends in Sn and Pb simultaneously. Theoretical values of the inertia overestimate the experimental ones in the light Sn region around $N=60$ and in the Pb isotopes around $N=110$-$120$. Thus none of the density dependence can reproduce the experimental trend precisely.

The pair-density derivative terms improve this situation. The main effect of the pair-density derivative terms in the neutron pairing rotational moment of inertia in Sn and Pb isotopes is the decrease near $N=60$ and $N=78$ in Sn isotopes, and $N=110$-$120$ in Pb isotopes. Irrespective of the density dependence of the pairing functional, the pair-density derivative terms improve the quantitative agreement of the neutron pairing rotational moment of inertia of both Sn and Pb isotopes.

Next we analyze the pairing gaps. We define two kinds of the average pairing gaps by extending equations (3) and (4) as

[TABLE]

and

[TABLE]

In the case of the present pairing EDF and choice of the phase (normally real representation is used for the pair density), the pairing gap $\Delta_{n}(\tilde{\rho})$ is related to the total pairing energy. Figures 5 and 6 show the average pairing gap $\Delta_{n}(\rho)$ and $\Delta_{n}(\tilde{\rho})$, respectively. The pairing gap $\Delta_{n}(\rho)$ increases with the strength of the pair-density derivative term $\tilde{C}^{\Delta\rho}_{n}$. The correspondence with the experimental OES is the best when pair-density derivative term is absent for the volume and mixed pairing, and at around $\tilde{C}_{n}^{\Delta\rho}=-50$ MeV fm5 for the surface pairing. The pairing gap $\Delta_{n}(\tilde{\rho})$ is more sensitive to the strength of the pair-density derivative term, and it overestimates the experimental OES when $\tilde{C}^{\Delta\rho}_{n}$ is larger than $-20$ MeV fm5.

Two problems arise from this analysis. First, the sensitivity to the pair-density derivative term is different between the two average pairing gaps $\Delta_{n}(\rho)$ and $\Delta_{n}(\tilde{\rho})$, thus it is not possible to determine this coupling constant uniquely using the pairing gaps and experimental OES. Second, in the neutron pairing rotational moment of inertia, the agreement with the experimental data is improved for large $\tilde{C}_{n}^{\Delta\rho}$ such as $-60$ and $-80$ MeV fm5, while the correspondence between the average pairing gap and OES is totally lost for large $\tilde{C}_{n}^{\Delta\rho}$. This indicates that the conventional pairing gaps in the presence of the pair-density derivative term do not correspond to the experimental OES due to the kinetic pairing contribution. On the other hand, the pairing rotational moment of inertia is well defined with the second derivative of the energy with respect to the particle number (or double binding-energy difference of even-even nuclei), and is free from this problem.

As described in equations (42) and (43), the pairing gaps are composed of several terms. To analyze the origin of the differences between the two pairing gaps $\Delta(\rho)$ and $\Delta(\tilde{\rho})$, and deviation of the pairing gaps from the experimental OES, each term in the pairing gaps is plotted separately in figure 7 in the case of the volume pairing. The contribution to the pairing gap from the second term, which includes the pair density coupling with $\Delta\tilde{\rho}_{n}$ is one order smaller than the contributions from the other terms. The conventional first term changes much with the pair-density derivative coupling constant $\tilde{C}^{\Delta\rho}_{n}$. It increases and reaches its maximum at around $\tilde{C}^{\Delta\rho}_{n}=-40$ MeV fm5, and then decreases for larger values of $\tilde{C}^{\Delta\rho}_{n}$. The increase of the pairing gap up to $\tilde{C}^{\Delta\rho}_{n}=-40$ MeV fm5 is due to the increase of the first term. Then for $\tilde{C}^{\Delta\rho}_{n}=-60$ MeV fm5 and $-80$ MeV fm5, the third term which includes the pair kinetic density $\tilde{\tau}_{n}$ suddenly increases, and provides the main contribution to the pairing gap. From the pair kinetic density for 120Sn in figure 3, it can be related to the change of the sign in pair kinetic density.

In $\Delta(\tilde{\rho})$, the third term is proportional to $\tilde{C}^{\tau}_{n}$ with a factor of two, while in $\Delta(\rho)$, the corresponding term is divided into the third and fourth terms. The third terms in $\Delta(\rho)$ and $\Delta(\tilde{\rho})$ are very similar, except for the difference in the factor of two in $\Delta(\tilde{\rho})$. The fourth term that include the particle-hole kinetic density $\tau_{n}$ has opposite sign in $\Delta(\rho)$, and this cancels the contribution from the third term partly in $\Delta(\rho)$. We can say that the difference between the two pairing gaps are mainly from the differences between the particle-hole and pair kinetic densities $\tau_{n}[$ and $\tilde{\tau}_{n}$.

6 Conclusions

We have extended the standard pairing EDF by including the pair-density derivative term that corresponds to the momentum-dependent term in the Skyrme effective interaction. Assuming the local gauge invariance, one coupling constant is introduced for neutrons and protons, respectively. The pairing coupling constants are fitted to reproduce the pairing rotational moments of inertia, which are microscopically computed as the Thouless-Valatin inertia using the finite-amplitude method for nuclear DFT. It is shown that the pair-density derivative term systematically improves the behavior of the neutron pairing rotational moments of inertia in Sn and Pb isotopes. The results indicate that the conventional theoretical average pairing gaps do not serve as measures of the experimental odd-even mass staggering when the pair-density derivative term is present.

The results demonstrate that the pairing rotational moment of inertia is effective for constraining the pairing EDF including the pair-density derivative term. Exploring the pairing tensor terms that are not taken into account in the present analysis is another interesting direction to investigate. Future work should include the systematic and global optimization of the pairing EDF using the pairing rotational moments of inertia.

Acknowledgments

Useful discussions with Witold Nazarewicz are gratefully acknowledged. This work is supported by the JSPS KAKENHI Grand Numbers 16K17680 and 17H05194. Numerical calculation was performed in the resources of High Performance Computing Center, Institute for Cyber-Enables Research, Michigan State University, and the COMA (PACS-IX) System at the Center for Computational Sciences, University of Tsukuba, through the Interdisciplinary Computational Science Program and HPCI System Research Project (hp170144).