Effective Coulomb interaction in actinides from linear response approach

Computational Materials Science 171 (2020) 109270

Contents lists available at ScienceDirect

Computational Materials Science journal homepage: www.elsevier.com/locate/commatsci

Effective Coulomb interaction in actinides from linear response approach a

b,⁎

Ruizhi Qiu , Bingyun Ao a b

, Li Huang

a,⁎

T

Science and Technology on Surface Physics and Chemistry Laboratory, Mianyang 621908, Sichuan, PR China Institute of Materials, China Academy of Engineering Physics, Mianyang 621907, Sichuan, PR China

A R T I C LE I N FO

A B S T R A C T

Keywords: Actinides Effective Coulomb interaction Linear response DFT + U

The effective on-site Coulomb interaction (Hubbard U) between 5f electrons in actinide metals (Th-Cf) is calculated with the framework of density-functional theory (DFT) using linear response approach. The value of Hubbard U seldom relies on the exchange-correlation functional, spin-orbital coupling, and magnetic state, but depends on the neighbor distance and actinide element. In particular, a pronounced dependence of Hubbard U on the Wyckoff position is shown for β -U. We find an abrupt increase of Hubbard U from α -Pu to α -Am along the actinide series, which is similar to that of atomic volume. Within local-density approximation or generalized gradient approximation of Perdew-Burke-Ernzerhof revised for solids, the calculated structural properties from DFT + U with U from linear response calculation are in good agreement with the experimental data.

1. Introduction Actinide elements produce a variety of fascinating physical behaviors, such as three charge-density wave phases of U metal [1], the unique phonon dispersion of α -U [2,3] and δ -Pu [4–6], unconventional heavy-fermion superconductivity in PuCoGa5 [7], hidden order phase in URu2Si2 [8], and high-rank multipolar order in NpO2 [9]. Along the actinide series, an atomic volume anomaly occurs near Pu, as shown in Fig. 1. As the atomic number increases, the equilibrium atomic volume first decreases parabolically, then increases abruptly about 43% from Pu to Am, followed by an almost constant evolution [10]. For the actinide metals, the 7s and 6d electronic bands are broad and s/d electrons are strongly delocalized, which do not change much along the series. With the increasing number of 5f electrons, the parabolic-like behavior of atomic volume for the light actinides (Th-Pu) is indicative of a system with itinerant 5f electrons that are strongly bonding and delocalized. This is similar to that of 5d transition metal series, in which the atomic volume first decreases due to filling of the 5d bonding states and then increases owing to the filling of the antibonding states [10]. While for the transplutonium elements (Am-Cf), the little volume change with increasing f electrons implies that the 5f electrons are localized. This is similar to the 4f rare-earth metals, in which the volume changes little along the 4f series [10]. To describe the volume anomaly and the delocalized-localized change, the effective on-site Coulomb interaction U between 5f electrons is extremely important. One may expect the increase of Coulomb parameter U with increasing f-electron count and an abrupt change of U

⁎

in the vicinity of Pu. Experimentally, U was defined as the energy cost of moving one 5f electron from one atom to another and the energy difference among different electron configurations is derived from the atomic spectral data [14]. The estimated U’s from the reactions 2( f n d 2s ) → f n − 1 d3s + f n + 1 ds [15] and 2( f n ds ) → f n − 1 d 2s + f n + 1 s [16] are listed in Table 1. The dependence of U on the valence is also considered [17] and actinides can exhibit different U. The values of trivalent actinide metals, as presented in Table 1, approximately equal to those in Ref. [16] except Pu and Am. The atomic limit was used above and the case of real interest is the metallic situation. By truncating the free-atom wave function at the Wigner-Seitz radius of the actinide metals and performing the Hartree-Fock calculations, Herbst et al. [18] used the energy difference between free-atom and Wigner-Seitz cell calculation to correct the U values. The corrected values of U, as shown in Table 1, are much larger than those from the atomic limit. The U of Cm is extraordinary large, which was attributed to the effect of Hund’s rule [18]: the stability of f 7 configuration. The Hartree-Fock atomic calculations were also directly performed for the Coulomb interaction strengths of actinide elements [19] and the values of U0 are given in Table 1. All the above data show a general increase in U from Th to Cm but no abrupt change of U between Pu and Am. Nowadays it is possible to compute the Coulomb parameter U from first-principles and the methods include constrained random-phase approximation (cRPA) [23,24], constrained density-functional theory (DFT) [25,26], and linear response approach [27]. In the cRPA, the effective interaction U is obtained as the expectation value of Coulomb interaction on the wave functions of the localized basis set. As this

Corresponding authors. E-mail addresses: [email protected] (B. Ao), [email protected] (L. Huang).

https://doi.org/10.1016/j.commatsci.2019.109270 Received 9 July 2019; Received in revised form 5 September 2019; Accepted 5 September 2019 0927-0256/ © 2019 Elsevier B.V. All rights reserved.

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

actinide except uranium [21]. In the work we apply the linear response approach to systematically calculate the Hubbard U of actinide metals. Different allotropic phases under ambient pressure, as shown in Fig. 1, are considered. Three approximations of Exc functionals including local-density approximation (LDA) [37,38], GGA of Perdew-Burke-Ernzerhof (PBE) [33] and PBE revised for solids (PBEsol) [39], are taken into account. For the calculated U of α -phase as a function of atomic number, there is indeed a jump near Pu. This trend of U along the actinide series will improve our understanding of the 5f electron behavior. Since the structural properties of Pu and the atomic volumes of actinides (U-Cm) were well described by GGA + U with U from cRPA calculation [22], it is straightforward here to test the validity of DFT + ULR with ULR from linear response calculation in the description of the structural properties of actinide metals. The rest of this paper is organized as follows. Computational method and details are given in Section 2. The computational models, i.e., various phases of the actinide metals, are introduced in the A. Section 3 presents the results of Coulomb parameters U and discussion. The assessment of DFT + ULR on the description of the structural properties is shown in Section 4. The last section summarizes the main achievement of this work.

Fig. 1. Atomic volume of each actinide metal (blue circle and line) calculated from the crystal structure of room temperature phase. The square block represents the crystallographic volume of the other solid-state allotropic phases of actinides at high temperature. All the data are obtained from Pearson’s Handbook [11] except β -Pa [12] and β -Cf [13].

2. Computational method

Table 1 Coulomb parameters U of actinides from the literature. Metal

Th

Pa

U

Np

Pu

Am

Cm

Ref. Ref. Ref. Ref. Ref. Ref. Ref. Ref.

1.7 1.5 1.53 3.0 2.72 2.3 – –

2.0 1.6 1.67 3.5 2.71 – – –

2.6 2.3 2.30 3.7 2.89 – 1.87–2.1 0.8

2.7 2.6 2.64 3.7 3.07 – – 1.0

3.2 3.5(± 1) 2.91 4.0 3.24 – – 0.95

3.0 5(± 1) 2.72 4.2 3.34 – – 1.5

– –

[15] [16] [17] [18] [19] [20] [21] [22]

2.1. DFT + U formalism In order to account explicitly for the on-site Coulomb interaction, a corrective functional inspired to the Hubbard model was introduced into the DFT + U formalism [26,40,41]. The localized orbitals in which “Hubbard U” will operate are introduced by defining the corresponding occupation number matrix (nI ) as the projection of Kohn-Sham orbitals σ (ψkν ) into the states of localized basis set of choice (ϕmI , e.g., 5f atomic orbitals in this work):

9.5 3.57 – – 3.4

Iσ nmm ′ =

∑

f kσν 〈ϕmI ∣ψkσν 〉〈ψkσν ∣ϕmI ′〉.

kν

approach allows for the calculations of the effective interaction matrix U and its energy dependence, it becomes particularly popular within the DFT + DMFT (dynamical mean field theory [28]) community. Recently Amadon [22,29] use the self-consistent cRPA to estimate the U of actinide metals, which are also listed in Table 1. The U values for U-Am are quite small, which is contrary to ordinary belief. Specifically, Hubbard U of Pu and Am are much smaller than the typical width of 5f bands (2.0 eV [10]), which indicates that the electronic correlations of Pu and Am are weak, being in contradiction to the photoemission experiments [30,31]. In addition, the effect of exchange-correlation (Exc) functional and allotropic phase on U was not evaluated. Nevertheless, using these small value of U, the structural properties of Pu and the volume jump in actinides could be well described by DFT + DMFT [28] and DFT + U [32] within the generalized gradient approximation (GGA) [33] of Exc functional [22,29,34]. For constrained DFT, the Hubbard U is obtained from the total-energy variation with respect to the occupation number of the localized orbitals, which is identified as the shift in the Kohn-Sham eigenvalues by virtue of Janak theorem [35]. The implementation of DFT using the localized basis sets makes it possible to change the occupation of the localized orbitals but fails to perform the screening of other delocalized orbitals. A further improvement of constrained DFT was the linear response approach [27], in which one can utilize the pseudopotential methods using plane wave basis sets and not just the localized basis sets. The localized orbitals are perturbed by a single-particle potential and the Hubbard U is determined by using the density response functions of the system with respect to these localized perturbations. The scheme is internally self-consistent and widely used in the first-principles calculation community [36]. To our best knowledge, computation of U using linear response approach has not been performed to pure

(1)

Here σ , k , ν , I , and m are spin, Bloch wavevector, band index, Hubbard σ atom index, and index of localized orbitals, respectively; f kν is the Fermi–Dirac distribution of the Kohn-Sham states. Note that the trace of σ occupation number matrix, nI = Tr[nI ] = ∑σm nmm , is the total occupation number of localized orbital in Hubbard atom I. If the localized orbitals are set as the f atomic orbitals, the index m might be understood as the magnetic quantum number, m = −3, −2, …, 3. In this work the DFT + U corrective energy functional is chosen as the simplified rotationally invariant form [41,27], i.e.,

EU [{nI }] =

∑ I

UI Tr[nI (1 − nI )]. 2

(2)

Note that the trace also includes the summation of spin σ . This simplified version is equivalent to the fully localized limit of fully rotationally invariant form [40] with the exchange parameter being zero. 2.2. Linear response approach From the corrective functional (2), it is realized that the Hubbard U I could be computed as the second derivative of the ground state total energy with respect to the occupation number nI . The linear response approach aimed at computing this second derivative [27]. First of all, the single-body potential is perturbed by an external potential that only acts on the localized orbitals of a Hubbard atom I, α I ∣ϕmI 〉〈ϕmI ∣, in which α I is the amplitude of the perturbation on Hubbard atom I. Solving the modified Kohn–Sham equations with the perturbed single-body potential yields an α I -dependent ground state total energy

E [{α I }] = min {EDFT [ρ (r )] + α I nI }, ρ (r )

2

(3)

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

7 × 7 × 7, 7 × 7 × 9, 5 × 7 × 3, 3 × 3 × 3, and 9 × 5 × 3 M-P k point meshes are chosen, respectively. Spin–orbit coupling (SOC) is considered in some linear response calculations for the discussion of the effect of SOC on the Hubbard U. For the structural properties, all the calculations are performed by including SOC. The magnetic orders considered here include diamagnetic (DM), ferromagnetic (FM) and non-collinear antiferromagnetic (AFM) and will be taken into account case by case. For the calculation of Hubbard U, the response function matrix χ and χ0 are derived by setting α parameter to ± 0.1, ± 0.2, ± 0.3, ± 0.4 , and ± 0.5 eV. The parameters are so chosen to remove the numerical noise. Only NWyckoff column of χ and χ0 is extracted from the calculation with NWyckoff being the number of Wyckoff positions in the cell. All the other matrix elements are reconstructed by symmetry. Considering the periodicity of DFT calculation in solid, the local perturbation should not be overlapped and thus a supercell approach is adopted here. For fcc, bcc, bct phase, α -U, β -Np, we used a 2 × 2 × 2 supercell. For dhcp phase and γ -Pu, a 2 × 2 × 1 supercell is used. A 1 × 2 × 2 supercell is constructed for α -Np and the unit cell is used for β -U, α -Pu and β -Pu. The U difference between the calculation using the above supercell and a larger supercell are within 0.05 eV. For example, for the calculated U of fcc Th within LDA, the value from a supercell of 64 atoms is 0.045 eV smaller than that from a supercell of 32 atoms. For bcc U, the U value from a supercell of 27 atoms is 0.02 eV smaller than that from a supercell of 16 atoms. For α -U, the U value from a supercell of 64 atoms is only 0.002 eV larger than that from a supercell of 32 atoms. Thus the result of linear response calculation will be presented with a precision of 0.1 eV.

with ρ (r ) being the electron density and EDFT [·] being the total energy function of pure DFT. Then Legendre-Fenchel transformation is used to transform this optimization problem into the corresponding dual problem. That is to say, an nI -dependent ground state total energy could be recovered,

E [{nI }]=E [{α I }] − α I

(

∂E [{α I }] ∂α I

)

=E [{α I }] − α I nI .

(4)

Based on this definition, the first derivative of the total energy with respect to nI is given by

∂E [{nI }] = −α I , ∂nI

(5)

and the second derivative,

∂2E [{nI }] ∂α I = − I. I 2 ∂ (n ) ∂n

(6)

{nJ }

In actual numerical calculation, all vary in response to the change of α I and the evaluated quality should be the response function χ whose matrix element is evaluated from the finite differences,

χI , J =

δnI . δα J

(7)

In addition, the electronic wave function of non-interacting electron systems would be reorganized in response to the perturbation and this response should be eliminated since it is not related to the electronelectron interaction. This response function χ0 could be evaluated by performing a non-self-consistent electronic structure calculation with the charge density being kept constant and collecting the response of this non-interacting system in terms of variation of all {nJ } . Finally the Hubbard U I is given by

U I = (χ0−1 − χ −1 )I , I .

3. Effective Coulomb interaction For actinide metals at various solid-state phases, the calculated values of effective Coulomb interaction U within different computational formulation are listed in Table 2. Before examining the trend of U along the series, it is instructive to consider its dependence on Exc functional, SOC and magnetic order. For most phases, the difference of U among different approximations of Exc functional (LDA, PBE, and PBEsol) is negligible, being smaller than the precision of 0.1 eV. The largest difference lies in Cf, whose pseudopotential was constructed recently and deserved to be carefully examined. The calculation related to Cf may be unreliable and is shown here for completeness. Even for the most complex low-symmetry metals, β -U, α -Pu, and β -Pu, the results from LDA/PBE/PBEsol are close. This consistent behavior of U illustrates the reliability of linear response approach and our calculation. We also perform linear response calculation of α -Th, α -Pa, α -U, and δ -Pu within LDA + SOC and find that the effect of SOC is also negligible. This could be understood by noting that the Hamiltonian of SOC is a single body term to employ the scalar relativistic correction. In the presence of SOC, the angular and spin part of single-electron wave function for a free atom will be altered but the radial part, which mainly contributes to the effective on-site Coulomb interaction, is unchanged. Most of the lowest-energy magnetic states are DM and thus the comparison among different magnetic states could not be made. For some magnetic systems from standard DFT calculation, the calculated U from AFM, FM, and DM states are compared. The difference of U among AFM, FM, and DM is smaller than 0.2 eV for most phases. The differences mainly lie in the transplutonium metals, in which the presence of the magnetic moment increase the localization of 5f electrons. As a whole, the computational formulation has limited effect on the value of Hubbard U. Note that the value of U is different for the inequivalent atom and the difference is pronounced for the structures with several Wyckoff positions such as β -U, α -Pu, and β -Pu. The most striking is β -U, in which the value of Hubbard U ranges from 1.3 to 2.9 eV. To understand this wide range of Hubbard U, it is important to note the coordination environment of uranium atoms. In Fig. 2 we plot the radial distribution

(8)

2.3. Computational details All calculations have been performed using the Vienna Ab-initio Simulation Package (VASP) [42]. The ion-electron interaction was described using the projector augmented wave (PAW) formalism [43,44], which has the accuracy of all electron methods because it defines an explicit transformation between the all-electron and pseudopotential wave functions by means of additional partial-wave basis functions [45]. We used the official PAW pseudopotentials for Th, Pa, U, Np, Pu, Am, Cm, and Cf, which were generated with the following reference valence configurations, 6s 2 6 p6 7s 2 5 f 16d1, 6s 2 6 p6 7s 2 5 f 2 6d1, 6s 2 6 p6 7s 2 5 f 3 6d1, 6s 2 6 p6 7s 2 5 f 4 6d1, 6s 2 6 p6 7s 2 5 f 56d1, 6s 2 6 p6 7s 2 5 f 6 6d1, 6s 2 6 p6 7s 2 5 f 7 6d1, and 6s 2 6 p6 7s 2 5 f 8 6d 2 , respectively. The energy levels of 5f and 6d are very close to each other and their occupations will be redistributed in the actual numerical calculation. The pseudopotential of Bk is not constructed and thus Bk is not included in our calculation. The linear response calculations of Hubbard U and DFT + U calculations are performed on the solid phases of actinide metals (Th-Cm, Cf) at ambient pressure, which are introduced in the A. Unless otherwise stated, the Hubbard U’s are calculated using the structures with the experimental lattice parameters, which are obtained from Pearson’s handbook [11] except β -Pa [12] and β -Cf [13]. Here we compared three approximations of Exc functionals: LDA, GGA of PBE and PBEsol. The partial occupancies for each Kohn-Sham wave function are determined by the tetrahedron method with Blöchl correction. For the convergence tests, the cutoff energy for plane wave basis and the Monkhorst-Pack (M-P) [46] k point meshes are determined. In this work, all the cutoff energy is chosen as 450 eV except Cf, which is set as 600 eV. A 9 × 9 × 9 M-P k-point meshes are used for fcc and bcc phases. For dhcp-phase, bct-Pa, α -U, β -U, α -Np, β -Np, α -Pu, β -Pu, and γ -Pu, the 9 × 9 × 3, 9 × 9 × 11, 11 × 7 × 7, 3 × 3 × 7, 3

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

Table 2 Calculated Hubbard U (in units of eV) of actinides from different computational schemes. Note that NWyckoff values are presented for each phase with NWyckoff being the number of Wyckoff positions. The prime in the symbol of Wyckoff position is added to distinguish from the other Wyckoff positions and has no special meaning. Here the reduced coordinates of 8i′, 2e′, and 4i′, are (0.3655, 0.0391, 0.0), (0.869, 0.25, 0.894), and (0.613, 0.0, 0.759), respectively [11]. Struct.

Mag.

Pos.

LDA

PBE

PBEsol

LDA + SOC

DM DM DM DM DM DM DMa

γ -U α -Np

bcc Pnma

DM DM

β -Np

P4212

DM

γ -Np

bcc

α -Pu

P 21/ m

AFM FM DM DMa

β -Pu

C 2/ m

FM

γ -Pu δ -Pu

P 63/ m fcc

∊-Pu α -Am

bcc dhcp

DM AFM FM DM DM AFM

4a 2a 2a 2a 4a 4c 2a 4f 8i 8i′ 8j 2a 4c 4c 2a 2c 2a 2a 2a 2e 2e 2e 2e 2e 2e 2e 2e′ 2a 4h 4i′ 4i 4i 8j 8j 8a 4a 4a 4a 2a 2a 2c 2a 2c 4a 4a 2a 2c 2a 2c 4a 4a 2a 2c 2a 2c 4a 4a

2.6 2.6 2.6 2.7 3.0 2.4 2.9 2.3 1.3 2.5 1.4 2.5 2.2 2.2 2.4 2.4 2.5 2.5 2.4 2.1 2.1 2.1 2.3 2.3 2.4 2.2 2.5 2.5 2.4 2.7 2.4 2.6 2.9 2.3 2.7 2.8 2.7 2.8 2.6 3.2 3.2 3.3 3.3 3.3 3.4 3.6 3.6 3.3 3.5 3.6 3.3

2.5 2.5 2.5 2.6 2.9 2.4 2.9 2.3 1.3 2.5 1.4 2.5 2.2 2.3 2.4 2.4 2.5 2.5 2.4 2.0 2.1 1.7 2.2 2.1 2.2 2.2 2.6 2.5 2.5 2.7 2.6 2.6 2.3 2.0 2.6 2.8 2.7 2.7 2.6 3.2 3.2 3.3 3.3 3.2 3.3 3.4 3.4 3.3 3.4 3.3 3.2 4.3 4.7 4.8 4.4 4.2 4.4

2.5 2.5 2.5 2.6 3.0 2.4 2.9 2.4 1.3 2.5 1.4 2.5 2.2 2.2 2.4 2.4 2.4 2.5 2.5 2.1 2.1 1.7 2.2 2.0 2.3 2.1 2.6 2.5 2.6 2.8 2.5 2.7 2.5 2.1 2.7 2.8 2.7 2.6 2.6 3.2 3.2 3.2 3.3 3.2 3.3 3.4 3.4 3.3 3.4 3.3 3.4 4.5 4.7 4.1 3.9 4.3 4.3

2.6

α -U β -U

fcc bcc bct bcc fcc Cmcm P 4 2/ mnm

α -Th β -Th α -Pa β -Pa

FM

β -Am

fcc

α -Cm

dhcp

AFM FM AFM FM

β -Cm

fcc

α -Cf

dhcp

AFM FM FM AFM

β -Cf

fcc

FM AFM

2.6

2.4

Fig. 2. The radial distribution function of α -, β -, and γ -U with respect to the atom at different Wyckoff positions.

that of Hubbard U. The subtle difference between U(β -U 8i) and U(β -U 8j) could be explained by virtue of the difference of the second-nearestneighbor distance. Apparently, the neighbor distance of actinide atoms is a dominant factor that determines the value of U. This also applies to α -Pu and β -Pu, in which the difference of U and neighbor distance is not that prominent. For α -Pu, the atoms with the largest U in the Wyckoff position 2e′ have the longest nearest-neighbor distance. For β -Pu, the atoms with the shortest nearest-neighbor distance in the Wyckoff position 4i′ have the smallest U. For the allotropes of actinides, U is also closely related to the neighbor distance. Usually, the neighbor distance of high-temperature phase is larger than that of room-temperature α -phase due to the thermal expansion. Accordingly, the calculated U of high-temperature phases is larger than that of α -phase. The larger the neighbor distance is, the larger the effective Coulomb interaction is. For Th, the nearestneighbor distance of α -Th (3.60 Å) is approximately equal to that of β -Th (3.56 Å), which implies the negligible difference of U between the two phases. For Pa, the nearest-neighbor distance of β -Pa (3.30 Å) is a bit larger than that of α -Pa (3.22 Å), indicating a relatively larger U of β -Pa with respect to α -Pa, as shown in Table 2. For U, the relationship between neighbor distance and U has been elaborated above. For Np, the U difference for the inequivalent atoms can be neglected and the order of U is U(α -Np) < U(β -Np) < U(γ -Np), which is same as that of nearest-neighbor distance. For Pu, the effect of inequivalent atoms on U cannot be neglected but here we also take the average value. The order of effective interaction, U(α -Pu) < U(β -Pu) < U(∊-Pu) < U(γ -Pu) < U(δ -Pu), is same as that of crystallographic volume which is shown in Fig. 1 and can be reckoned as the average neighbor distance. For Am and Cm, the approximately equivalent nearest-neighbor distances between α - and β -phases yield the approximately equivalent values of U. To clearly illustrate the evolution of effective interaction along the series, we plot the values of U within PBE and atomic volumes for α -phases of actinide in Fig. 3. As we expected, there is a clear jump of U near Pu, which is analogous to the jump of atomic volume. The 42% increase of U from α -Pu to α -Am is an indication of electronic delocalized-localized change because the on-site Coulomb interaction could be understood as the one-center Coulomb integral between the localized wave functions in the Hubbard model and the increasing localization means the increasing effective interaction. This behavior of U could in turn explain the evolution of atomic volume along the actinide series. For the lighter actinides, the effective on-site Coulomb interaction among 5f electrons is small, yielding delocalized 5f electrons which participate the bond and the parabolic-like behavior of atomic volume.

2.7

a The calculated magnetic moment in these cases is smaller than 0.3 μB per atom. For other cases the magnetic moment is negligible.

function (RDF) of β -U for five atoms at different Wyckoff positions in Fig. 2. The RDFs of α -U and γ -U are also presented for comparison. The nearest-neighbor distance d ranges from 2.09 to 3.11 Å, which is responsible for the variation of Hubbard U. In particular, the order of first neighbor distance is d(β -U 8i) = d(β -U 8j) < d(β -U 4f) < d(α -U 4c) < d(β -U 8i′) < d(γ -U 2a) < d(β -U 2a), which is exactly same as 4

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

4. Structural properties The structural properties of actinide metals were already well described by standard DFT within GGA [57–65,56,66–69,71]. For example, most of the characteristics of the unique Pu phase diagram could be well reproduced by GGA [61]. But GGA also has some disadvantages such as the inclusion of artificial magnetic moment [70] to increase the localization of 5f electrons [10]. Considering the order of the typical bandwidth of f bands in actinide metals (2 eV [10]), the values of previously evaluated U’s in Table 1 and our calculated U’s in Table 2 equal to or are greater than the bandwidth, implying that the actinides are electronic systems with strong or intermediate correlation. Thus a consistent theoretical framework of actinides should include the explicit description of the effect Coulomb interaction such as DMFT [47,5,48–50,22], exact diagonalization [72], and Gutzwiller approximation [73]. For δ -Pu, the phonon spectra predicted by DMFT [5] were shown to be surprisingly accurate by later inelastic X-ray scattering measurement [6] and the predicted photoemission spectra [47,49,50] are also in good agreement with the experimental results [74]. Because of the computational cost, the above theoretical calculations cannot be easily performed for a large system [34]. Thus a static treatment of effective interaction, DFT + U, is widely used in the description of the structural and phonon properties of actinide metals [51–54,21,20,55,34]. These DFT + U calculations often rely on empirical values of U [18] except the above-mentioned theoretical U [21,20,34]. Since our calculated U’s are different from those used in the literature, it is necessary to assess the validity of DFT + ULR (ULR from linear response calculation) in the description of structural properties of actinide metals. The theoretical lattice parameters from DFT and DFT + ULR will be compared with the experiments [11,12] and the previous theoretical values [57–64,54,65,56,66–69,34]. To obtain the optimal lattice parameters of actinide metals, we first perform the structural relaxation using the conjugate-gradient and quasi-Newton algorithm. As mentioned above, SOC is considered in all calculations and the relaxations keep the experimental structural symmetry. Within the framework of DFT + ULR , the lattice should be relaxed using structurally consistent Hubbard U [75] since ULR is dependent on the lattice spacing. However, the difference of ULR between the optimized and experimental structure is very small and then the deviation of lattice parameters using constant ULR from that using structurally consistent U is negligible. Thus for simplification, the computational scheme is chosen as DFT + ULR with constant ULR from Table 2. The values of ULR we used for α -Np, β -Np, and dhcp structure are set as the average of calculated U’s for different Wyckoff positions. For β -U, α -Pu, and β -Pu, it is not appropriate to choose the average value of U but after many tests, we are not able to obtain the converged results for a calculation with different values of U for the various Wyckoff positions. Thus their results are not given in this work. The presented optimal lattice parameters and bulk moduli are obtained by fitting the energy-volume data with the third-order BirchMurnaghan equation of state (EOS) [76]. This procedure has to keep the continuity of the energy-volume curve, which partially rules out the presence of metastable states in the magnetic and DFT + ULR calculation. The problem of metastable states could be fixed by controlling the occupation number matrix [54,77]. For light actinides and Am, the constraint of zero magnetic moment, i.e., DM, is imposed. It should be emphasized that for Th and Pa, the lowest-energy magnetic states are DM but for U, Np, Pu, and Am, the lowest-energy magnetic states from DFT + ULR are not DM. For δ -Pu and Am, the magnetic order of the lowest-energy states is AFM and the corresponding lattice parameters are also listed. For Cm, the lattice parameters of AFM α -phase and FM β -phase are shown. The calculated lattice parameters and bulk moduli from LDA/PBE/ PBEsol and LDA/PBE/PBEsol + ULR are shown in Table 3. These results are in reasonable agreement with the previous theoretical calculation. For LDA, the lattice parameters are underestimated due to the

Fig. 3. The calculated Hubbard U using PBE (left) and atomic volume (right) of α -phase for actinides. The two background colours reflect the delocalized-localized change of 5f electrons along the series.

For the transplutonium, the interaction is large, giving rise to localized 5f electrons and the almost constant behavior in atomic volume. Note that the slight variation of U for lighter actinides is attributed to the parabolical behavior of atomic volume. To clarify this, we perform the linear response calculation of U for fcc An (An = Pa, U, Np, Pu) with the same lattice parameter 5.0 Å and find that all the values are about 2.9 eV. Finally let us compare our calculated U’s with the Coulomb parameters of actinides in the literature, as shown in Table 1. For the elements between Th and Cm, the values of our calculated U’s range between 2.0 and 3.4 eV, which is a compromise between those estimate using the atomic limit and the Hartree-Fock atomic calculation. The difference suggests that it may be not optimal to use empirical U directly in DFT + U calculations. The empirical U is still useful and can provide guidelines for the calculated U from first-principles. For fcc Th, the value U in this work is close to that in Ref. [20] which is calculated using constrained DFT [25]. The small difference may result from the absence of the screening of other delocalized s and d orbitals in the constrained DFT. For U, our PBE results are consistent to that in Ref. [21] which also use the linear response approach. The small difference between our PBE results and that from Ref. [21] is due to the different computational parameters. For β -U, the effect of different Wyckoff positions on the calculation of U may be not considered in Ref. [21]. Our calculations differ strongly from the U’s from self-consistent cRPA calculations [22] for Th-Am but are close to that for Cm. We suspect that the metallic feature of actinide influenced the localization of the constructed Wannier orbitals and then underestimated the on-site Coulomb repulsion. In addition, the self-consistent calculation of U may be affected by the choice of Exc. Nevertheless, using these small value of U, the structural properties of Pu and the volume jump in actinides could be described by DFT + DMFT [28] and DFT + U [32] within the generalized gradient approximation (GGA) [33] of Exc functional [22,29,34]. For δ -Pu and transplutoniums, the commonly used value of U in the DFT + DMFT [47,5,48–50] and DFT + U [51–55] calculations is 4.0 eV, which is greater than our calculated U’s. Those calculations could provide a different electronic structure by capturing the correlated nature of 5f electrons [56] but for the structural properties, standard GGA could well reproduce the experimental data [57–65,56,66–69]. Take δ -Pu as an example, one can only use an AFM order to obtain equilibrium properties in good agreement with experiments within GGA and LDA + U (U = 4.0 eV) [54]. The problem of including an artificial magnetic moment [70] is ubiquitous [10,34].

5

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

Table 3 Calculated structural parameters and bulk moduli from different computational schemes. The experimental values of the lattice parameters and bulk moduli are obtained from Refs. [11,12] and Refs. [78–84], respectively. All our calculations take into account SOC. LDA

LDA + ULR

PBE

PBE + ULR

PBEsol

PBEsol + ULR

Exp.

a

4.8883 4.946a

4.9958

5.1147

4.9410

5.0343

5.0861

B

65.1 78.6a

70.7

62.2

62.6

67.5

58

a a

3.9014 3.8323

3.9646 3.9389

4.0613 4.0200

3.9364 3.8608

3.9933 3.9587

4.11 3.932

c

3.815g 3.1119

3.1920

3.2632

3.1434

3.2192

3.238

B

3.104g 115.6

115.4

101.4

110.8

112.2

157

a a

117.9g 3.5678 2.7187

3.6716 2.8399

3.7513 3.1470

3.5998 2.7470

3.6948 2.8670

3.81 2.8537

b

5.7541

5.8051

5.9279

5.7615

5.8229

5.8695

c

4.8079

4.9633

5.3661

4.8398

5.0035

4.9548

y

0.0973

0.1020

0.1121

0.0986

0.1036

0.1024

B

186.4

140.6

53.6

171.4

130.7

135.5

γ -U

a

3.3488

3.4892

3.7088

3.3708

3.5101

3.5320

α -Np

a

6.4848

6.6079

6.8689

6.5058

6.6460

6.663

b

4.5756

4.6814

4.8624

4.5966

4.7089

4.723

c

4.7148

4.8562

5.0688

4.7407

4.8907

4.887

x1

0.0404

0.0373

0.0329

0.0395

0.0357

0.036

z1

0.2131

0.2141

0.2109

0.2133

0.2132

0.208

x2

0.3197

0.3158

0.3099

0.3186

0.3144

0.319

z2

0.8534

0.8517

0.8342

0.8519

0.8479

0.842

B

247.8

151.0

53.0

232.9

132.1

118

a c z a a b c B

4.5648 3.3672 0.6777 3.2291 3.3807 5.5422 7.1726 270.0

4.7228 3.4996 0.6559 3.4506 3.2178 5.5706 9.3972 56.2

5.0327 5.056a 5.046b 5.055c 5.024d 5.036e 5.045f 58.9 73.1a 56d 58e,f 4.0058 3.9256 3.942d 3.927e 3.907g 3.1967 3.193d 3.212e 3.185g 99.2 92.0d 92.2e 98.2g 3.6658 2.8094 2.845h 2.797i 2.829d 5.8475 5.818h 5.867i 5.770d 4.9123 4.996h 4.893i 4.950d 0.0992 0.1025h 0.098i 0.10d 140.3 133.0h 147.0i 138.0d 3.4307 3.46h 3.43i 6.6037 6.588d 4.6659 4.706d 4.8231 4.800d 0.0410 0.04d 0.2154 0.22d 0.3173 0.32d 0.8553 0.86d 193.0 150d 4.6504 3.4406 0.6739 3.3054 3.4508 5.6597 7.3107 201.9 33k

4.9802 3.4673 0.6400 3.7250 3.2699 5.9390 10.8379 31.1

4.5898 3.3707 0.6725 3.2465 3.3803 5.5619 7.2392 251.4

4.7746 3.5273 0.6428 3.5064 3.2171 5.6596 9.8000 44.4

4.897 3.388 0.625 3.52 3.1587 5.7682 10.162 25.7

α -Th

β -Th α -Pa

β -Pa α -U

β -Np

γ -Np γ -Pu

(continued on next page) 6

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

Table 3 (continued) LDA

LDA + ULR

PBE

PBEsol

PBEsol + ULR

Exp.

4.1527

4.5940

4.6347

121.3

47.7

29–30

4.3612

4.8079

4.6347

4.604l 34.4

73.3

39.8

29–30

3.8001

3.1883

3.5802

3.6357

3.530l 3.5075 11.3433 39.4 3.6114 12.4303 17.8 4.9484 5.1989

2.9644 9.4425 107.7 3.3140 10.2282 42.9 4.1688 4.5699

3.3980 10.9804 45.4 3.4436 11.6729 25.6 4.7789 4.8867

3.463 11.231 29.7 3.463 11.231 29.7 4.894 4.894

3.5629 11.6281 39.2 5.0654

3.3710 10.6627 45.0 4.7140

3.4653 11.2842 45.1 4.9161

3.502 11.32 36.5 4.93

PBE + ULR l

δ -Pu

a B

4.1159 4.032j 141.6

4.5551 53.8

4.2473 4.146a 4.290j

4.643j

BAFM

72.5

46.9

∊-Pu

a

3.1707

3.5190

α -Am

a c B aAFM cAFM BAFM a aAFM

2.9275 9.3769 127.3 3.2545 9.5820 54.4 4.1221 4.4546

3.3453 10.7985 50.0 3.3768 11.2143 34.2 4.7099 4.7720

aAFM cAFM BAFM aFM

3.2772 10.3450 48.8 4.5603

3.4067 11.0565 50.8 4.8130

aAFM

β -Am

α -Cm

β -Cm

4.7021

33 4.7956

4.2770 4.142j 94.5

36.4 37l 4.9585

4.5301 4.292a 4.549j 4.635k

4.995j

54.5 41k 3.2553 3.664k 3.0831 9.7025 64.6 3.4435 11.1240 29.8 4.3299 4.7884 4.821m 4.874n 3.4888 11.1104 38.1 4.9150

a

LAPW method with LDA + SOC and GGA + SOC [57]. LAPW method with PW91 + SOC [59]. c Norm-conserving pseudopotential method with SOC [63]. d Norm-conserving pseudopotential method without SOC [85]. e Mixed-basis pseudopotential method with SOC [69]. f PAW pseudopotential method with SOC [65]. g PAW pseudopotential method with SOC [66]. h LAPW method with SOC [60]. i PAW pseudopotential method with SOC and GGA of PW91 [86]. j PAW pseudopotential method without SOC [54]. AFM order and Liechtenstein’s DFT + U scheme [40] with U = 4.0 and J = 0.7 eV are used. k LAPW method with SOC [61]. l PAW pseudopotential method with SOC [34]. The parameters U = 0.94 eV and J = 0.58 eV are used. m LAPW method with SOC [62]. n LAPW method with SOC [64]. b

overbinding. As pointed in Ref. [87] and confirmed here, the performance of PBEsol is similar to LDA for actinide metals. For systems with slowly varying electronic densities, LDA and PBEsol are a good choice but not for the actinide metals in which the electronic density varies rapidly. This discrepancy could be remedied by GGA of PBE [57–65,56,66–69] or the inclusion of Hubbard U [51,54,34]. This point is confirmed by our calculation in Table 3. Now it is necessary to compare LDA/PBEsol + ULR with PBE. To clearly illustrate their performance, we use the data of lattice parameters in Table 3 to plot the deviations of PBE and PBEsol + ULR from the experimental results in Fig. 4. It can be seen that both PBE and LDA/ PBEsol + ULR could well reproduce the experimental lattice parameters for Th, Pa, U, Cm, AFM δ -Pu and Am. But for the internal parameter x1 of α -Np, the lattice constant c of β -Np, and a of γ -Np and ∊-Pu, PBE results deviate from the experimental values while LDA/PBEsol + ULR performs well. Clearly, LDA/PBEsol + U performs better than PBE on the lattice parameters of low symmetry structures of actinide metals. Even for α -U, there is also little improvement on the prediction of internal parameter y. In addition, for the lattice parameters of DM Pu and Am, PBE results differ strongly from the experiment while LDA/ PBEsol + ULR preforms much better. Artificial magnetic moment [70] was included in the calculation [10,67,34] to match the experiment, which can also be seen from our calculation of AFM δ -Pu and Am in Table 3. Here for DFT + U, the inclusion of artificial magnetic moment

Fig. 4. The deviation of lattice parameters calculated using PBE + SOC, LDA + ULR +SOC, and PBEsol + ULR +SOC from the experimental values. The data are presented in Table 3. The bar plots in the retangles represents the data from AFM systems.

is not necessary. This improvement of DFT + U over PBE in the DM system was already pointed out by Amadon and Dorado [34] and confirmed here. Compared to the large deviation of PBE (0–28%), the largest deviation of LDA + ULR and PBEsol + ULR is only 7.5% and 4.1%, respectively. This manifests the improvement of LDA/ PBEsol + ULR with respect to PBE in the prediction of lattice parameters 7

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

of structural properties and it is found that LDA/PBEsol + U performs better than PBE for diamagnetic actinide metals. Nevertheless, DFT + U approach has many inherent limits such as the static character and the calculated U is not frequency-dependent. In addition, the constraint of zero magnetic moment is questionable since diamagnetic order is not the lowest-energy magnetic states with DFT + U for U, Np, Pu, and Am. But for Pu and Am, the exclusion of artificial magnetic moment indicates a small improvement. And due to its internal self-consistency, DFT + U with U calculated from linear response approach may be a useful computational tool to model actinide materials.

of actinide metals. As concern the bulk moduli B, both PBE and LDA/PBEsol + ULR could reproduce the experimental values except the controversial data for α -Pa [79,88]. For α -Th and α -Cm, PBE performs better than LDA/ PBEsol + ULR and for DM α -Np, γ -Pu, and δ -Pu, LDA/PBEsol + ULR is better. For other systems, the deviation of LDA/PBEsol + ULR is almost the same as that of PBE. However, at this stage one cannot conclude that LDA/PBEsol + ULR could provide a better description of actinide metals than PBE. On the one hand, the structural properties are calculated under the constraint of DM which is not the lowest-energy magnetic order. The complex structures β -U, α - and β -Pu are not considered here. On the other hand, the physical properties of actinide metals were already well reproduced by PBE in many theoretical calculation [57–65,56,66–69,71] and a thorough comparison of the calculated magnetic, spectral, phonon properties, and high-pressure behavior using LDA/PBEsol + ULR with the literature remains to be done. These studies of each actinide metals will be pursued in the near future.

CRediT authorship contribution statement Ruizhi Qiu: Conceptualization, Investigation, Software, Validation, Writing - original draft. Bingyun Ao: Funding acquisition, Writing review & editing. Li Huang: Funding acquisition, Conceptualization, Writing - review & editing. Data availability

5. Conclusion The raw/processed data required to reproduce these findings are available to download from [ https://doi.org/10.17632/784tzvwwx8. 2].

In summary, we perform the linear response calculation of the effective Coulomb interaction U between 5f electrons of pure actinides. The calculated U is found to seldom rely on the exchange-correlation functional, spin-orbital coupling, and magnetic state, but depend on the neighbor distance and element. In particular, the values of U vary from atom to atom for β -U, which is interpreted using the variation of the nearest-neighbor distance. The jump of U from Pu to Am is similar to that of atomic volume, implying the itinerant-localized change along this series. These values of U have been used in the DFT + U calculation

Acknowledgments We would like to acknowledge the financial support from the Science Challenge Project of China (Grant No. TZ2016004), National Natural Science Foundation of China (Grant Nos. 21771167, 11874329, and 21601167), and the CAEP project (Grant No. TCGH0708).

Appendix A. Actinide metals In this appendix we will show our computational models, i.e., the solid phases of actinide metals at ambient pressure. All the experimental lattice parameters are obtained from Pearson’s handbook [11] except β -Pa [12] and β -Cf [13]. Thorium (Th) is the first element in the actinides series with empty 5f orbitals for free atom but a substantial occupation of 5f orbitals in its metallic condensed phase. Under ambient condition, Th has a close-packed face-centered cubic (fcc) structure (α -Th) and at elevated temperature, Th transforms from fcc to body-centered cubic (bcc) at approximately 1673 K and bcc Th (β -Th) melts at approximately 2023 K [89]. The narrow 5f bands, which could induce a Peierls distortion [90], don’t take effect at ambient pressure but play a role under high pressure. That is, fcc Th transforms toward a body centered tetragonal (bct) structure at about 60 GPa [91]. Protactinium (Pa) is the second element of actinides series with the electron configuration of free atom being [Rn]7s 2 5 f 2 6d1. The 5f electrons begin to play an important role and participate in the metallic bonding. A low-symmetry bct structure is adopted by the solid-state phase of Pa under ambient condition (α -Pa). The melting point of Pa is about 1845 K and when approaching the melting point, there is another solid-state phase, β -Pa. The crystal structure of β -Pa is controversial [88]. Marples prepared Pa metal by reducing the tetrafluoride with calcium and predicted a bcc form of β -Pa from the extrapolation of thermal expansion data [12]. Asprey et al. found a new fcc structure in a quenched arc-melted sample besides the bct structure of α -Pu [92]. This fcc form of β -Pa was confirmed by the reversible transition between bct and fcc phases above 1473 K using both X-ray and impurity analyses [93]. For more details and discussion, one can refer to Ref. [88]. Since the experimental lattice parameter of fcc Pa (5.018 Å) is much larger than any theoretical values even within DFT + U, bcc is considered as the structure of β -Pa here. As the third element in the actinide series, uranium (U) has many very unique features that result from its 5f electrons. For example, the roomtemperature orthorhombic crystal structure and the low-temperature charge density wave (CDW) transition of U are unique for an element at ambient pressure. Neglecting the CDW phase, U exists in three solid-state phases at ambient pressure, which are labeled as α, β , γ -U and shown in Fig. A.5. α -U is orthorhombic with space group Cmcm (No. 63) and all atoms are located at 4c Wyckoff positions (0, y, 1/4). The structure is parameterized by the three lattice constants a, b, c , and the internal parameter y. β -U is tetragonal with space group P 4 2/ mnm (No. 136) and all atoms are located at five Wyckoff positions with seven internal parameters. There are thirty atoms in the unit cell of β -U. γ -U is bcc structure, which is adopted by most of metals when approaching the melting point for reasons of lattice stability. At elevated temperatures, U transforms from α to β at approximately 941 K and β transforms to γ at approximately 1048 K. Fig. A.5. Crystal structures of uranium.

8

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

Fig. A.6. Crystal structures of neptunium.

Fig. A.7. Crystal structures of plutonium.

Neptunium (Np) is the fourth element in the series, which also exists in three solid-state phases at ambient pressure: orthorhombic α , tetragonal β , bcc γ . The structures are shown in Fig. A.6. The room-temperature α -Np structure has space group Pnma (No. 62), which is the subgroup of Cmcm. All the Np atoms are also located at two 4c Wyckoff positions ( x1,2 , 1/4, z1,2 ). Thus a Peierls distortion along a-direction could make α -U structure transition toward α -Np structure, which has been theoretically demonstrated by our previous research [68]. The α -Np structure is parameterized with the three lattice constants a, b, c , and four internal parameters. The space group of β -Np is P4212 (No. 90), which also has a lower symmetry than β -U. But the structure of β -Np is much simpler. The Np atoms are located at 2a Wyckoff positions (0, 0, 0) and 2c Wyckoff positions (0, 1/2, z). Only two lattice constants a, c and one internal parameter z are enough to characterize the β -Np structure. At about 551 K, α -Np transforms to β -Np and β -Np transforms to bcc Np (γ -Np) at about 823 K [94]. As the fifth element of actinides, plutonium (Pu) is the most complex element in the periodic table. The electron configuration of Pu atom is [Rn] 7s 2 5 f 6 , which implies that only one electron is needed to reach a stable electron configuration, i.e., half-filled f-shell. On the other hand, there is spin–orbit splitting in the 5f states and the j = 5/2 states are filled for Pu atom. Thus it has many states close to each other in energy but dramatically different in crystal structure, which is adopted in response to minor changes in its surroundings. Before it melts at about 913 K, Pu undergoes six different phases, which is more than any other element. See Fig. A.7 for a sense of complexity. α -Pu is monoclinic with space group P 21/ m (No. 11) and all atoms are located at eight 2e Wyckoff positions ( x1 … 8 , 1/4, z1 … 8 ). α -Pu transforms to β -Pu at about 395 K and β -Pu is also monoclinic with space group C 2/ m (No. 12). β -Pu has thirty-four atoms per unit cell, i.e., the largest unit cell of pure metals. β -Pu transforms to γ -Pu at about 479 K and γ -Pu is face-centered orthorhombic with space group Fddd (No. 70). Pu transforms from γ phase to fcc δ phase at about 585 K and δ -Pu has the lowest density. δ -Pu transforms to bct δ ′-Pu at about 724 K and δ ′-Pu transforms to bcc ∊-Pu at about 758 K. For the transplutonium elements, i.e., americium (Am), curium (Cm), berkelium (Bk), and californium (Cf), all the metals are known to have a double hexagonal closed-packed (dhcp) α form and a high temperature fcc β form at normal pressure. There is a possible hexagonal closed-packed (hcp) phase of Cf, which is controversial [13] and not considered here. Research interests of transplutonium metals are focused on their highpressure behavior since the pressure could induce the transition from localization to delocalization of 5f electrons. As can be seen from Fig. 1, the atomic volume of α -phase is nearly equal to or even greater than that of β -phase. This behavior results from their localized 5f electrons at normal pressure.

with a transition temperature above 18 K, Nature 420 (6913) (2002) 297–299. [8] J.A. Mydosh, P.M. Oppeneer, Colloquium: hidden order, superconductivity, and magnetism: The unsolved case of URu2Si2, Rev. Mod. Phys. 83 (2011) 1301–1322. [9] P. Santini, S. Carretta, G. Amoretti, R. Caciuffo, N. Magnani, G.H. Lander, Multipolar interactions in f-electron systems: the paradigm of actinide dioxides, Rev. Mod. Phys. 81 (2009) 807–863. [10] K.T. Moore, G. van der Laan, Nature of the 5f states in actinide metals, Rev. Mod. Phys. 81 (2009) 235–298. [11] P. Villars, Pearson’s Handbook of Crystallographic Data for Intermediate Phases, 1997. [12] J.A.C. Marples, On the thermal expansion of protactinium metal, Acta Crystallogr. 18 (4) (1965) 815–817. [13] S. Heathman, T. Le Bihan, S. Yagoubi, B. Johansson, R. Ahuja, Structural investigation of californium under pressure, Phys. Rev. B 87 (2013) 214111 . [14] L. Brewer, Energies of the electronic configurations of the lanthanide and actinide neutral atoms, J. Opt. Soc. Am. 61 (8) (1971) 1101–1111. [15] C. Herring, Magnetism: Exchange Interactions Among Itinerant Electrons, vol. IV, 1967. [16] B. Johansson, Nature of the 5f electrons in the actinide series, Phys. Rev. B 11 (1975) 2740–2743. [17] D. van der Marel, G.A. Sawatzky, Electron-electron interaction and localization in d

References [1] G. Lander, E. Fisher, S. Bader, The solid-state properties of uranium A historical perspective and review, Adv. Phys. 43 (1) (1994) 1–111. [2] M.E. Manley, B. Fultz, R.J. McQueeney, C.M. Brown, W.L. Hults, J.L. Smith, D.J. Thoma, R. Osborn, J.L. Robertson, Large harmonic softening of the phonon density of states of uranium, Phys. Rev. Lett. 86 (2001) 3076–3079. [3] M.E. Manley, M. Yethiraj, H. Sinn, H.M. Volz, A. Alatas, J.C. Lashley, W.L. Hults, G.H. Lander, J.L. Smith, Formation of a new dynamical mode in α-uranium observed by inelastic X-ray and neutron scattering, Phys. Rev. Lett. 96 (2006) 125501 . [4] G.H. Lander, Sensing electrons on the edge, Science 301 (5636) (2003) 1057–1059. [5] X. Dai, S.Y. Savrasov, G. Kotliar, A. Migliori, H. Ledbetter, E. Abrahams, Calculated phonon spectra of plutonium at high temperature, Science 300 (5621) (2003) 953–955. [6] J. Wong, M. Krisch, D.L. Farber, F. Occelli, A.J. Schwartz, T.-C. Chiang, M. Wall, C. Boro, R. Xu, Phonon dispersions of fcc δ-plutonium-gallium by inelastic X-ray scattering, Science 301 (5636) (2003) 1078–1080. [7] J.L. Sarrao, L.A. Morales, J.D. Thompson, B.L. Scott, G.R. Stewart, F. Wastin, J. Rebizant, P. Boulet, E. Colineau, G.H. Lander, Plutonium-based superconductivity

9

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

[49] J.H. Shim, K. Haule, G. Kotliar, Fluctuating valence in a correlated solid and the anomalous properties of δ-plutonium, Nature 446 (7135) (2007) 513–516. [50] C.A. Marianetti, K. Haule, G. Kotliar, M.J. Fluss, Electronic coherence in δ-Pu: a dynamical mean-field theory study, Phys. Rev. Lett. 101 (2008) 056403 . [51] S.Y. Savrasov, G. Kotliar, Ground state theory of δ-Pu, Phys. Rev. Lett. 84 (2000) 3670–3673. [52] J. Bouchet, B. Siberchicot, F. Jollet, A. Pasturel, Equilibrium properties of delta-Pu: LDA+ U calculations (LDA = local density approximation), J. Phys.: Condens. Matter 12 (8) (2000) 1723. [53] A.B. Shick, V. Drchal, L. Havela, Coulomb-U and magnetic-moment collapse in δ-Pu, Europhys. Lett. 69 (4) (2005) 588. [54] G. Jomard, B. Amadon, F. Bottin, M. Torrent, Structural, thermodynamic, and electronic properties of plutonium oxides from first principles, Phys. Rev. B 78 (2008) 075125 . [55] B. Dorado, F.M.C. Bottin, J. Bouchet, Phonon spectra of plutonium at high temperatures, Phys. Rev. B 95 (2017) 104303 . [56] P. Söderlind, G. Kotliar, K. Haule, P.M. Oppeneer, D. Guillaumont, Computational modeling of actinide materials and complexes, MRS Bull. 35 (2010) 883–888. [57] M.D. Jones, J.C. Boettger, R.C. Albers, D.J. Singh, Theoretical atomic volumes of the light actinides, Phys. Rev. B 61 (2000) 4644–4650. [58] Y. Wang, Y. Sun, First-principles thermodynamic calculations for δ-Pu and ∊-Pu, J. Phys.: Condens. Matter 12 (2000) L311–L316. [59] J. Kuneš, P. Novák, R. Schmid, P. Blaha, K. Schwarz, Electronic structure of fcc Th: spin-orbit calculation with 6p1/2 local orbital extension, Phys. Rev. B 64 (2001) 153102 . [60] P. Söderlind, First-principles elastic and structural properties of uranium metal, Phys. Rev. B 66 (2002) 085113 . [61] P. Söderlind, B. Sadigh, Density-functional calculations of α, β, γ , δ , δ′, and ∊ plutonium, Phys. Rev. Lett. 92 (2004) 185702 . [62] M. Pénicaud, Localization of 5f electrons and phase transitions in americium, J. Phys.: Condens. Matter 17 (2) (2005) 257–267. [63] J. Bouchet, F. Jollet, G. Zérah, High-pressure lattice dynamics and thermodynamic properties of Th: an ab initio study of phonon dispersion curves, Phys. Rev. B 74 (2006) 134304 . [64] D. Gao, A.K. Ray, Relativistic density functional study of fcc americium and the (111) surface, Phys. Rev. B 77 (2008) 035123 . [65] P. Souvatzis, T. Björkman, O. Eriksson, P. Andersson, M.I. Katsnelson, S.P. Rudin, Dynamical stabilization of the body centered cubic phase in lanthanum and thorium by phonon-phonon interaction, J. Phys.: Condens. Matter 21 (17) (2009) 175402. [66] K.O. Obodo, N. Chetty, First principles LDA+U and GGA+U study of protactinium and protactinium oxides: dependence on the effective U parameter, J. Phys.: Condens. Matter 25 (14) (2013) 145603. [67] P. Söderlind, First-principles phase stability, bonding, and electronic structure of actinide metals, J. Electron Spectrosc. Relat. Phenom. 194 (2014) 2–7. [68] R. Qiu, H. Lu, B. Ao, T. Tang, P. Chen, First-principles studies on the charge density wave in uranium, Modell. Simul. Mater. Sci. Eng. 24 (5) (2016) 055011 . [69] P. González-Castelazo, R. de Coss-Martínez, O. De la Peña Seaman, R. Heid, K.P. Bohnen, Electron-phonon coupling and superconductivity in the light actinides: a first-principles study, Phys. Rev. B 93 (2016) 104512 . [70] J.C. Lashley, A. Lawson, R.J. McQueeney, G.H. Lander, Absence of magnetic moments in plutonium, Phys. Rev. B 72 (2005) 054416 . [71] P. Söderlind, A. Landa, B. Sadigh, Density-functional theory for plutonium, Adv. Phys. 68 (1) (2019) 1–47. [72] A.B. Shick, J. Kolorenc, J. Rusz, P.M. Oppeneer, A.I. Lichtenstein, M.I. Katsnelson, R. Caciuffo, Unified character of correlation effects in unconventional Pu-based superconductors and δ-Pu, Phys. Rev. B 87 (2013) 020505 . [73] N. Lanatà, Y. Yao, C.-Z. Wang, K.-M. Ho, G. Kotliar, Phase diagram and electronic structure of praseodymium and plutonium, Phys. Rev. X 5 (2015) 011008 . [74] A.J. Arko, J.J. Joyce, L. Morales, J. Wills, J. Lashley, F. Wastin, J. Rebizant, Electronic structure of α- and δ-Pu from photoelectron spectroscopy, Phys. Rev. B 62 (2000) 1773–1779. [75] H. Hsu, K. Umemoto, M. Cococcioni, R. Wentzcovitch, First-principles study for low-spin LaCoO3 with a structurally consistent Hubbard U, Phys. Rev. B 79 (2009) 125124 . [76] F. Birch, Finite elastic strain of cubic crystals, Phys. Rev. 71 (1947) 809–824. [77] B. Dorado, B. Amadon, M. Freyss, M. Bertolus, DFT+U calculations of the ground state and metastable states of uranium dioxide, Phys. Rev. B 79 (2009) 235125 . [78] G. Bellussi, U. Benedict, W. Holzapfel, High pressure X-ray diffraction of thorium to 30 GPa, J. Less-Common Met. 78 (1) (1981) 147–153. [79] U. Benedict, J. Spirlet, C. Dufour, I. Birkel, W. Holzapfel, J. Peterson, X-ray diffraction study of protactinium metal to 53 GPa, J. Magn. Magn. Mater. 29 (1) (1982) 287–290. [80] C.-S. Yoo, H. Cynn, P. Söderlind, Phase diagram of uranium at high pressures and temperatures, Phys. Rev. B 57 (1998) 10359–10362. [81] S. Dabos, C. Dufour, U. Benedict, M. Pags, Bulk modulus and P-V relationship up to 52 GPa of neptunium metal at room temperature, J. Magn. Magn. Mater. 63–64 (1987) 661–663. [82] P. Söderlind, A. Landa, J.E. Klepeis, Y. Suzuki, A. Migliori, Elastic properties of Pu metal and Pu-Ga alloys, Phys. Rev. B 81 (2010) 224110 . [83] S. Heathman, R.G. Haire, T. Le Bihan, A. Lindbaum, K. Litfin, Y. Méresse, H. Libotte, Pressure induces major changes in the nature of americium’s 5f electrons, Phys. Rev. Lett. 85 (2000) 2961–2964. [84] S. Heathman, R.G. Haire, T. Le Bihan, A. Lindbaum, M. Idiri, P. Normile, S. Li, R. Ahuja, B. Johansson, G.H. Lander, A high-pressure structure in curium linked to magnetism, Science 309 (5731) (2005) 110–113. [85] J. Bouchet, R.C. Albers, Elastic properties of the light actinides at high pressure, J.

and f transition metals, Phys. Rev. B 37 (1988) 10674–10684. [18] J.F. Herbst, R.E. Watson, I. Lindgren, Coulomb term U and 5f electron excitation energies for the metals actinium to berkelium, Phys. Rev. B 14 (1976) 3265–3272. [19] T. Bandyopadhyay, D.D. Sarma, Calculation of Coulomb interaction strengths for 3d transition metals and actinides, Phys. Rev. B 39 (1989) 3517–3521. [20] A.V. Lukoyanov, M.O. Zaminev, V.I. Anisimov, Pressure-induced modification of the electron structure of metallic thorium, J. Exp. Theor. Phys. 118 (1) (2014) 148–152. [21] W. Xie, W. Xiong, C.A. Marianetti, D. Morgan, Correlation and relativistic effects in U metal and U-Zr alloy: validation of ab initio approaches, Phys. Rev. B 88 (2013) 235128 . [22] B. Amadon, First-principles DFT+DMFT calculations of structural properties of actinides: role of Hund’s exchange, spin-orbit coupling, and crystal structure, Phys. Rev. B 94 (2016) 115148 . [23] F. Aryasetiawan, M. Imada, A. Georges, G. Kotliar, S. Biermann, A.I. Lichtenstein, Frequency-dependent local interactions and low-energy effective models from electronic structure calculations, Phys. Rev. B 70 (2004) 195104 . [24] F. Aryasetiawan, K. Karlsson, O. Jepsen, U. Schönberger, Calculations of Hubbard U from first-principles, Phys. Rev. B 74 (2006) 125106 . [25] O. Gunnarsson, O.K. Andersen, O. Jepsen, J. Zaanen, Density-functional calculation of the parameters in the Anderson model: application to Mn in CdTe, Phys. Rev. B 39 (1989) 1708–1722. [26] V.I. Anisimov, O. Gunnarsson, Density-functional calculation of effective Coulomb interactions in metals, Phys. Rev. B 43 (1991) 7570–7574. [27] M. Cococcioni, S. de Gironcoli, Linear response approach to the calculation of the effective interaction parameters in the LDA+U method, Phys. Rev. B 71 (2005) 035105 . [28] G. Kotliar, S.Y. Savrasov, K. Haule, V.S. Oudovenko, O. Parcollet, C.A. Marianetti, Electronic structure calculations with dynamical mean-field theory, Rev. Mod. Phys. 78 (2006) 865–951. [29] B. Amadon, Erratum: first-principles DFT+DMFT calculations of structural properties of actinides: role of Hund’s exchange, spin-orbit coupling, and crystal structure [Phys. Rev. B 94, 115148 (2016)], Phys. Rev. B 97 (2018) 039903. [30] T. Gouder, L. Havela, F. Wastin, J. Rebizant, Evidence for the 5f localisation in thin Pu layers, Europhys. Lett. 55 (5) (2001) 705. [31] J.R. Naegele, L. Manes, J.C. Spirlet, W. Müller, Localization of 5f electrons in americium: a photoemission study, Phys. Rev. Lett. 52 (1984) 1834–1837. [32] B. Amadon, F. Jollet, M. Torrent, γ and β cerium: LDA+U calculations of groundstate parameters, Phys. Rev. B 77 (2008) 155104. [33] J.P. Perdew, K. Burke, M. Ernzerhof, Generalized gradient approximation made simple, Phys. Rev. Lett. 77 (1996) 3865–3868. [34] B. Amadon, B. Dorado, A unified and efficient theory for the structural properties of actinides and phases of plutonium, J. Phys.: Condens. Matter 30 (40) (2018) 405603. ∂E = ∊ in density-functional theory, Phys. Rev. B 18 (1978) [35] J.F. Janak, Proof that ∂ni 7165–7168. [36] B. Himmetoglu, A. Floris, S. Gironcoli, M. Cococcioni, Hubbard-corrected DFT energy functionals: the LDA+U description of correlated systems, Int. J. Quantum Chem. 114 (1) (2014) 14–49. [37] D.M. Ceperley, B.J. Alder, Ground state of the electron gas by a stochastic method, Phys. Rev. Lett. 45 (1980) 566–569. [38] J.P. Perdew, A. Zunger, Self-interaction correction to density-functional approximations for many-electron systems, Phys. Rev. B 23 (1981) 5048–5079. [39] J.P. Perdew, A. Ruzsinszky, G.I. Csonka, O.A. Vydrov, G.E. Scuseria, L.A. Constantin, X. Zhou, K. Burke, Restoring the density-gradient expansion for exchange in solids and surfaces, Phys. Rev. Lett. 100 (2008) 136406 . [40] A.I. Liechtenstein, V.I. Anisimov, J. Zaanen, Density-functional theory and strong interactions: orbital ordering in Mott-Hubbard insulators, Phys. Rev. B 52 (1995) R5467–R5470. [41] S.L. Dudarev, G.A. Botton, S.Y. Savrasov, C.J. Humphreys, A.P. Sutton, Electronenergy-loss spectra and the structural stability of nickel oxide: an LSDA+U study, Phys. Rev. B 57 (1998) 1505–1509. [42] G. Kresse, J. Furthmüller, Efficient iterative schemes for ab initio total-energy calculations using a plane-wave basis set, Phys. Rev. B 54 (1996) 11169–11186. [43] P.E. Blöchl, Projector augmented-wave method, Phys. Rev. B 50 (1994) 17953–17979. [44] G. Kresse, D. Joubert, From ultrasoft pseudopotentials to the projector augmentedwave method, Phys. Rev. B 59 (1999) 1758–1775. [45] K. Lejaeghere, G. Bihlmayer, T. Björkman, P. Blaha, S. Blügel, V. Blum, D. Caliste, I. E. Castelli, S.J. Clark, A. Dal Corso, S. de Gironcoli, T. Deutsch, J.K. Dewhurst, I. Di Marco, C. Draxl, M. Dułak, O. Eriksson, J.A. Flores-Livas, K.F. Garrity, L. Genovese, P. Giannozzi, M. Giantomassi, S. Goedecker, X. Gonze, O. Grånäs, E.K.U. Gross, A. Gulans, F. Gygi, D.R. Hamann, P.J. Hasnip, N.A.W. Holzwarth, D. Iuşan, D.B. Jochym, F. Jollet, D. Jones, G. Kresse, K. Koepernik, E. Küçükbenli, Y.O. Kvashnin, I.L.M. Locht, S. Lubeck, M. Marsman, N. Marzari, U. Nitzsche, L. Nordström, T. Ozaki, L. Paulatto, C.J. Pickard, W. Poelmans, M.I.J. Probert, K. Refson, M. Richter, G.-M. Rignanese, S. Saha, M. Scheffler, M. Schlipf, K. Schwarz, S. Sharma, F. Tavazza, P. Thunström, A. Tkatchenko, M. Torrent, D. Vanderbilt, M.J. van Setten, V. Van Speybroeck, J.M. Wills, J.R. Yates, G.-X. Zhang, S. Cottenier, Reproducibility in density functional theory calculations of solids, Science 351 (2016) 6280. [46] H.J. Monkhorst, J.D. Pack, Special points for Brillouin-zone integrations, Phys. Rev. B 13 (1976) 5188–5192. [47] S.Y. Savrasov, G. Kotliar, E. Abrahams, Correlated electrons in δ-plutonium within a dynamical mean-field picture, Nature 410 (6830) (2001) 793–795. [48] S.Y. Savrasov, K. Haule, G. Kotliar, Many-body electronic structure of americium metal, Phys. Rev. Lett. 96 (2006) 036404 .

10

Computational Materials Science 171 (2020) 109270

R. Qiu, et al.

the crystal structures of metals, Nature 374 (1995) 524–525. [91] Y.K. Vohra, J. Akella, 5f bonding in thorium metal at extreme compressions: phase transitions to 300 GPa, Phys. Rev. Lett. 67 (1991) 3563–3566. [92] L. Asprey, R. Fowler, J. Lindsay, R. White, B. Cunningham, A new high-temperature phase of protactinium metal, Inorg. Nucl. Chem. Lett. 7 (10) (1971) 977–980. [93] J. Bohet, W. Müller, Preparation and structure studies of Van Arkel protactinium, J. Less-Common Met. 57 (2) (1978) 185–199. [94] W.H. Zachariasen, Crystal chemical studies of the 5f-series of elements. XVIII. Crystal structure studies of neptunium metal at elevated temperatures, Acta Crystallogr. 5 (5) (1952) 664–667.

Phys.: Condens. Matter 23 23 (21) (2011) 215402 . [86] C.D. Taylor, Evaluation of first-principles techniques for obtaining materials parameters of α-uranium and the (001) α-uranium surface, Phys. Rev. B 77 (2008) 094119 . [87] P. Söderlind, A. Gonis, Assessing a solids-biased density-gradient functional for actinide metals, Phys. Rev. B 82 (2010) 033102 . [88] H. Blank, Experimental results on properties of Pa revisited, J. Alloys Compd. 343 (1–2) (2002) 108–115. [89] P. Chiotti, High temperature crystal structure of thorium, J. Electrochem. Soc. 101 (11) (1954) 567–570. [90] P. Söderlind, O. Eriksson, B. Johansson, J.M. Wills, A.M. Boring, A unified picture of

11