Nuclear response in a finite-temperature relativistic framework

TL;DR

This paper develops a finite-temperature relativistic framework for nuclear response using the Bethe-Salpeter equation and Matsubara Green's functions, enabling analysis of nuclear oscillations at various temperatures.

Contribution

It introduces a temperature-dependent projection operator to simplify the Bethe-Salpeter equation at finite temperature within a relativistic quantum hadrodynamics model.

Findings

Successfully applied to monopole, dipole, and quadrupole responses of $^{48}$Ca.

Analyzed the evolution of nuclear collective oscillations with temperature.

Demonstrated the method's effectiveness in finite-temperature nuclear response studies.

Abstract

A thermal extension of the relativistic nuclear field theory is formulated for the nuclear response. The Bethe-Salpeter equation (BSE) with the time-dependent kernel for the particle-hole response is treated within the Matsubara Green's function formalism. We show that, with the help of a temperature-dependent projection operator on the subspace of the imaginary time (time blocking), it is possible to reduce the BSE for the nuclear response function to a single frequency variable equation also at finite temperature. The approach is implemented self-consistently in the framework of quantum hadrodynamics based on the meson-nucleon Lagrangian. The method is applied to the monopole, dipole and quadrupole response of $^{48}$Ca and to the dipole response of the tin isotopes $^{100,120,132}$Sn, in particular, to a study of the evolution of nuclear collective oscillations with temperature. The article is dedicated to the memory of Pier Francesco Bortignon and devoted to the developments related to his pioneering ideas.

Equations

$$\Omega (\lambda ,T)=E-TS-\lambda N$$

$$S=-k\text{mbox}Tr(\hat{\rho }\text{mbox}ln\hat{\rho }),N=\text{mbox}Tr(\hat{\rho }\hat{\mathcal{N}}).$$

$$E[\hat{\rho },{\varphi }_m]$$

$$\hat{\rho }=\frac{e^{-(\hat{\mathcal{H}}-\lambda \hat{\mathcal{N}})/kT}}{\text{mbox}Tr[e^{-(\hat{\mathcal{H}}-\lambda \hat{\mathcal{N}})/kT}]},$$

$$n_1(T)=n({\epsilon }_1,T)=\frac{1}{1+\text{mbox}exp\{{\epsilon }_1/T\}},$$

$$N_{\mu }(T)\equiv N({\Omega }_{\mu },T)=\frac{1}{\text{mbox}exp\{{\Omega }_{\mu }/T\}-1},$$

$$\mathcal{R}(12,34)={\tilde{\mathcal{R}}}^0(12,34)+$$

$$\displaystyle+\sum_{5678}^{\tau}{{\widetilde{\mathcal{R}}}}^{0}(12,56)\Big{[}{\widetilde{V}}(56,78)+{\mathcal{W}}(56,78)\Big{]}{\mathcal{R}}(78,34),$$

$$\sum _{\mathrm{12..}}^{\tau }=\sum _{k_1{k}_2...}{\int }_0^{1/T}d{\tau }_{1}d{\tau }_2\cdots$$

$$\tilde{\mathcal{G}}(2,1)=\sum _{\sigma }{\tilde{\mathcal{G}}}^{\sigma }(2,1),$$

$${\tilde{\mathcal{G}}}^{\sigma }(2,1)=-\sigma {\delta }_{12}n[-\sigma {\epsilon }_1,T]e^{-{\epsilon }_1{\tau }_{21}}\theta (\sigma {\tau }_{21}),$$

$${\tilde{\mathcal{G}}}_{21}({\epsilon }_{\ell })=\frac{1}{2}{\int }_{-1/T}^{1/T}d{\tau }_{21}{e}^{i{\epsilon }_{\ell }{\tau }_{21}}\tilde{\mathcal{G}}(2,1)$$

$$=\frac{{\delta }_{12}}{i{\epsilon }_{\ell }-{\epsilon }_1+\mu }={\delta }_{12}{\tilde{\mathcal{G}}}_1({\epsilon }_{\ell }),$$

$${\mathcal{R}}^e={\tilde{\mathcal{R}}}^0+{\tilde{\mathcal{R}}}^0\mathcal{W}{\mathcal{R}}^e$$

$$\mathcal{R}={\mathcal{R}}^e+{\mathcal{R}}^e\tilde{V}\mathcal{R}$$

$$\mathcal{W}(12,34)={\tilde{\mathcal{G}}}^{-1}(3,1){\Sigma }^e(2,4)$$

$${\Sigma }_{12}^e({\epsilon }_{\ell }+{\omega }_n)$$

$${\mathcal{R}}^e$$

$$\Theta (12,34)$$

$$\Theta (12,34;T)={\delta }_{{\sigma }_1,-{\sigma }_2}\theta ({\sigma }_1{\tau }_{41})\theta ({\sigma }_1{\tau }_{32})\times$$

$$\times [n({\sigma }_1{\epsilon }_2,T)\theta ({\sigma }_1{\tau }_{12})+n({\sigma }_2{\epsilon }_1,T)\theta ({\sigma }_2{\tau }_{12})],$$

$$\tilde{\mathcal{D}}(12,34)=\Theta (12,34;T){\tilde{\mathcal{R}}}^{0({\sigma }_1{\sigma }_2)}(12,34)=$$

$$={\delta }_{{\sigma }_1,-{\sigma }_2}{\tilde{\mathcal{G}}}^{{\sigma }_1}(3,1){\tilde{\mathcal{G}}}^{{\sigma }_2}(2,4)\theta ({\sigma }_1{\tau }_{41})\theta ({\sigma }_1{\tau }_{32})\times$$

$$\times [n({\sigma }_1{\epsilon }_2,T)\theta ({\sigma }_1{\tau }_{12})+n({\sigma }_2{\epsilon }_1,T)\theta ({\sigma }_2{\tau }_{12})]$$

$${\mathcal{R}}_{12,34}^e$$

$${\mathcal{R}}_{12,34}^e({\omega }_n)=T^2\sum _{\ell }\sum _{{\ell }^{\prime }}{\mathcal{R}}_{12,34}^e({\omega }_n,{\epsilon }_{\ell },{\epsilon }_{{\ell }^{\prime }}),$$

$${\mathcal{R}}_{12,34}^e(\omega )={\tilde{\mathcal{R}}}_{12,34}^0(\omega )+\sum _{56,78}{\tilde{\mathcal{R}}}_{12,56}^0(\omega ){\Phi }_{56,78}(\omega ){\mathcal{R}}_{78,34}^e(\omega ),$$

$${\tilde{\mathcal{R}}}_{12,34}^0(\omega )={\delta }_{13}{\delta }_{24}\frac{n({\epsilon }_2,T)-n({\epsilon }_1,T)}{\omega -{\epsilon }_1+{\epsilon }_2}$$

$${\Phi }_{12,34}(\omega )=\frac{{\delta }_{{\sigma }_1,-{\sigma }_2}{\sigma }_1}{n({\epsilon }_4,T)-n({\epsilon }_3,T)}\times$$

$$\times \sum _{56;\mu }\sum _{{\eta }_{\mu }=\pm 1}{\eta }_{\mu }{\zeta }_{12,56}^{\mu {\eta }_{\mu }}{\zeta }_{34,56}^{\mu {\eta }_{\mu }\ast }\times$$

$$\times \frac{[N({\eta }_{\mu }{\Omega }_{\mu })+n({\epsilon }_6,T)][n({\epsilon }_6-{\eta }_{\mu }{\Omega }_{\mu },T)-n({\epsilon }_5,T)]}{\omega -{\epsilon }_5+{\epsilon }_6-{\eta }_{\mu }{\Omega }_{\mu }}.$$

$${\zeta }_{12,56}^{\mu {\eta }_{\mu }}$$

$$g_{12}^{\mu (\sigma )}={\delta }_{\sigma ,+1}{g}_{12}^{\mu }+{\delta }_{\sigma ,-1}{g}_{21}^{\mu \ast }.$$

$$g_{12}^{\mu }=\sum _{34}{\tilde{V}}_{12,34}{\rho }_{34}^{\mu },$$

Full text

1 Introduction

Response to external probes is one of the major characteristics of strongly-correlated systems, which provides complete information about them. In nuclear systems, the spectral properties extracted from various nuclear reactions serve also for understanding the underlying nuclear forces. Being extremely complicated already for a two-nucleon system, the nucleon-nucleon interaction is, in addition, considerably modified by numerous processes in the strongly-correlated medium, which are commonly called medium polarization. It was recognized early by A. Bohr and B.R. Mottelson BohrMottelson1969 ; BohrMottelson1975 that the major contribution to the medium polarization originates from the coupling between nucleonic and vibrational degrees of freedom, at least in medium-mass and heavy nuclei. Since then, this idea has advanced many directions in both theoretical and experimental nuclear physics. Pier Francesco Bortignon was one of the pioneers of the nuclear field theory (NFT) which developed the concept of the interaction between fermionic and bosonic (vibrational) degrees of freedom in atomic nuclei Broglia1976 ; BortignonBrogliaBesEtAl1977 ; Bortignon1978 ; Bortignon1981 ; BertschBortignonBroglia1983 ; MahauxBortignonBrogliaEtAl1985 ; Bortignon1986 ; Bortignon1997 ; ColoBortignon2001 . These developments lead to a very successful description of a large amount of experimental data and advanced tremendously the understanding of many nuclear structure phenomena, especially those related to collective effects. They also inspired experimental groups in various fields of research, such as giant resonances at finite temperature and multi-phonon excitations. This has become possible after a number of works on single-particle and collective excitations in hot nuclei within the thermal NFT Bortignon1986 ; BrogliaBortignonBracco1992 ; Giovanardi1996 ; Donati1996 ; Giovanardi1998 as well as on the theory of thermal shape fluctuations OrmandBortignonBrogliaEtAl1990 ; OrmandBortignonBroglia1996 . A large corpus of theoretical knowledge and experimental data was summarized in the volume ”Giant Resonances: Nuclear Structure at Finite Temperature” by Pier Francesco Bortignon and his co-authors Angela Bracco and Ricardo A. Broglia Bortignon1998 . In his late years, he was involved in new advancements of the NFT Bortignon2016 ; Bortignon2016a ; Broglia2016a ; Colo2017 , which also open new horizons for future research.

In this article dedicated to his memory we would like to focus on one of our recent developments which are closely related to the pioneering ideas of Pier Francesco Bortignon on the nuclear response at finite temperature. Inspired by those, the temperature dependence of the nuclear response, mostly in the dipole channel and mostly of its high-energy part associated with the giant dipole resonance (GDR), was extensively studied experimentally in the past Gaardhoje1984a ; GaardhojeEllegaardHerskindEtAl1986 ; Bracco1989 ; Ramakrishnan1996 ; Mattiuzzi1997 ; Heckman2003 , see also a relatively recent review Santonocito2006 and a newer study of Refs. Wieland2006 ; Bracco2007 . In more recent experiments on the dipole response, a concentration of strength has been observed in the low-energy region SavranAumannZilges2013 ; Isaak2019 , being most prominent in neutron-rich nuclei. Its origin was attributed to the coherent oscillation of the neutron skin around the isospin-saturated core, which is called pygmy dipole resonance (PDR). In the gamma decay of quasicontinuum excitations another enhancement of the dipole strength has been found at very low gamma energies VoinovAlginAgvaanluvsanEtAl2004 ; WiedekingBernsteinKrtickaEtAl2012 . These new features of the low-energy strength functions become especially interesting in the context of astrophysical modeling, because they can, in particular, significantly enhance the reaction rates of the rapid neutron-capture nucleosynthesis (r-process) MumpowerSurmanMcLaughlinEtAl2016 and modify the nuclear matter equation of state. Similar enhancements of the low-lying spin-isospin response are found in the theoretical finite-temperature studies of Refs. Minato2009 ; Dzhioev2010 ; Niu2011 , which indicate that the electron capture and beta decay rates in the astrophysical environments are also affected by the sensitivity of nuclear transitions to finite temperature. Thus, both phases of the r-process (neutron capture and beta decay), which occurs in the neutron star mergers (NSM), and electron capture, which is crucial for the core-collapse supernovae (CCSN) simulations, require precise knowledge about low-energy nuclear transitions at finite temperature. The nature and progress of the r-process depend on the specific astrophysical scenario, which can vary in rather broad limits, however, a consistent description of merging neutron stars require the information about the CCSN which precedes the neutron star formation. As discussed in Ref. Aprahamian2018 , ultimately, identical simulation frameworks are to be used for both the NSM and the CCSN.

As the complete description of the r-process requires information about many nuclei which are beyond current experimental capabilities, a reliable theory of the finite-temperature nuclear response is highly desirable. An accurate theory for the response of nuclei at finite temperature, or compound nuclei, which is required for applications, is very challenging. The most common practice for the microscopic approaches to the finite-temperature nuclear response is to confine the framework by the simplest one-loop approximation represented by the thermal random phase approximation (TRPA) RingRobledoEgidoEtAl1984 ; NiuPaarVretenarEtAl2009 , thermal quasiparticle RPA (TQRPA) RingRobledoEgidoEtAl1984 ; YuekselColoKhanEtAl2017 , and by the continuum TQRPA LitvinovaKamerdzhievTselyaev2003 ; KhanVanGiaiGrasso2004 ; LitvinovaBelov2013 . Some extensions beyond T(Q)RPA have been proposed over the years, for instance, the finite- temperature self-consistent RPA DukelskyRoepkeSchuck1998 , the collision-integral approach LacroixChomazAyik1998 , the finite-temperature second RPA (SRPA) AdachiSchuck1989 ; Yannouleas1986 , the thermal NFT with a coupling of nucleons to collective surface vibrations Bortignon1986 , and the quasiparticle-phonon model formulated in terms of the thermofield dynamics Storozhenko2004 . These microscopic approaches elaborated on different aspects of the finite-temperature nuclear response and, together with the theories of thermal shape fluctuations Kusnezov1998 ; OrmandBortignonBroglia1996 advanced considerably the understanding of spectral properties of hot nuclei.

In this article we discuss another microscopic approach to the finite-temperature nuclear response, which was developed recently in a relativistic framework. The aim of this approach is to advance the relativistic version of the nuclear field theory (RNFT) LitvinovaRing2006 ; LitvinovaRingTselyaev2007 ; LitvinovaRingTselyaev2008 ; LitvinovaRingTselyaev2013 ; Litvinova2015 ; Litvinova2016 ; RobinLitvinova2016 ; RobinLitvinova2018 to the finite- temperature case. The RNFT, which is based on the covariant energy density functional of quantum hadrodynamics, extends the relativistic RPA (RRPA) by the (quasi)particle-vibration coupling in a parameter-free way. This approach performs very well in the description of nuclear transitions from ground to excited states LitvinovaRingTselyaevEtAl2009 ; LitvinovaLoensLangankeEtAl2009 ; LitvinovaRingTselyaev2010 ; EndresLitvinovaSavranEtAl2010 ; TamiiPoltoratskaNeumannEtAl2011 ; MassarczykSchwengnerDoenauEtAl2012 ; LitvinovaRingTselyaev2013 ; SavranAumannZilges2013 ; LanzaVitturiLitvinovaEtAl2014 ; PoltoratskaFearickKrumbholzEtAl2014 ; Oezel-TashenovEndersLenskeEtAl2014 ; EgorovaLitvinova2016 . It has been extended recently in Refs. LitvinovaWibowo2018 ; WibowoLitvinova2018 to finite temperatures and applied to the dipole response of selected medium-light and medium-heavy nuclei. In the present article we discuss some more details of the dipole response as well as further aspects of the finite-temperature theory, such as the monopole and quadrupole excitations and its broader impacts on the r-process and on the nuclear matter equation of state. We will consider both the low-temperature range (1-2 MeV) relevant for the astrophysical applications, such as NSM and CCSN, and higher temperatures which can be achieved in heavy-ion collisions.

2 Method

We begin with modeling a compound nucleus with the finite-temperature relativistic mean-field (RMF) theory. The grand potential $\Omega(\lambda,T)$

[TABLE]

is minimized with the Lagrange multipliers $\lambda$ and $T$ determined by the average energy $E$, particle number $N$, and the entropy $S$ Sommermann1983 . The quantities $N$ and $S$ are thermal averages with the one-body nucleonic density operator $\hat{\rho}$ :

[TABLE]

Here $\hat{\cal N}$ is the particle number operator, $k$ is the Boltzmann constant, and the energy $E$ is a covariant functional of the nucleonic density and classical meson and photon fields $\phi_{m}$ VretenarAfanasjevLalazissisEtAl2005 :

[TABLE]

The energy is thus determined by the nucleon mass $M$ and non-linear sigma-meson potentials $U(\phi_{m})$ containing meson masses and meson-nucleon coupling vertices adjusted to bulk nuclear characteristics Lalazissis1997 ; VretenarAfanasjevLalazissisEtAl2005 . In Eq. (3) the sign ”+” is associated with the scalar $\sigma$-meson, ”-” with the vector $\omega$- and $\rho$-mesons, and the index ”$m$” runs over the mesons, photon and Lorentz indices VretenarAfanasjevLalazissisEtAl2005 . The variation of Eq. (1) determines the operator of the nucleonic density with the unity trace

[TABLE]

where ${\hat{\cal H}}=\delta E[{\hat{\rho}},\phi_{m}]/\delta{\hat{\rho}}$. The eigenvalues of the density matrix are the Fermi-Dirac occupation factors:

[TABLE]

where the subscript ”1” runs over the complete set of the single-particle quantum numbers in the Dirac-Hartree basis, which diagonalizes ${\hat{\cal H}}$ with the eigenvalues of the single-particle energies $\varepsilon_{1}={\tilde{\varepsilon}_{1}}-\lambda$ measured from the chemical potential $\lambda$. In Eq. (5) and in the following we adopt the units where the value of the Boltzmann constant is $k=1$. In the following, we will deal also with the counterpart of Eq. (5) in the bosonic sector:

[TABLE]

where $\Omega_{\mu}$ are the frequencies of the phonon modes composed of correlated particle-hole excitations, which possess the quantum numbers $\mu$.

As in Refs. LitvinovaWibowo2018 ; WibowoLitvinova2018 , in this work we will consider non-superfluid nuclear systems, such as doubly-magic nuclei and nuclei at temperatures above the critical temperature, when pairing correlations vanish. In this case, the particle-hole finite-temperature four-point response function ${\mathcal{R}}(12,34)$ can be conveniently formulated in the Dirac-Hartree basis of the single-particle states defined above and in the imaginary-time representation, so that the number multi-indices $1=\{k_{1},\tau_{1}\}$ with the $k_{1}$ characterizing the single-nucleon quantum numbers and $\tau_{1}$ being the imaginary time. The response function ${\mathcal{R}}(12,34)$ obeys the Bethe-Salpeter equation (BSE)

[TABLE]

which is formally similar to the one at zero temperature, but in the imaginary-time representation. In Eq. (2) the summations imply integrations over the time arguments:

[TABLE]

and the uncorrelated particle-hole response ${\widetilde{\mathcal{R}}}^{0}(12,34)=\\ {\widetilde{\mathcal{G}}}(3,1){\widetilde{\mathcal{G}}}(2,4)$ is defined via the one-body Matsubara Green’s function ${\widetilde{\mathcal{G}}}$ in the thermal mean-field Abrikosov1965 :

[TABLE]

where $\tau_{21}=\tau_{2}-\tau_{1}$ ($-1/T<\tau_{21}<1/T$), $\theta(\tau)$ is the Heaviside step-function and the index $\sigma=+1(-1)$ denotes the retarded (advanced) component of $\widetilde{\mathcal{G}}$. The Fourier transformation of Eqs. (9,10) with respect to the imaginary time gives the mean-field Matsubara Green’s function in the domain of the discrete imaginary energy variable:

[TABLE]

where $\varepsilon_{\ell}=(2\ell+1)\pi T$ with $\ell$ being integer numbers.

The interaction kernel of the BSE (2) is divided into two parts $\widetilde{V}$ and $\cal W$, which are responsible for the short and long-range correlations, respectively. This division is very general and can be obtained, for instance, by the equation of motion method for a model-independent Hamiltonian with a time-independent two-body bare interaction AdachiSchuck1989 ; DukelskyRoepkeSchuck1998 ; SchuckTohyama2016 . As in our previous works LitvinovaRingTselyaev2007 ; LitvinovaRingTselyaev2008 ; LitvinovaRingTselyaev2010 ; LitvinovaRingTselyaev2013 , the short-range contribution $\widetilde{V}$ is associated with the exchange by effective mesons and photon parametrized in accordance with Eq. (3) Lalazissis1997 , and the long-range term $\cal W$ is induced by the medium polarization effects. In practice, the short-range term is almost instantaneous, and its temporal locality translates to the energy independence of its Fourier transform. In case of neglecting the second long-range term, where the retardation effects are important, the BSE reduces to the response-theory form of the random phase approximation (RPA). If the latter term is retained, it translates to an energy-dependent kernel while the time integrations transform to the integrations over the energy variables of the internal fermionic loops. Because of a clear separation of these two terms, it is convenient to divide the problem of solving Eq. (2) into two parts, namely: (i) to calculate the correlated propagator $\mathcal{R}^{e}$ from the equation (in the operator form):

[TABLE]

and (ii) to solve the remaining equation

[TABLE]

for obtaining the full response $\mathcal{R}$. The problem of finding ${\mathcal{R}}^{e}$ is, thus, the central problem of accounting for the complex nuclear long-range correlations. The exact solution of this problem is hardly possible at both finite and zero temperatures and, thus, requires some model assumptions and simplifications which, at the same time, would retain the leading contribution of correlations beyond RPA. A truncation scheme based on retaining only two-body correlation functions leads to the following leading approximation for the long-range part of the interaction kernel:

[TABLE]

where $\Sigma^{e}(3,1)$ and $\mathcal{U}^{e}(12,34)$ are the time-dependent parts of the fermionic self-energy and induced interaction, respectively. The diagrammatic expression of $\Sigma^{e}(3,1)$ is given in Fig. 1 and the induced interaction $\mathcal{U}^{e}(12,34)$ is related to it by the dynamical consistency condition

[TABLE]

which we assume as a finite-temperature generalization of the $T=0$ one KamerdzhievTertychnyiTselyaev1997 using $\omega_{n}=2n\pi T$ with integer $n$.

Thus, the equation for the correlated propagator takes the form:

[TABLE]

At $T=0$ it can be solved in the time blocking approximation Tselyaev1989 ; KamerdzhievTertychnyiTselyaev1997 which is obtained by introducing a time projection operator into the integral part of Eq. (2). The analogous imaginary-time projection operator for the finite temperature case would look as follows:

[TABLE]

With such time projection, the triple Fourier transform onto the energy domain and integrating over two energy variables leads to the equation of the Dyson type with a single frequency at $T=0$. This allows for obtaining a clearly positive definite spectral (strength) functions and stable numerical schemes performing a sort of time ordering. However, it turns out that at $T>0$ this operator does not lead to a similar result.

We found that, in order to reach the desired form of the kernel, the imaginary-time projection operator has to be modified as follows:

[TABLE]

i.e. it needs to contain an additional multiplier with the dependence on the Fermi-Dirac occupation numbers, which turns to unity in the $T=0$ limit at the condition $\sigma_{1}=-\sigma_{2}$. Acting by the projection operator $\Theta(12,34;T)$ on the components of $\widetilde{\mathcal{R}}^{0}(12,34)$, we construct an operator $\widetilde{\mathcal{D}}$, which reads:

[TABLE]

and make a replacement $\widetilde{\mathcal{R}}^{0}\rightarrow\widetilde{\mathcal{D}}$ in the second term of Eq. (2). In Eq. (2), $\sigma_{k}=+1(-1)$ for particle (hole) components. The multiplier $\delta_{\sigma_{1},-\sigma_{2}}$ constrains the possible combinations of $(\sigma_{1},\sigma_{2})$ to be ($+1,-1$) and ($-1,+1$). A pair of states $\{12\}$ is regarded as a $ph$ ($hp$) pair if the energy difference $\varepsilon_{1}-\varepsilon_{2}$ is larger (smaller) than zero. As at $T=0$, the replacement of $\widetilde{\mathcal{R}}^{0}$ by $\widetilde{\mathcal{D}}$ corresponds to the elimination of configurations, which are more complex than $ph\otimes phonon$ ones. Thus, besides confining by only two-fermion correlation functions, in the leading approximation we keep the terms with only $ph$ and $ph\otimes phonon$ configurations and neglect higher-order terms.

After the triple Fourier transformation with respect to the imaginary-time,

[TABLE]

the summation over the fermionic discrete variables $\ell$ and $\ell^{\prime}$,

[TABLE]

and the analytical continuation to the real frequencies, we obtain

[TABLE]

where the uncorrelated propagator is given by:

[TABLE]

and the particle-vibration coupling amplitude reads:

[TABLE]

Here we defined the phonon vertex matrices $\zeta^{\mu\eta_{\mu}}$ according to:

[TABLE]

and the matrix elements of the phonon vertices $g^{\mu(\eta_{\mu})}$ Tselyaev1989 ; LitvinovaRingTselyaev2007 as:

[TABLE]

The phonon vertices $g_{12}^{\mu}$ are formally related to the transition densities $\rho_{12}^{\mu}$

[TABLE]

which can be extracted from the full response function as described in Ref. WibowoLitvinova2018 , however, in practice it is sufficient to calculate them in the random phase approximation. The static interaction, in the framework based on the covariant energy density functional (3), obeys the well-known relation: ${\hat{\widetilde{V}}}=\delta^{2}E[{\hat{\rho},\phi_{m}}]/(\delta{\hat{\rho}}\delta{\hat{\rho}})$. In the case of effective mesons, whose masses and coupling vertices are adjusted to bulk properties of finite nuclei, the double counting of the particle-vibration coupling (PVC) is removed by correcting the PVC amplitude as follows:

[TABLE]

i.e. by the subtraction of itself at $\omega=0$. This subtraction procedure and its role in the density functional based calculations is discussed in more detail in Refs. Tselyaev2013 ; LitvinovaTselyaev2007 ; LitvinovaRingTselyaev2007 for the case of $T=0$. At finite temperature we perform the subtraction in the complete analogy to the zero-temperature case.

Finally, the nuclear strength function $\tilde{S}(E)$ is defined as the response to a certain external probe associated with a one-body operator ${\hat{V}}^{0}$:

[TABLE]

The exponential factor $\left[1-\exp(-E/T)\right]^{-1}$ appears due to the detailed balance between the absorption and emission parts of the strength function WibowoLitvinova2018 ; Sommermann1983 ; RingRobledoEgidoEtAl1984 . This factor brings a new feature to the finite-temperature strength function, for instance, it changes the low-energy behavior and the zero-energy limit of $\tilde{S}(E)$. In addition to the appearance of new poles in the response function due to the thermal unblocking, the exponential factor leads to a finite strength function at $E=0$, in contrast to the spectral density $S(E)$ which has zero limit at $E\to 0$.

3 Numerical details, results and discussion

In this work, we illustrate the performance of the developed approach named finite-temperature relativistic time-blocking approximation (FT-RTBA) for the response of the even-even spherical nuclei 48Ca and 100,120,132Sn. The calculation scheme consists of the following steps. First, the closed set of the RMF equations is solved using the NL3 parametrization Lalazissis1997a of the non-linear sigma-model with the thermal fermionic occupancies (5). The obtained temperature-dependent single-particle Dirac spinors and the corresponding single-nucleon energies form the basis for subsequent calculations. Second, the finite-temperature relativistic random phase approximation (FT-RRPA) equations, which are equivalent to Eq. (2) without the time-dependent kernel $\cal W$, are solved to obtain the phonon vertices $g^{\mu}$ and their frequencies $\Omega_{\mu}$. The set of phonons, together with the RMF single-particle basis, forms the $1p1h\otimes$phonon configurations for the particle-phonon coupling amplitude ${\Phi}(\omega)$. Then, Eq. (2) is solved in the truncated configuration space which includes excitations below 25 MeV and the total response function is computed by solving Eq. (13) in either momentum or configuration space as described in LitvinovaRingTselyaev2007 ; LitvinovaRingTselyaev2008 for the $T=0$ case. Finally, the strength function is computed according to Eq. (2) with external fields of the isoscalar (L=0,2) and electromagnetic (L=1) character:

[TABLE]

The particle-hole basis was limited by $\varepsilon_{ph}\leq 100$ MeV and $\varepsilon_{\alpha h}\geq-1800$ MeV with respect to the positive-energy continuum. A direct verification with $\varepsilon_{ph}\leq 300$ MeV eliminating the dipole spurious translational mode completely showed that the physical states of the excitation spectra converge reasonably well with $\varepsilon_{ph}\leq 100$ MeV. The values of the smearing parameter $\Delta=500$ keV and $\Delta=200$ keV were adopted for the calculations of the strength functions for 48Ca and 100,120,132Sn, respectively. The phonon space was formed of the vibrations with quantum numbers of spin and parity $J^{\pi}=2^{+},\;3^{-},\;4^{+},\;5^{-},\;6^{+}$ below the energy cutoff, which amounts to 15 and 20 MeV for heavy and medium-mass nuclei, respectively. This cutoff is justified by our previous calculations at $T=0$ with the subtraction (27) which minimizes the the high-energy phonon contribution. An additional truncation was done according to the values of the reduced transition probabilities of the corresponding electromagnetic transitions. The modes with the values of the reduced transition probabilities $B(EL)$ less than 5% of the maximal one (for each $J^{\pi}$) were neglected. We have put a condition that the same truncation criteria on the phonon energy, $J^{\pi}$ and the reduced transition probability should be kept for all temperature regimes in order to make a fair comparison of the calculated strength distributions. At zero and low temperatures varying these limits would almost not change the resulting strength distributions and, therefore, we typically adopt these cutoff values for $T=0$ calculations. At high temperatures, as a consequence of the significant thermal unblocking, we see the appearance of many additional phonon modes. This effect is illustrated below in Figs. 4-6 for $J^{\pi}=0^{+},\;1^{-},\;2^{+}$ and it is similar for all $J^{\pi}$ values. On average, at $T\sim 5-6$ MeV the number of phonon modes in the truncated model space becomes by an order of magnitude larger than at $T=0$. We note that at high temperatures the saturation of the results with respect to the $B(EL)$ cutoff is rather slow, but the 5% cutoff, which is at the limit of our current computational capabilities, is still very reasonable. We have checked that by varying the cutoff value between 10% and 5%. Another truncation was made on the absolute values of the numerator of Eq. (22), so that the contributions with $|n(\varepsilon_{2},T)-n(\varepsilon_{1},T)|\leq 0.01$ were excluded from the BSE. FT-RTBA calculations at high temperatures are otherwise very prohibitive, while at low and moderate temperatures this truncation does not affect the results. We found, however, that in some cases at high temperature values this truncation can cause a noticeable contribution of the spurious translational mode in the lowest-energy parts of the dipole spectra, therefore, the truncation on the occupation numbers was released as much as possible to remove this contribution where it appeared.

An illustration of the thermal RMF calculations for the compound nuclei 48Ca and 132Sn at $J^{\pi}=0^{+}$ is given in Fig. 2. As it follows from Eqs. (3,5), the effect of finite temperature on the total (binding) energy of a thermally excited nucleus is mainly induced by the change of the fermionic occupation numbers. Their direct effect is caused by promoting the fermions to higher-energy orbits, which obviously increases the total energy. Another, indirect, contribution is related to the change of the nucleonic densities which play the role of sources for the meson and photon fields WibowoLitvinova2018 . The changes in the bosonic fields, in turn, translate to the nucleons via the self-consistent set of the thermal RMF equations. Thus, the thermodynamical equilibrium is achieved throughout the iteration procedure. As one can observe in Fig. 2, the dependence $E^{\ast}(T)$ is nearly parabolic, as it is expected for the non-interacting Fermi gas, except for the lowest temperature values. The flat behavior of the excitation energy at these temperatures is a result of the presence of the large shell gaps right above the Fermi surface in the doubly-magic nuclei: small temperatures are insufficient to promote the nucleons over the shell gaps. At $T\geq$ 1 MeV the dependence $E^{\ast}(T)$ is very well approximated by the parabolic fits, which are also shown in Fig. 2. The corresponding coefficient at the quadratic term is found to be in a very good agreement with the empirical Fermi gas level density parameter $a=A/k$, where $8\leq k\leq 12$ Harakeh2001 .

Fig. 3 explains qualitatively the new mechanism of formation of nuclear excitations, which emerges with the temperature growth. Due to the thermal occupation factors of the single-nucleon orbits (5), a sharp borderline between particle and hole states becomes diffuse as all states acquire a non-zero probability to be occupied. Therefore, the notion of particle-hole pairs adopts a conditional character. As far as the response is concerned, in the simplest approximation it is governed by the uncorrelated propagator of Eq. (22), which has formally the same ansatz as at $T=0$. At zero temperature, when the occupation factors take the values of zero above the Fermi energy and unity below the Fermi energy, the only non-vanishing contributions come, obviously, from the pairs of single-nucleon states located on different sides from the Fermi surface (particle-hole, or $ph$ pairs). At finite temperature, due to the fractional occupancies, the uncorrelated propagator of Eq. (22) acquires contributions from all pairs of states, except for those of the identical states. This fact allows for an extended notion of the particle-hole pairs. As mentioned in the previous section, we regard a pair of states $\{12\}$ as a ’thermal’ $ph$ ($hp$) pair if the energy difference $\varepsilon_{1}-\varepsilon_{2}$ is larger (smaller) than zero. Remarkably, in this context a single-particle state alone does not carry a particle or a hole character, but adopts it only being paired with another state. In Fig. 3 we have denoted as $ph$ the pairs which have a particle-hole character already at $T=0$ and as $\widetilde{ph}$ the new, or ’thermal’ $ph$-pairs which are classified as particle-particle or hole-hole ones at $T=0$ and do not contribute to the zero-temperature uncorrelated propagator and, thus, to the zero-temperature response.

The isoscalar monopole, electromagnetic dipole, and isoscalar quadrupole strength distributions were computed for the doubly-magic 48Ca in a wide range of temperatures, according to Eqs. (2 - 2,13) with the value of the imaginary part of the energy variable (smearing parameter) $\Delta=500$ keV. The results of the FT-RTBA calculations are displayed in the right panels of Fig. 4 and can be compared to the FT-RRPA results given in the left panels. At $T=0$ RTBA provides a very reasonable description of the strength distributions in 48Ca, in agreement with data, see, for instance, a direct comparison for the dipole response in our recent work EgorovaLitvinova2016 . This can serve as a good benchmark of the theory and, together with the thermal mean-field calculations presented above, justify the reliability of our approach to the strength distributions at finite temperature. The isoscalar monopole response (ISMR) shown in the top panels of Fig. 4 represents the so-called breathing mode, or the compression mode, and serves as one of the most important nuclear characteristics. It is directly related to the nuclear compressibility which is studied very extensively, by both theory and experiment. One of the main motivations for these studies is that the characteristics of response of finite nuclei to compression can be translated or extrapolated to the compressibility of the uniform nuclear matter. Its incompressibility parameter, $K_{\infty}$, is one of the most important parameters of the nuclear matter equation of state (EOS) and, thus, plays the key role in modeling astrophysical objects, from core collapse supernovae to neutron star mergers. The incompressibility of a finite nucleus with the mass number $A$, $K_{A}$, can be determined from the centroid of ISMR and related to the $K_{\infty}$ through the well-known leptodermous expansion. The possibility of such extraction of $K_{\infty}$ uniquely is, however, very much debated as it leads to different results if different nuclei are used for this procedure. In particular, nuclear superfluidity seems to play a role in nuclear compressibility, because the results of analyses based on the monopole responses of closed-shell and open-shell nuclei vary, see, for instance, a recent review Garg2018 . As far as the astrophysical aspect of nuclear compressibility is concerned, its temperature dependence could be, in principle, another factor to be taken into account. One can see, comparing the two top panels of Fig. 4, that PVC effects are not strongly manifested in the case of monopole response at $T=0$. This is a well-known phenomenon related to a significant cancellations of terms with the self-energy and induced interaction in Eq. (2) in the monopole channel LitvinovaRingTselyaev2007 . At $T=1$ MeV both FT-RRPA and FT-RTBA strength distributions demonstrate the appearance of a new soft mode at $E\approx$ 9 MeV. This is already a clear manifestation of the thermal unblocking mechanism discussed above. At $T=3$ MeV we observe first drastic changes in the monopole strength distribution: the soft mode at 9 MeV strengthens considerably while a new prominent mode emerges at $E\approx$ 5 MeV. The locations of these states reveal the presence of $0\hbar\omega$ excitations which are barely possible at $T=0$. Another new effect at this temperature is enforcing the fragmentation of the main high-energy peak whose origin stems from $2\hbar\omega$ excitations. The latter two features – the emergence of strong low-lying states and the reinforcement of fragmentation of the high-energy peak - become even more amplified at $T=5$ MeV, moreover, the low-energy strength begins to dominate. In addition, the entire strength distribution shifts toward lower energies. Such changes in the isoscalar monopole strength distribution obviously affect its centroid, that, in turn, means a change of the nuclear compressibility. We do not go into more details of the nuclear compressibility here, however, point out that its sensitivity to the temperature increase might be studied elsewhere more extensively for the temperature regimes which are relevant for astrophysical processes.

It can be seen from the middle panels of Fig. 4 that the main thermal effects which are observed for the case of ISMR, are also visible in the dipole channel. At zero and low temperatures, the latter is dominated by the giant dipole resonance (GDR), which is formed by collective oscillations of protons and neutrons against each other and centered at about $E_{GDR}=18A^{-1/3}+25A^{-1/6}$ MeV Bortignon1998 . In the context of its temperature dependence, the dipole response was studied most extensively in experiments with heavy ion collisions Gaardhoje1984a ; GaardhojeEllegaardHerskindEtAl1986 ; Bracco1989 ; Ramakrishnan1996 ; Mattiuzzi1997 ; Heckman2003 ; Santonocito2006 ; Wieland2006 ; Bracco2007 . The general conclusions from these studies are (i) a fast growth of the GDR’s width with temperature and (ii) disappearance of the high-frequency collective oscillation at high temperatures.

As discussed in detail in Ref. Bortignon1998 , the latter phenomenon occurs largely due to the particle emission which is not taken into account consistently in our approach and should be included with the exact continuum as it is done, for instance, in Ref. LitvinovaBelov2013 . However, even in the present calculations we observe a reduction of the coherence of the GDR with temperature. Indeed, one can see that, in the case of 48Ca, at $T=5$ MeV the soft mode starts to dominate over the high-frequency mode while the latter becomes more fragmented. At higher temperatures these effects become even more pronounced, see, for instance, an example of our calculations in Refs. LitvinovaWibowo2018 ; WibowoLitvinova2018 .

Another enhancement of the low-energy strength comes with the exponential factor in Eq. (2) which is analyzed below. Experimental and theoretical studies of the nuclear dipole response also have far reaching broader impacts. The electric dipole polarizability, which is proportional to the $m_{-1}$ moment of the dipole strength distribution, provides another very important aspect of the nuclear EOS, namely the slope of its symmetry energy, see, for instance, a recent review Roca-Maza2018 . As the $m_{-1}$ moment of the dipole strength distribution is particularly sensitive to the low-energy fraction of the strength, it should be also sensitive to the temperature increase, which can become crucial for the symmetry-energy sector of the EOS at high-temperature regimes.

The specific feature of a typical isoscalar quadrupole response (ISQR) at $T=0$ is the presence of two main collective structures: one at low (few MeV) and one at high energy. The properties of the ISQR are also of a broader interest as it is sensitive to the nuclear matter incompressibility, provides a constraint on the effective mass, and contributes to determining the nuclear symmetry energy and the neutron skin thickness Roca-Maza2018 . In the context of theoretical models for nuclear response beyond the random phase approximation, the low-lying quadrupole mode is of a particular interest, because, together with its octupole counterpart, it provides the most important contribution to the PVC mechanism. The low-energy collective quadrupole and octupole phonons couple most strongly to the single-particle degrees of freedom and produce the largest fraction of the fragmentation of both single-particle and collective nuclear states. The ISQR in 48Ca at various temperatures is displayed in the bottom panels of Fig. 4, where FT-RTBA results can be compared to those of FT-RRPA. The first thermal unblocking effect can be seen already at $T=1$ MeV as the formation of the new low-energy mode. With the temperature increase, i.e. at $0\leq T\leq 3$ MeV the fragmentation of the high-energy quadrupole peak enforces gradually while at $T=5$ MeV it becomes visibly stronger. The latter is common for all $0^{+},1^{-}$ and $2^{+}$ strength distributions and related to the tremendously enhanced low-energy strength of all multipolarities $J^{\pi}=2^{+},\;3^{-},\;4^{+},\;5^{-},\;6^{+}$ included into the PVC amplitude (2). These new strong low-lying modes tend to also couple strongly to the single-particle degrees of freedom and to contribute significantly to the $\Phi(\omega)$ amplitude, while the quadrupole contribution is one of the dominant ones.

Here it should be noted that in Fig. 4 we plot not the final strength functions, but their zero-temperature analogs $S(E)$ of Eq. (2) called spectral densities. The strength function at finite temperature ${\tilde{S}}(E)$ has to obey the detailed balance between absorption and emission and, thus, to be corrected by the exponential factor $\left[1-\exp(-E/T)\right]^{-1}$ which is singular at $E=0$. In order to illustrate the role of this new factor, Fig. 5 shows the FT-RTBA dipole strength distributions in 48Ca with and without this factor. For a better assessment of its role at low temperatures, the right panel displays the same strength functions as shown in the left panel, but on the logarithmic scale. One can notice that the high-frequency giant dipole resonance remains almost unaffected at all temperatures. The situation is, however, different for the low-energy part. It can be seen that already at $T=1$ MeV the strength function acquires a finite value at $E=0$ because in the present case $S(0)=0$ (no physical poles appear at zero energy) and the influence of the exponential factor on the low-energy strength grows with temperature. At $T=5$ MeV, for instance, the low-energy peak receives an enhancement of a factor of two. Thus, although such high temperatures are not relevant for typical astrophysical situations, temperatures of $\sim 1-2$ MeV can be reached in many cases, so that the exponential factor should not be neglected in the finite-temperature astrophysical modeling of quantities which are sensitive to the zero-energy limit of the strength functions.

In this work we mostly focus on the spectral density $S(E)$ because it has a well-defined zero-energy limit. Self-consistent calculations for the dipole strength should demonstrate, in particular, the absence of the admixture of the spurious translational mode. The presence of the singular factor in the full strength function ${\tilde{S}}(E)$ which makes it finite at $E=0$ conceals the low-energy behavior of the distribution $S(E)$ and makes an assessment of the quality of the numerical implementation difficult. We note that the exponential factor gives a very small contribution at temperatures below 4 MeV and can be easily included for a comparison with experimental data when they are available for the low-lying strength. By these reasons we found showing $S(E)$ more informative than the full strength function and illustrated the influence of the exponential factor in Fig. 5 as an example of its generic behavior. For a reliable prediction of the low-energy strength distributions at finite temperatures at least one more aspect, namely the (single-particle) continuum should be included into the description, because it can further modify the trend of the low-energy behavior of the strength distributions of all multipoles LitvinovaBelov2013 .

The picture of the temperature evolution of the nuclear response would not be complete without considering heavier nuclear systems. Therefore, we have performed calculations for some isotopes of the tin chain. The nucleus 120Sn is of special interest because of the availability of experimental data for its giant dipole resonance in terms of the temperature-dependent width Fultz1969 ; Heckman2003 ; Ramakrishnan1996 . We have also considered the neutron-rich doubly magic 132Sn because of its particular importance for the astrophysical r-process nucleosynthesis, so that the precise knowledge about its strength functions is required for the r-process modeling MumpowerSurmanMcLaughlinEtAl2016 . The neutron-deficient 100Sn was chosen in order to find out whether the temperature evolution of the nuclear response is different in the neutron-rich and neutron-deficient cases. The results of our calculations for the dipole response of these three nuclei with the smearing parameter $\Delta=200$ keV are displayed in Fig. 6. One can observe immediately that the doubly-magic 100,132Sn isotopes reveal a trend in the temperature evolution of the dipole strength, which is similar to one already seen in 48Ca. The situation with the open-shell 120Sn is somewhat different in the low-temperature regime. This nucleus is superfluid at $T=0$, and the superfluidity effects typically produce a more spread strength distribution because of the presence of more poles in the response function. However, above the critical temperature of $T_{c}\approx$ 0.66 MeV the superfluidity, in its traditional understanding within the Bardeen-Cooper-Schrieffer or Hartree-(Fock)-Bogoliubov models, vanishes, while at $T=1$ MeV the thermal unblocking is not yet well pronounced. As a result, the dipole strength in 120Sn at $T=1$ MeV is less broad than at $T=0$, in contrast to the case of the closed-shell nuclei. Otherwise, at $T\geq$ 1 MeV the strength distributions in all three nuclei evolve in a similar manner. As in the case of 48Ca, we observe a gradual increase of the fragmentation and broadening of the high-energy peak, due to both Landau damping and PVC, and the formation of the enhanced strength at low energies. Comparing the neutron-balanced 120Sn, the neutron-rich 132Sn and the neutron-deficient 100Sn nuclei, one can see that the proton-neutron composition does not play the major role in the temperature evolution of the dipole strength. In particular, the formation of the new pronounced low-energy structure at high temperature is very similar for the three considered tin isotopes. This means that this structure is not associated with the pygmy resonance as a neutron skin oscillation, but has a purely thermal unblocking origin. Indeed, from the Fig. 3 it is clear that in the calculations within a sufficiently large single-particle basis the number of $\widetilde{ph}$ pairs above the Fermi surface can grow tremendously with the temperature increase and form a coherent low-lying mode, because many of such pairs can have close values of the energy differences which enter the denominator of Eq. (22). Accounting for the single-particle continuum is expected to provide a more accurate description of the low-energy strength and will be addressed in future work. Another interesting aspect of this new low-energy structure is a relatively weak fragmentation due to the PVC mechanism, which is obtained for the strength functions of all multipolarities. As it was discussed in Refs. LitvinovaWibowo2018 ; WibowoLitvinova2018 , this can possibly occur due to missing contributions of the PVC-induced ground state correlations KamerdzhievTertychnyiTselyaev1997 , which should be considered among other future directions.

The width of the strength distribution and the energy weighted sum rule (EWSR) are the most important integral characteristics of the GDR which are often addressed in theoretical and experimental studies. For example, they help constraining theoretical approaches because of their almost model-independent character. In theoretical calculations, verification of the sum rules usually helps in testing the consistency of numerical implementations.

In order to illustrate the present case, the evolution of GDR’s width $\Gamma(T)$ with temperature obtained in FT-RTBA for 120Sn and 132Sn nuclei is shown together with experimental data in the left panel of Fig. 7. Experimental data which are available only for 120Sn are given as well. The theoretical values for the widths at $T=0$ are provided by our earlier calculations LitvinovaRingTselyaev2008 ; LitvinovaRingTselyaev2007 , respectively. As follows from Fig. 2, the thermal unblocking effects hardly appear at $T\leq 1$ MeV in 132Sn because of its specific shell structure, namely, the presence of the large shell gap in the vicinity of the Fermi energy in both proton and neutron subsystems. In particular, for protons, which form the $Z=50$ closed shell and have the next available orbitals only in the next major shell, temperatures below 1 MeV are not sufficient to promote particles over the shell gap with a sizable occupancy. For neutrons in 132Sn as well as for protons in 120Sn the situation is similar. The lowest orbit available for neutrons in 120Sn is the intruder $1h_{11/2}$ state where particles get promoted relatively easily, but the next shell gap occurs right after this orbit. As a result, at $T=1$ MeV there are still not many options for the $\widetilde{ph}$ pair formation and for the thermal unblocking. Together with the disappearance of superfluidity above the critical temperature $T_{c}\approx 0.66$ MeV, this explains, for instance, the unexpectedly small GDR’s width at $T=1$ MeV in 120Sn reported in Ref. Heckman2003 . After $T=1$ MeV in 132Sn and $T=2$ MeV in 120Sn we observe a quick growth of $\Gamma(T)$ due to the formation of the low-energy shoulder of the dipole strength distribution and due to the reinforced fragmentation of the high-energy peak emerging from the finite-temperature effects in the PVC amplitude $\Phi(\omega)$. In general, at all temperatures $T\geq$ 1 MeV the low-energy shoulder of 132Sn is stronger, which leads to a larger overall width in 132Sn, compared to 120Sn. This occurs because of the more neutron-rich character of 132Sn. Notice here that, since our standard Lorentzian fit of the microscopic strength distribution fails in recognizing the distribution as a single peak structure at high temperatures, the GDR’s widths for T$>$3 MeV in 132Sn and for T$>$4 MeV in 120Sn are not presented in Fig. 7.

At $T=0,T=1$ MeV and $T=3$ MeV the agreement on the GDR’s width in 120Sn between the FT-RTBA calculations and data is very reasonable while at $T=2$ MeV FT-RTBA underestimates the experimental value of the width by a factor of two. Here we should note that effects of thermal shape fluctuations are not included into the present approach. In order to do this in a microscopic way, we would need to generalize our method to the deformed case and to consider a superposition of shapes and orientations, as discussed, for instance, in Ref. Bortignon1998 . As a very large number of shapes and orientations would be required, the authors of Refs. OrmandBortignonBrogliaEtAl1990 ; OrmandBortignonBroglia1996 , for instance, did that in an effective way using phenomenological relations for the GDR. However, as it is shown, for instance, in Refs. Bonche.1984 ; Lisboa.2016 ; Zhang.2017 , nuclear deformation tends to disappear at some critical temperature. These studies agree on the upper bound of the critical temperature $T_{c}\approx 2.0-4.0$ MeV. This means that above those temperatures medium-mass and heavy nuclei reveal the tendency to take the spherical shape. On the other hand, one can see from Refs. OrmandBortignonBrogliaEtAl1990 ; OrmandBortignonBroglia1996 that, although the theory of thermal shape fluctuations describes the GDR’s width at temperatures $1.5\leq T\leq 3$ MeV well, it does not reproduce $\Gamma(T)$ in the low temperature regime at $T<1.5$ MeV in 120Sn. Our microscopic theory, however, shows a reasonable agreement with data for $T=0$ and $T=1$ MeV. From these arguments we can conclude that the thermal shape fluctuations in medium-heavy nuclei should be important in a limited temperature range between $T\approx 1$ MeV and $T\approx 3$ MeV. They could be, perhaps, included on top of our model in the way proposed in Refs. OrmandBortignonBrogliaEtAl1990 ; OrmandBortignonBroglia1996 , however, we leave this beyond the scope of our paper admitting that thermal shape fluctuations have to be included for a correct comparison to experimental data. Apart from that, our results for $\Gamma(T)$ in the entire range of temperatures under study show a nearly quadratic dependence, which is in agreement with the Fermi liquid theory Landau1957 . The FT-RTBA results are also consistent with those of the microscopic approach of Ref. Bortignon1986 , which are available for the GDR energy region at $T\leq$ 3 MeV.

The energy-weighted sum rule for 48Ca and 132Sn nuclei as a function of temperature within FT-RRPA and FT-RTBA is given in the right panel of Fig. 7, in the percentage with respect to the Thomas-Reiche-Kuhn (TRK) sum rule. The EWSR at $T>0$ can be calculated in analogy with the zero-temperature case Sommermann1983 ; Barranco1985a . In our framework the meson-exchange interaction is velocity-dependent, so that already in RRPA and RQRPA at $T=0$ we observe up to 40% enhancement of the TRK sum rule within the energy regions which are typically studied in experiments LitvinovaRingTselyaev2007 ; LitvinovaRingTselyaev2008 , in agreement with data. In the resonant time blocking approximation without the ground state correlations associated with the PVC mechanism the EWSR is known to have exactly the same value as in RPA Tselyaev2007 . A small violation comes only with the subtraction procedure Tselyaev2007 ; LitvinovaTselyaev2007 . At $T=0$ in the RTBA corrected by the subtraction we find a few percent less EWSR in finite energy intervals below 25-30 MeV than in RRPA, however, this difference becomes smaller for larger intervals. As in RTBA the strength distributions are more spread, one finds more strength outside finite intervals. The situation is very similar at T$>$0. In Fig. 7 (b) one can see that the EWSR stays nearly flat with a slow decrease with the temperature growth due to the fact that the entire strength distribution shifts toward lower energies. The EWSR values from FT-RRPA and FT-RTBA become nearly equal at $T=6$ MeV in both nuclei, because apparently at this temperature the high-energy tails of the strength distributions become less important.

As it follows from Eq. (2), calculations of the nuclear response with small smearing parameter $\Delta$ are capable of resolving individual states in the low-energy region. Such calculations can reveal fine details of the evolution of the thermally emergent strength, as it was illustrated in Ref. WibowoLitvinova2018 . In addition, such calculations allow for extracting the transition densities from the response function, which provide a valuable information about the underlying structure of individual states. Fig. 8 gives an illustration of such kind with the radial dependencies of the neutron and proton transition densities in the neutron-rich 48Ca (top panels) and in the neutron-deficient 100Sn (bottom panels). The proton and neutron transition densities are shown for the most prominent peaks below 10 MeV, which are selected separately for each temperature value. In 48Ca proton and neutron transition densities show in-phase oscillations inside the nucleus while neutron oscillations become weakly dominant outside at $T=0$. At $T=1$ MeV the major oscillations shift toward the nuclear surface while in the outer area one observes the typical pattern of the neutron skin oscillation known as pygmy dipole resonance. At $T=3$ MeV the general picture changes very little, but at $T=5$ MeV the oscillations of proton and neutron subsystems in the outer region start to develop the pattern of the out-of-phase oscillation. The latter is more typical for the giant dipole resonance and reproduced in numerous theoretical calculations at $T=0$. Our previous analysis of the thermal transition densities in a heavier nuclear system, such as 68Ni, has shown that the low-energy peak at high temperature has indeed some features of collective nature WibowoLitvinova2018 . In order to clarify whether this new feature is related to a considerable neutron excess, we have performed a similar analysis for the neutron-deficient 100Sn nucleus. As it can be seen from the bottom panels of Fig. 8, at $T=0$, it also exhibits the in-phase oscillations of protons and neutrons inside the nucleus, but with a little dominance of proton oscillations in the surface area. Only little changes are observed at $T=1$ MeV, and at $T=3$ MeV one can see a clear proton dominance in the outer area. As in the case of 48Ca, at $T=5$ MeV one starts to distinguish a GDR-like pattern of the out-of-phase oscillation in the low-lying state at $E=$3.92 MeV while the in-phase oscillations are still present inside the nucleus and on the surface. Notice that at $3\leq T\leq 5$ MeV the oscillations extend to far distances from the nuclear central region, which indicates that these oscillations are mainly associated with nucleons in the continuum. The latter is again consistent with the picture of Fig. 3 which illustrates the possibility to create a large amount of $\widetilde{ph}$ transitions between the states located above $\varepsilon_{F}$.

In the spectra of compound nuclei at finite temperature, especially at the high temperature values, the thermal unblocking allows for the appearance of the states with arbitrary small energy values. Such states can be indeed observed in Figs. 4 - 6 and, in principle, they would form a smooth continuum if the calculations were done in the complete basis LitvinovaBelov2013 . We have analyzed the underlying structure of the lowest dipole peaks at $E=$1.44 MeV in 48Ca and at $E=$0.76 MeV in 100Sn at $T=6$ MeV that is displayed in Fig. 9. As at these transition energies the PVC effects included in the present model are rather weak, in this figure we show only the FT-RRPA calculations for the proton and neutron transition densities. In both neutron-rich 48Ca and neutron-deficient 100Sn we have obtained very similar pictures for the behavior of the proton and neutron transition densities in the lowest dipole states: they oscillate in phase and have their major peaks inside the nucleus while protons in 48Ca and neutrons in 100Sn show a pronounced minor peak in the outer area. The lowest dipole state at $E=$1.44 MeV in 48Ca is composed predominantly by the proton $\widetilde{ph}$ configuration $3p_{3/2}-2d_{5/2}$ and the analogous state at $E=$0.76 MeV in 100Sn is dominated by the neutron $\widetilde{ph}$-configuration $1i_{13/2}-2h_{11/2}$. Other states in the lowest-energy domain are similar and their single $\widetilde{ph}$-configuration character is in agreement with their relatively low transition probabilities.

4 Summary

In this article, which is dedicated to the memory and inspiring ideas of Pier Francesco Bortignon, we discuss a finite- temperature extension of the nuclear response theory beyond the relativistic RPA. We review the formalism of the thermal relativistic mean field and the time blocking method generalized for the finite-temperature case as the leading-order approach to the time-dependent part of the in-medium nucleon-nucleon interaction. Technically, the temperature-dependent imaginary-time projection operator is used to reduce the Bethe-Salpeter equation to an equation with one energy (frequency) variable in the energy domain.

The finite-temperature relativistic time blocking approximation was implemented on the base of the meson-nucleon Lagrangian of quantum hadrodynamics with the NL3 parametrization. We investigated the temperature dependence of the monopole, dipole and quadrupole response in the closed-shell medium-light 48Ca and of the dipole response of the medium-heavy 100,120,132Sn isotopes. It was found that the temperature dependence of the nuclear response is very generic and exhibits similar features for all multipoles and all nuclei under study. The most remarkable effects of the temperature increase on the spectra are (i) the reinforcement of both the Landau damping and the particle-vibration coupling, (ii) the formation of the low-energy strength due to the thermal unblocking, and (iii) the shift of the entire strength distribution toward lower energies. We have investigated numerically the integral characteristics of the obtained dipole spectra, such as the energy-weighted sum rule and the width of the giant dipole resonance, in finite energy intervals up to high-temperature regimes. The EWSR was found to be very robust up to very high temperatures and the width exhibits a nearly parabolic growth with temperature, in accordance with the Landau theory. The obtained temperature dependence of the GDR’s width is consistent with the available data, except for the temperature regime where thermal shape fluctuations become important. The results obtained for the dipole and quadrupole response functions are also consistent with the NFT calculations of Ref. Bortignon1986 .

We have discussed some broader impacts of the presented developments, from the astrophysical r-process nucleosynthesis to the nuclear matter equation of state. The proposed treatment of the many-body correlations at finite temperature is of a general character and can be widely applied to the response of strongly-correlated systems other than finite nuclei. Future developments may consider the inclusion of continuum effects and ground state correlations associated with the particle-vibration coupling.

Acknowledgements

The authors greatly appreciate discussions with Peter Schuck and Jian Li. This work is supported by the NSF Career Grant PHY-1654379.