Internuclear potentials from the Sky3D code

The Sky3d code is a published (v1.0 [1] and 1.1 [2]) time-dependent nuclear Hartree–Fock code which is widely used to study nuclear dynamics for process such as fusion [3], fission [4, 5], neutron star matter [6], and giant resonances [7–12].

In studies of fusion, the Frozen Hartree–Fock (FHF) approximation [13–17] is sometimes used as a method to obtain details of the ion-ion potential of the fusing reactants. It is a microscopic method working at the level of wave functions and effective nuclear interactions. As such, one may learn about role of the microscopic internucleon interactions giving rise to internuclear potentials using FHF. Such input may come from, for example, the neutron skin as predicted by a particular effective interaction, which will directly affect the internuclear potential as calculated by FHF.

Here we present instructions to make minimal adjustments to the distributed Sky3d code v1.1 [2] in order to calculate FHF potentials, and document how to run the resulting code and process the output. In the next section we give a brief description of the FHF method, then follow with a section containing the necessary instructions for adapting the code, followed by instructions for running the code to extract a FHF potential.

Frozen Hartree–Fock's name derives from the Hartree–Fock link to the density functional (EDF) formalism in which the starting quantity is the EDF which gives the energy of a nucleus as a function of the densities involved [18] as

$$\begin{eqnarray}&&E=\int d\vec{r}{ \mathcal H }[\rho (\vec{r})].\end{eqnarray} \tag{ 1 }$$

Here, the Hamiltonian density ${ \mathcal H }$ is given as a functional of the particle density $\rho (\vec{r})$. This is a shorthand for all relevant densities and currents. The Sky3d code uses the Skyrme interaction [19, 20] and full details of the link between the interaction and the functional can be see in the original derivations [21], further developments [22, 23], and the documentation of the Sky3d code [2].

The FHF potential comes from the placement of two nuclei made up from their ground-state densities, at a chosen relative position $\vec{R}$. The interaction potential between them is then given as the EDF energy of the combined system, from which is subtracted the energies of the two participating nuclear ground states of the nuclei in isolation:

$$\begin{eqnarray}&&{V}_{\mathrm{FHF}}(\vec{R})=\int d\vec{r}{ \mathcal H }[{\rho }_{1}(\vec{r})+{\rho }_{2}(\vec{r}-\vec{R})]-E[{\rho }_{1}]-E[{\rho }_{2}],\end{eqnarray} \tag{ 2 }$$

where ρ₁ and ρ₂ refer to the ground state densities of the two interacting nuclei.

Note that there is no accounting for the Pauli principle when adding the densities together. As the relative vector between the nuclei, $\vec{R}$, becomes small, the FHF approximation becomes worse because of this. More sophisticated approaches are available in which the deficiency is removed [24, 25]. FHF remains a useful first-order microscopic approach to assess barrier positions, systematic effects of nuclear interactions, and as a reference for more sophisticated methods.

In principle, any static Hartree–Fock solver should be able to evaluate equation (2). In practice, the ability to place multiple nuclei at an arbitrary position relative to each other is not available in most codes, but is present in codes in which the static HF calculation is a preamble to time-dependent calculations of collisions, for which the initial placement of the nuclei at time t = 0 is exactly what is needed for a FHF calculation.

In this section the prescription for amending the Sky3d v1.1 code is given. It is assumed that the reader has been able to obtain the code from the repository given in the published article concerning the code [2], successfully compiled and run it.

The minimal change necessary to the code to enable a FHF calculation is to turn off the re-orthogonalization that is performed after the placement of multiple nuclei on the coordinate space grid as the initial condition to a time-dependent calculation. The call to the relevant routine occurs in line 173 of main3d.f90. This line should be commented out, along with the IF-block in which it sits. The code block in main3d.f90 should read, with lines of context:

Here, the only changes are the four exclamation marks at the beginning of the lines 172–175 removing the statements that lead to calling the schmid routine. With this change, a recompiled code is able to calculate the FHF potential. The amended code may still be used in static mode to calculate the ground states successfully.

A Frozen Hartree–Fock potential between two nuclei is calculated by running the Sky3d code as per its documentation to obtain ground states for the desired nuclei, and then constructing an input file as if for a time-dependent calculation, but with the number of time-steps set to zero. The placement of the nuclei in the input file gives the separation for the FHF calculation. One must also ensure that no initial boost is given to the nuclei, as would usually be the case for setting up a collision. As an example, we assume that the sample static calculation (for ¹⁶O) has been run from the Test directory in the distributed code, to give the file O16 containing the ground state. Here we present a sample for005 input file which calculates the $\int d\vec{r}{ \mathcal H }[{\rho }_{1}(\vec{r})+{\rho }_{2}(\vec{r}-\vec{R})]$ term in (2) for two ¹⁶O nuclei for a single value or $\vec{R}$:

The input file consists of a series of Fortran namelists, with values as follows: wffile should be set to NONE to suppress the writing of the (large) output file containing the full set of wave functions, which is not needed for the FHF calculation. In the &force namelist, the same force name and pairing choice should be used as in the static calculation. Likewise Coulomb should be switched on with tcoul if it was used in the static calculation. imode must be 2 to indicate a time-dependent calculation, so that the placement of multiple nuclei on the grid is allowed. nof=2 indicates that two nuclei are to be initialised.

The &grid namelist defines the coordinate space grid. The grid spacing dx should match that of the static calculation. The number of grid points in each direction nx, ny, and nz should be such that the entire static calculation can be placed offset from the origin by a sufficient amount to provide the desired range in $\vec{R}$. The variable periodic should be set to F so that that two nuclei are isolated and not in a lattice of repeated charges.

nt=0 is set to ensure that no actual time-steps are solved in the TDHF equations. dt should be given some value, but it is inconsequential since no time-steps take place. In &fragments, filename gives the name of the two files containing the ground state wave functions. In this case, we show a repeated file name since we are dealing with the potential for two identical nuclei; ¹⁶O + ¹⁶O. fix_boost is set to T so that initial velocity vectors for each nucleus are specified (rather than the alternative setup of giving a centre of mass energy and impact parameter). These initial velocity boost vectors default to zero and so do not need to be otherwise specified. Finally, the two fcent vectors give the centres of the two nuclei on the (x,y,z) grid. The values should be chosen to give the desired value of $\vec{R}$. Note that the values should be such that the centres fit exactly on defined grid points. Though the code will interpolate if needed, the interpolation errors are significant at the level of the size of the FHF potential.

If the amended Sky3d code is run with this input file, the output file energies.res will contain the required energy in the fifth column. In order to obtain ${V}_{\mathrm{FHF}}(\vec{R})$ as in equation (2) the ground state energies of each fragment must be subtracted. This is found in the file conver.res from the static run(s), in the second column. By repeatedly running the FHF calculation for different fcent values, a full curve for ${V}_{\mathrm{FHF}}(\vec{R})$ may be obtained. We present the result for ¹⁶O+¹⁶O with the Skyrme force SV-bas [26] in figure 1. For repeated calculations to give a full set of points at different $\vec{R}$, simple scripts should be used to drive the code and generate a set of input files, as well as processing the outputs.

Zoom In Zoom Out Reset image size

We have presented a simple modification to a published nuclear time-dependent Hartree–Fock code (Sky3d) which allows the calculation of nucleus-nucleus potentials using the Frozen Hartree–Fock approximation—an application not originally built-in to the running of the code.

The authors gratefully acknowledge support from the UK Science and Technology Facilities Council (STFC) via awards ST/P005314/1 and ST/L005743/1 as well as computing time on the STFC DiRAC computer.