Number conserving particle-hole RPA for superfluid nuclei

TL;DR

This paper introduces a number conserving particle-hole RPA method for superfluid nuclei, which smoothly describes collective excitations across the superfluid transition, improving upon traditional RPA approaches.

Contribution

It develops a novel number conserving particle-hole RPA formalism based on antisymmetric geminal powers, bridging the normal and superfluid phases seamlessly.

Findings

Excitations evolve smoothly across the phase transition.

The method aligns closely with exact results in a minimal model.

It avoids the basis change required by traditional RPA methods.

Abstract

We present a number conserving particle-hole RPA theory for collective excitations in the transition from normal to superfluid nuclei. The method derives from an RPA theory developed long ago in quantum chemistry using antisymmetric geminal powers, or equivalently number projected HFB states, as reference states. We show within a minimal model of pairing plus monopole interactions that the number conserving particle-hole RPA excitations evolve smoothly across the superfluid phase transition close to the exact results, contrary to particle-hole RPA in the normal phase and quasiparticle RPA in the superfluid phase that require a change of basis at the broken symmetry point. The new formalism can be applied in a straightforward manner to study particle-hole excitations on top of a number projected HFB state.

Equations

$$|PHFB⟩={\Gamma }^{†M}|0⟩\text{, }{\Gamma }^{†}=\sum _{i=1}^L{z}_i{a}_i^{†}{a}_{\overline{i}}^{†},$$

$$C_{ij}=z_i{a}_i^{†}{a}_j-z_j{a}_{\overline{j}}^{†}{a}_{\overline{i}},\forall i\ne j.$$

$$\left.\begin{array}[c]{c}S_{ij}^{+}=z_{i}a_{i}^{\dagger}a_{\overline{j}}+z_{j}a_{j}^{\dagger}a_{\overline{i}}\\ S_{ij}^{-}=z_{i}a_{\overline{i}}^{\dagger}a_{j}+z_{j}a_{\overline{j}}^{\dagger}a_{i}\end{array}\right\}~{}\forall i>j.$$

$$C_{ij}|PHFB⟩=S_{ij}^{+}|PHFB⟩=S_{ij}^{-}|PHFB⟩=0.$$

$$\left[a_{\alpha }^{†}{a}_{\beta },a_{\gamma }^{†}{a}_{\delta }\right]={\delta }_{\beta \gamma }{a}_{\alpha }^{†}{a}_{\delta }-{\delta }_{\alpha \delta }{a}_{\gamma }^{†}{a}_{\beta },$$

$$a_i^{†}{a}_j$$

$$a_{\overline{i}}^{†}{a}_{\overline{j}}$$

$$a_i^{†}{a}_{\overline{j}}$$

$$[C_{ij},S_{pq}^{+}]$$

$$[C_{ij},S_{pq}^{-}]$$

$$[C_{ij},C_{pq}]$$

$$[S_{ij}^{+},S_{pq}^{-}]$$

$$+z_q(C_{jp}{\delta }_{iq}+C_{ip}{\delta }_{jq}),$$

$$[S_{ij}^{+},S_{pq}^{+}]$$

$$Q_{\mu }^{†}=\sum _{i>j}{X}_{ij}^{\mu }{D}_{ij}^{†}-Y_{ij}^{\mu }{D}_{ij}\text{, }{D}_{ij}=\frac{C_{ij}}{\sqrt{(z_i^2-z_j^2)(n_i-n_j)}},$$

$$Q_{\mu }|RPA⟩=0,$$

$$⟨PHFB|[Q_{\nu },Q_{\mu }^{†}]|PHFB⟩=\sum _{i>j}{X}_{ij}^{\nu }{X}_{ij}^{\mu }-Y_{ij}^{\nu }{Y}_{ij}^{\mu }={\delta }_{\mu \nu }.$$

$$⟨PHFB|[Q_{\nu },Q_{\mu }]|PHFB⟩=\sum _{i>j}{X}_{ij}^{\nu }{Y}_{ij}^{\mu }-X_{ij}^{\nu }{Y}_{ij}^{\mu }=0.$$

$$\sum _{\mu }\left(X_{ij}^{\mu }{X}_{kl}^{\mu }-Y_{ij}^{\mu }{Y}_{kl}^{\mu }\right)={\delta }_{ik}{\delta }_{jl},$$

$$\sum _{\mu }\left(Y_{ij}^{\mu }{X}_{kl}^{\mu }-X_{ij}^{\mu }{Y}_{kl}^{\mu }\right)=0,$$

$$D_{ij}=\sum _{\mu }{X}_{ij}^{\mu }{Q}_{\mu }+Y_{ij}^{\mu }{Q}_{\nu }^{†},$$

$$D_{ij}^{†}=\sum _{\mu }{X}_{ij}^{\mu }{Q}_{\mu }^{†}+Y_{ij}^{\mu }{Q}_{\nu }.$$

$$\left(\begin{array}[]{cc}A&B\\ -B^{\ast}&-A^{\ast}\end{array}\right)\left(\begin{array}[]{c}X\\ Y\end{array}\right)=E_{\mu}\left(\begin{array}[]{c}X\\ Y\end{array}\right),$$

$$A_{ijkl}=⟨PHFB|\left[D_{ij},\left[H,D_{kl}^{†}\right]\right]|PHFB⟩,$$

$$B_{ijkl}=-⟨PHFB|\left[D_{ij}^{†},\left[H,D_{kl}^{†}\right]\right]|PHFB⟩.$$

$$T_{\mu }$$

$$H=J_0-\frac{\Sigma }{(2j-1)}\sum _{\sigma {\sigma }^{\prime }}{A}_{\sigma }^{†}{A}_{{\sigma }^{\prime }}-\frac{\chi }{2(2j-1)}[J_{+}^2+J_{-}^2],$$

$$A_{\sigma }^{†}=\sum _{m=1}^j{a}_{\sigma ,m}^{†}{a}_{\sigma ,-m}^{†},$$

$$A_0^{†}=\sum _{m=1}^j\left(a_{-1m}^{†}{a}_{1,-m}^{†}-a_{-1-m}^{†}{a}_{1,m}^{†}\right)$$

$$J_{+}=\sum _{m=-j}^j{a}_{1m}^{†}{a}_{-1m}=(J_{-}{)}^{†},$$

$$J_0=\frac{1}{2}\sum _{m=-j}^j\left(a_{1m}^{†}{a}_{1m}-a_{-1m}^{†}{a}_{-1m}\right).$$

$$|PHFB⟩={\Gamma }^{†j}|0⟩,{\Gamma }^{†}=z_{-}{A}_{-1}^{†}+z_{+}{A}_1^{†},$$

$$|PHFB⟩=\sum _{l=0}^j\frac{j!}{\left(j-l\right)!l!}{z}_{+}^l{z}_{-}^{\left(j-l\right)}{A}_{+1}^{†l}{A}_{-1}^{†\left(j-l\right)}|0⟩.$$

$$C=z_{-}{J}_{-}-z_{+}{J}_{+}.$$

Full text

††journal: Physics Letter B

1 Introduction

Excitation spectra of nuclei show a tremendous richness and diversity, going from low energy collective states to giant resonances. The basic approach for the description of those phenomena is time-dependent mean field with the use of more or less sophisticated energy density functionals (EDFs). In the small amplitude limit, this means the consistent solution of mean field and RPA equations. For superfluid nuclei this has to be generalized to Hartree-Fock-Bogoliubov (HFB) and quasiparticle RPA (QRPA) equations. Many works of this type based on Skyrme or Gogny EDFs as well as relativistic HFB exist in the literature, see for example [1, 2, 3, 4, 5, 6]. The weakness of the HFB approach is that it violates, among other symmetries, particle number conservation. This is no problem in macroscopic superfluid systems. However, in finite nuclei this becomes annoying if there is a spread of several units of particles in the calculations where one desires a sharp number for each nucleus. The restoration of good particle number is technically feasible with modern projection techniques [7, 8] but becomes rather cumbersome on the level of QRPA where one has to project two-quasiparticle states with a subsequent orthogonalization [9]. This might be the reason why only in relatively rare cases it has been applied mainly to $\beta$ and double-$\beta$ decays [10, 11]. Moreover, number projected QRPA, being a projection of two quasiparticle states, involves a mixture of particle-hole (ph) excitations of the $A$ nucleus, particle-particle excitations of the $A-2$ nucleus and hole-hole of the $A+2$ nucleus. At the end one, therefore, one has to filter out the most interesting excitations of the $A$ system which are then excitations of the ph type. These problems are particularly important for transitional regions since there the particle number fluctuations are strong and the transition from normal to superfluid is rounded due to the finiteness of the system. Crossing the superfluid transition requires a change in the reference state from a HF Slater determinant to a HFB quasiparticle vacuum. The collective ph excitations described by RPA in the former case and by QRPA in the latter yield a sharp transition with a significant kink [12].

Collective states with projections, in particular particle number projection, have been investigated with the Generator Coordinate Method (GCM). However, GCM is generally tailored for large amplitude collective motion and only a few collective coordinates are considered which are obtained from constrained and symmetry projected HFB calculations. RPA theory which is considered in the present work, on the contrary, treats collective motion in the small amplitude limit where collective and non-collective states are generally treated on the same footing. A recent paper with many references to the GCM-PHFB method can be found in [13]. In this situation, we thought it useful to adapt a formalism, known in quantum chemistry since many years [14], to the nuclear physics case. While in quantum chemistry this formalism has had little success due the absence of attractive pairing correlations, for the same reason we do believe that it could be extremely useful in nuclear structure. In chemical physics it was shown [15] that the so-called antisymmetrized geminal power (AGP) state, which is the same as the well known number projected HFB (PHFB) ground state, is the vacuum of a set of ph operators which kill the PHFB vacuum. This remarkable property can be exploited to construct a particle-hole RPA (phRPA) approach along the same lines as in the standard RPA theory. One adds to the particular set of ph killers which play the role of backward terms, their adjoints as the forward terms with the corresponding amplitudes and establishes with the equation of motion technique the RPA matrix. We call this theory number conserving particle-hole RPA (NCphRPA). This new particle number conserving technique may be implemented without too much difficulty and effort in the existing QRPA codes to collective modes in superfluid nuclei like, e.g., the chain of Sn isotopes. Since NCphRPA adds correlations on top of standard ph-RPA even at magicity, one could consider all Sn isotopes from A = 100 to A = 132 and beyond with the same approach.

We will first classify the complete set of ph operators in terms of killers of an arbitrary PHFB state, their adjoints and diagonal operators. We will then proceed to derive the NCphRPA based on the killers and adjoints in complete analogy to the RPA over a HF reference state. Finally, we will test numerically this approach in a minimal model of pairing plus monopole ph interactions.

2 PHFB killers

Let us start by defining the PHFB wave function of $M$ nucleon pairs as a general pair condensate in the canonical basis where the density matrix is diagonal

[TABLE]

where the operator $a_{i}^{\dagger}$ ($a_{\overline{i}}^{\dagger}$) create a particle in the state $i$ ($\overline{i}$). $\left\{i,\overline{i}\right\}$ is a single particle basis with $2L$ states that has been grouped into $L$ pairs $i$ and $\overline{i}$. The $\overline{i}$ states are conjugate orbitals to the ${i}$ single particle states. These pairs are not necessarily degenerated but have the same occupation probabilities $\langle a_{i}^{\dagger}a_{i}\rangle=\langle a_{\overline{i}}^{\dagger}a_{\overline{i}}\rangle$. The free parameters $z_{i}$ can eventually be determined variationally as for example in HFB with number projected before variation. Similar to a HF Slater determinant where it is possible to define a set of killers in terms of which the RPA operators are constructed, PHFB has a set of killers that we classify into particle killers and spin killers.

Particle killers:

[TABLE]

• 2.

Spin killers:

[TABLE]

Note that particle killers move nucleons without changing the spin, while spin killers move nucleons flipping the spin. It is easy to check that all killers commute with the pair operator (1), $[C_{ij},\Gamma^{\dagger}]=[S_{ij}^{+},\Gamma^{\dagger}]=[S_{ij}^{-},\Gamma^{\dagger}]=0$, and therefore, each one of the killers annihilates the $PHFB$ wave function (1)

[TABLE]

The set of killers (2-3), their adjoints and the diagonal operators $a_{i}^{\dagger}a_{i}$, $a_{\overline{i}}^{\dagger}a_{i}$, $a_{i}^{\dagger}a_{\overline{i}}$ and $a_{\overline{i}}^{\dagger}a_{\overline{i}}$ form a complete set of $4L^{2}$ ph operators that generate a $U(2L)$ algebra of particle-hole operators defined by the commutators

[TABLE]

where $\alpha,\beta,\gamma,\delta$ stand for any orbital $j$ and their conjugate, $\overline{j}$. Within this set there are a total of $L(L-1)$ particle killers and the same number of spin killers. The total number of adjoints that create particle-type or spin-type excitations is $2L(L-1)$. Adding up the killers, the adjoints and the $4L$ diagonals, we obtain the total of $4L^{2}$ ph operators. The completeness of this set guaranties that any ph operator can be expressed as a linear combination of the above generators

[TABLE]

The denominators that appear in the inversion relations could be the source of numerical instabilities in case of level crossings or for states deep in the Fermi sea with occupations close to 1. In the former case, the pair of approaching levels could be considered as effectively degenerated and redefine de killers taking into account this degeneracy. In the latter case, it might be possible to consider these states as an inert core and exclude them from the numerical calculation. On the other hand, degeneracies due to the conservation of symmetries like the rotational symmetry (spherical nuclei) can be explicitly taken into account in the definition of the killers.

Using the commutation algebra of the ph operators (5) it is easy to derive the commutation relations between killers

[TABLE]

The complete set of commutators including the adjoints and diagonals follows in a similar manner.

3 Number conserving particle-hole RPA

We will derive the NCphRPA equations following closely the equation of motion method of [16, 17] assuming a PHFB pair condensate (1) as the reference state. In order to ease the formalism we restrict to particle excitations, but the generalization to spin excitations is straightforward. In fact, if the PHFB state is axially symmetric, excitations with a given $k$ value of the $z$ component of the total angular momentum could be constructed as a mixture of particle and spin excitations with the same $k$.

As usual, we assume that excited states $\left|\mu\right\rangle=Q_{\mu}^{\dagger}\left|0\right\rangle$ are generated by collective RPA operators

[TABLE]

where the square root is a normalization factor and $n_{i}=\langle a_{i}^{\dagger}a_{i}\rangle$ is the occupation probability in the PHFB state. With this definition, the backward term annihilates the PHFB state in complete analogy to the RPA based on HF.

Moreover, the RPA assumes a killing condition

[TABLE]

for the correlated RPA ground state which later will be relaxed replacing it by the PHFB reference state.

We can now check that the RPA operators (8) are properly normalized in the PHFB state

[TABLE]

As well as

[TABLE]

Together with the closure relations

[TABLE]

[TABLE]

it allows us to invert the RPA operators

[TABLE]

[TABLE]

The equation of motion method leads to the RPA equation

[TABLE]

with

[TABLE]

[TABLE]

Given a transition operator $T$, its transition matrix element $T_{\mu}=\left\langle 0\right|T\left|\mu\right\rangle$ can be evaluated as

[TABLE]

In NCphRPA we evaluate expectation values of the $A$ and $B$ matrices and the transition matrices in the PHFB ground state in a similar way as the standard RPA replaces the correlated vacuum defined by the killing condition by the mean-field reference state, either a HF Slater determinant or a quasiparticle vacuum. However, the inversion relations (6) and (14-15) would allow us to go beyond NCphRPA towards self-consistent RPA [18] evaluating most of the expectation values in the true RPA vacuum using the killing condition.

4 Benchmark with the Agassi model

As a minimal model to benchmark the NCphRPA approximation, we will consider a two-level model with pairing and particle-hole monopole interactions. The model was introduced by Agassi [19] as a tool to test many-body theories that deal with the interplay between ph and superfluid correlations. It has been studied at the mean-field level by Agassi himself [19] and later by Davis and Heiss [20]. They derived the phase diagram of the model and the collective excitations in each of the phases by means of HFB and, phRPA and QRPA approximations respectively. More advanced many-body methods like the merging of Coupled Cluster with symmetry restored HFB theory [21] revived the model as an excellent test-bed. Recently, García-Ramos et al [22] generalized the model by the introduction of new interaction terms that give rise to an extremely rich phase diagram. We will use here the original Agassi model described by a Hamiltonian which is a superposition of the Lipkin model and the two-level pairing model

[TABLE]

where $\sigma=\pm 1$ labels each of the two single particle levels with degeneracy $2j$, and $\Sigma$ and $\chi$ are suitably scaled coupling parameters for pairing and ph interactions. The pair creations operators are

[TABLE]

and the ph operators are

[TABLE]

[TABLE]

Fig. 1 shows the phase diagram of the Agassi model at half filling [20] with a normal (spherical) phase for $\chi<1$ and $\Sigma<1$, a parity broken (deformed) phase for $\chi>1$ and $\chi>\Sigma$, and a superfluid phase for $\Sigma>1$ and $\Sigma>\chi$, for details see [20, 22]. Interestingly, the model has no stable phase with both symmetries broken. We will test and compare the different RPA approximations along the horizontal line $\chi=\frac{1}{2}$ which contains sizeable ph correlations. As a function of $\Sigma$, it crosses the transition from normal to superfluid requiring a change of the reference state from a Slater determinant (HF) to a quasiparticle vacuum (BCS). Since the line $\chi=\frac{1}{2}$ runs across the normal and superfluid regions that preserve parity, we use the following PHFB ansatz at half filling

[TABLE]

where the structure parameters $z_{\pm}$ can be obtained by a straightforward minimization of the energy or by number projected HFB calculation before variation. In our case, we expand (25) using a binomial

[TABLE]

Note that for $z_{+}=0$ we recover the HF Slater determinant with all particles filling the lower level.

PHFB improves over the mean field used to construct the phase diagram taking into account some of the correlations even in the normal phase.

Due to the simplicity of the model, we are left with a couple of ph operators $J_{+}$ and $J_{-}$, in terms of which the killer of the pair condensate (25) is

[TABLE]

It can be easily checked that $[C,\Gamma^{\dagger}]=0$ and therefore, $C$ annihilates the pair condensate (25). Following the general derivation of the NCphRPA, the ph RPA operator is

[TABLE]

The RPA matrices $A$ and $B$ as defined in (17) and (18) can be obtained from the double commutators

[TABLE]

and

[TABLE]

by taking expectation values in the PHFB state (26).

It can be easily checked that NCphRPA reduces to phRPA if we replace the PHFB reference state (26) by the Slater determinant obtained in the limit $z_{+}=0$ with all particles occupying the lower level.

In Fig. 2 we show the excitation energy

$\omega=\sqrt{A^{2}-B^{2}}$ along the line of $\chi=0.5$, as a function of $\Sigma$. It is compared with the phRPA for $0<\Sigma<1$ in the spherical regions and with QRPA for $\Sigma>1$ in the superfluid region [23]. We can see a clear improvement of NCphRPA over both RPA’s, specially around the phase transition which NCphRPA crosses smoothly, while in the conventional RPA one has to change from a number conserving HF basis to a number non-conserving quasiparticle vacuum with a consequent kink. This improvement could be traced back to the structure of the PHFB reference state shown in the inset. This parameter $z=\frac{z_{+}}{z_{-}}$ displays a smooth behavior across the transition due to the inclusion of pair correlations and number fluctuations absent at the mean-field level.

In order to test the RPA wave function we evaluate the transition matrix element of the operator $J_{x}=\frac{1}{2}(J_{+}+J_{-})$.

[TABLE]

In Fig. 3 we show this transition matrix element. Again we see a clear improvement of the present approach over the conventional RPA’s [23], mainly in the transitional region. The properties of the excited state as described by the NCphRPA reflect the net improvement seen in the ground state at the level of PHFB as seen in the inset of Fig. 2 as well as in the inset of Fig. 3 for the expectation value of $J_{0}$.

5 Conclusions

In this work we have introduced the NCphRPA following a similar approach presented long ago in quantum chemistry. NCphRPA starts with a linear transformation of the one particle operators such that they can be classified into killers of a PHFB state, their adjoints, and diagonals. Based on this new one particle operators, NCphRPA is derived in a standard way using the equations of motion method. The formalism is a straightforward extension of phRPA based on a HF determinant as a reference state, and it requires the sole information of the one and two-body density matrices in the optimized PHFB state. In this regard, it is much simpler than the number projected QRPA that requires the projection of a complete set of two quasiparticle states and the subsequent re-orthonormalization. The general formalism of NCphRPA developed and tested here is ready to be used with density-dependent forces treated with PHFB to describe collective motion in transitional and superfluid nuclei. Even for magic nuclei, NCphRPA with a PHFB reference state accounting for some of the pair fluctuations, would give non-trivial corrections to phRPA based on a Slater determinant. Therefore, NCphRPA can be used with benefit throughout the mass table. These kinds of corrections to phRPA can be seen in figures 2 and 3 for $0<\Sigma<1$. We want also to stress that the transition from non-superfluid to superfluid nuclei is always completely smooth since the reference state (PHFB) does not show any abrupt change in the crossover. The fact that the NCphRPA uses the equation of motion method to derive the RPA equations allows for an improved treatment of the reference in the direction of a number conserving self-consistent RPA, by using the RPA killing condition and the inversion of the RPA operators. Before closing we would like to emphasize the ph character of the NCphRPA. Within this formalism particle-particle correlation are only described at the level of the PHFB reference state. However, particle-particle and hole-hole fluctuations could be taken into account in a second RPA based on quadratic killers and quadratic adjoints.

6 Acknowledgements

This work has been partially supported by the Consejería de Economía, Conocimiento, Empresas y Universidad de la Junta de Andalucía (Spain) under Groups FQM-160 and FQM-370 and by European Regional Development Fund (ERDF), ref.  SOMM17/

6105/UGR. We acknowledge financial support from the Spanish Ministerio de Ciencia, Innovación y Universidades and the ERDF under Projects No. FIS2015-63770-P, FIS2017-88410-P and PGC2018-094180-B-I00. Resources supporting this work were provided by the CEAFMC and Universidad de Huelva High Performance Computer (HPC@UHU) funded by

FEDER/MINECO project UNHU-15CE-2848.