Coordinate-space solver for finite-temperature Hartree-Fock-Bogoliubov calculations using the shifted Krylov method

Background: In order to study structure of protoneutron stars and those in subsequent cooling stages, it is of great interest to calculate inhomogeneous hot and cold nuclear matter in a variety of phases. The finite-temperature Hartree-Fock-Bogoliubov (FT-HFB) theory is a primary choice for this purpose; however, its numerical calculation for superfluid (superconducting) many-fermion systems in three dimensions requires enormous computational costs.

Purpose: To study a variety of phases in the crust of hot and cold neutron stars, we propose an efficient method to perform the FT-HFB calculation with the three-dimensional (3D) coordinate-space representation.

Methods: Recently, an efficient method based on the contour integral of Green's function with the shifted conjugate-orthogonal conjugate-gradient method was proposed [Phys. Rev. C 95, 044302 (2017)]. We extend the method to finite temperature, using the shifted conjugate-orthogonal conjugate-residual method.

Results: We benchmark the 3D coordinate-space solver of the FT-HFB calculation for hot isolated nuclei and fcc phase in the inner crust of neutron stars at finite temperature. The computational performance of the present method is demonstrated. Different critical temperatures of the quadrupole and the octupole deformations are confirmed for ${}^{146}\mathrm{Ba}$. The robustness of the shape coexistence feature in ${}^{184}\mathrm{Hg}$ is examined. For the neutron-star crust, deformed neutron-rich Se nuclei embedded in the sea of superfluid low-density neutrons appear in the fcc phase at a nucleon density of 0.045 ${\mathrm{fm}}^{-3}$ and a temperature of $k_{B}T=200$ keV.

Conclusions: The efficiency of the developed solver is demonstrated for nuclei and inhomogeneous nuclear matter at finite temperature. It may provide a standard tool for nuclear physics, especially for the structure of hot and cold neutron-star matter.

Figure 1

Schematic illustration of the contour $C$ and poles of the Green's function $G\left(z\right)$ (closed circles) and of the Fermi-Dirac function $f_T\left(z\right)$ (open circles) in the complex plane.

Figure 2

Schematic illustration of the contour $C^{\prime }$ and poles of the Green's function $G\left(z\right)$ (closed circles).

Figure 3

The convergence behavior of the shifted-COCG and shifted-COCR methods for solutions of Eq. (28). Three typical points, $z(\theta =0)$ (COCG/COCR: purple/green) and $z(\pi /2)$ (orange/black) on the contour of Eq. (40) with $E_{\mathrm{cut}}=100$ MeV and $h=16\pi k_{B}T$, in addition to the reference point $z=0$ (red/blue), have been taken as examples. The norms of the residual vectors $|{\mathit{r}}_k^{\sigma }/{\mathit{r}}_0^{\sigma }|$ are shown as functions of the iteration number. This is a case of the FT-HFB calculation for the center of mass of a ${}^{146}\mathrm{Ba}$ nucleus (${\mathit{r}}^{\prime }={\mathit{r}}_c$) with a temperature of $k_{B}T=200$ keV. The mesh and box sizes are the same as those in Sec. 3c.

Figure 4

The convergence behavior of (a) the neutron spin-up density $\rho ({\mathit{r}}^{\prime }\uparrow ,{\mathit{r}}^{\prime }\uparrow )$ and (b) the neutron pair density $\nu \left({\mathit{r}}^{\prime }\right)=\kappa ({\mathit{r}}^{\prime }\uparrow ,{\mathit{r}}^{\prime }\downarrow )$ at the center of mass of a ${}^{146}\mathrm{Ba}$ nucleus (${\mathit{r}}^{\prime }={\mathit{r}}_c$) with $k_{B}T=200$ keV. The shifted-COCG method is shown by red lines and the shifted-COCR method by black lines. See text for details.

Figure 5

(a) Calculated neutron average paring gap, (b) quadrupole and octupole deformation parameters, and (c) specific heat as functions of temperature for ${}^{146}\mathrm{Ba}$. In panel (b), the quadrupole deformation of the dripped uniform neutrons is shown by the solid line.

Figure 6

Nucleon density profiles in the $z-x$ plane for ${}^{146}\mathrm{Ba}$ at different temperatures: (a) $T=0$, (b) $k_{B}T=0.8$ MeV, (c) $k_{B}T=1.2$ MeV, and (d) $k_{B}T=1.6$ MeV.

Figure 7

Potential energy surface, calculated with the constrained FT-HFB method, as a function of the quadrupole deformation, for ${}^{184}\mathrm{Hg}$. Panel (a) is the total energy $E$, while panel (b) is the free energy $F$. See text for details.

Figure 8

Nucleonic density distribution in the inner crust of neutron stars at beta equilibrium with the neutron chemical potential ${\mu }_n=10$ MeV and the proton chemical potential ${\mu }_p=-58.5$ MeV: (a) 3D plot of baryon density, (b) 2D contour plot of neutron density in the $z=0$ plane, and (c) the same as (b) but for protons.

Figure 9

Calculated neutron gap ${\mathrm{\Delta }}_n\left(\mathit{r}\right)$ on the $x$ axis. The protons are clustered near the both edges.