Pre-melting hcp to bcc Transition in Beryllium

0.73 MB·

by  Y. Lu

Pre-melting hcp to bcc Transition in Beryllium
Y. Lu,1 T. Sun,2 Ping Zhang,3 P. Zhang,4,1 D.-B. Zhang,1,∗ and R. M. Wentzcovitch5
1Beijing Computational Science Research Center, Beijing 100193, China
2Key Laboratory of Computational Geodynamics,
University of Chinese Academy of Sciences, Beijing 100049, China
3Institute of Applied Physics and Computational Mathematics, Beijing 100088, China
4Department of Physics, University at Buffalo, State University of New York, Buffalo, New York 14260, USA
5Department of Chemical Engineering and Materials Science,
7 and Minnesota Supercomputing Institute, University of Minnesota, Minneapolis, MN, 55455, USA
1
0 Beryllium (Be) is an important material with wide applications ranging from aerospace compo-
2 nents to X-ray equipments. Yet a precise understanding of its phase diagram remains elusive. We
have investigated the phase stability of Be using a recently developed hybrid free energy compu-
n
tation method that accounts for anharmonic effects by invoking phonon quasiparticles. We find
a
that the hcp → bcc transition occurs near the melting curve at 0 < P < 11 GPa with a positive
J
Clapeyronslopeof41±4K/GPa. Thebccphaseexistsinanarrowtemperaturerangethatshrinks
4
withincreasingpressure,explainingthedifficultyinobservingthisphaseexperimentally. Thiswork
1
alsodemonstratesthevalidityofthistheoreticalframework basedonphononquasiparticletostudy

structural stability and phase transitions in strongly anharmonic materials.
]
i
c
PACSnumbers: 63.20.Ry,81.30.Bx,61.50.Ks,63.20.D-,64.70.K-
s
-
l
r
t Elementalsolidsusuallyundergoaseriesofphasetran- used quasi-harmonic approximation (QHA) and Debye
m
sitionsfromambientconditionstoextremeconditions[1– modelarenotabletocapturesucheffect[23–29]. Forthis
. 3]. Knowledge of their phase diagrams is a prerequisite reason, bcc Be and the associated hcp/bcc phase transi-
t
a
for establishing their equations of state (EOS), a funda- tion remain poorly understood for P <11 GPa (density
m
mental relation for determining thermodynamics prop- < 2.1g/cc) where bcc Be is dynamically unstable at 0
- erties and processes at high pressures and temperatures K [24]. At P >11 GPa, bcc Be is dynamically stabilized
d
n (PT). However, resolving phase boundaries is challeng- bypressureandtheQHAmight,inprinciple,beapplied.
o ing for experiments given the uncertainties from several However,the hcp/bccboundary [23–25] predicted by the
c sources, especially at very high PT. Beryllium (Be) is QHA does not agree with experiments [13], suggesting
[
a typical system whose phase diagram remains an open that anharmonic effects still play an important role at

1 problem despite intense investigations. It assumes a higher pressures.
v hexagonal close-packed (hcp) structure at relatively low
7 In this Letter, we report a new investigation of the
T [1]. Near the melting temperature T ( 1,550 K
6 M ∼ phase stability of bcc Be and the associated hcp bcc
8 at 0 GPa), a competing phase with the body-centered phase transition boundary up to 30 GPa and tem→pera-
3 cubic symmetry (bcc) seems to emerge [4–6]. However,
tures up to 2,000 K. We have used a recently developed
0 not all experiments [7–13] have observed this bcc phase,
hybrid approach [30, 31] that combines first-principles
1. causing confusion and controversies. Be is important for
molecular dynamics (MD) and lattice dynamics calcula-
0 bothfundamentalresearch[14–17]andpracticalapplica-
tions to address anharmonic effects in the free energy.
7 tions. Being astrongandlight-weightmetal, ithas been
1 widelyusedinabroadrangeoftechnologicalapplications In this method, the concept of phonon quasiparticles of-
: fers a quantitative characterization of the effects of lat-
v in harsh environments and extreme PT conditions, e.g.,
tice anharmonicity [32, 33]. We show that Be exhibits
i up to T > 4,000 K and P > 200 GPa [18–22].
X pronounced anharmonic effects in both the bcc and hcp
The bcc phase of Be was directly observed [4, 5] only
r phases. Specifically,ourresultsrevealthedynamicalsta-
a at T > 1,500 K around ambient pressure before melt-
bilization of bcc Be with increasing T. The bcc phase,
ing. Measurements of the temperature dependent resis-
however,isfavorableonlyinanarrowtemperaturerange
tance suggestedthat bcc Be is a high pressurephase and
near T , with hcp bcc phase boundary having a pos-
M
the hcp/bcc phase boundary between 0 < P < 6 GPa →
itive Clapeyron slope of 41 4 K/GPa. The bcc sta-
has a negative Clapeyron slope ( 52 8 K/GPa) [6]. ±
− ± bility field shrinks with increasing pressure and eventu-
However, recent experiments have challenged this con-
ally disappears at around 11 GPa. This result seems
clusion [8, 9, 11–13]. For example, it was reported that
to agree overall with early experiments [4, 5] and dif-
bcc Be was not observed for 8 < P < 205 GPa and
fersfromother hcp/bccphase boundaries(e.g., Mg [34]),
300 < T < 4,000 K [13]. On the theory side, the
usually displaying a negative Clapeyron slope.
study of Be’s phase diagram using conventional meth-
ods encounters significant difficulties. The lattice dy- Inthepresentapproach,phononquasiparticlesarenu-
namics of bcc Be is highly anharmonic, and the widely merically characterized by the mode projected velocity
2
autocorrelationfunction [30, 31],
hV(0)·V(t)iq,s =t0l→im∞t10 Z0t0Vq∗,s(t′)·Vq,s(t′+t)dt′, (a )/K)V(t)k(q,sB4000 400K ( b )/K)V(t)k(q,sB10000 1,000K
where V (t) = N √M v(t) ˆǫi exp(iq R ) is t(h1e) V(0)-400 V(0)-1000
q,s i=1 i · q,s · i 0 1 2 0 1 2
mode projected aPnd mass weighted velocity for normal t (p s) t (ps)
mode (q,s) with wave vector q; v (t)(i = 1,...,N) is the (c) (d)
i 400K 1,000K
alattoimonisc wveitlohciNty aptroomdus,ceadndbyMfiirsatn-dpriRnciipalreestMheDatsoimmuic- Gq,s Gq,s
mass and coordinate of atom i, respectively. ǫˆi (i =
q,s
1,...,N) is the polarization vector of normal mode (q,s)
calculated using density functional perturbation theory 400 500 600 700 400 500 600 700
(cm-1) ( cm-1)
(DFPT) [35]. For a well-defined phonon quasiparticle,
FIG. 1. Mode projected velocity auto-correlation func-
the velocity autocorrelation function displays an oscilla-
tions of the TA1 acoustic mode (q,s) at q = N with a har-
torydecayingbehavioranditsFouriertransform,i.e.,the monic frequency of 599 cm−1 of bcc Be at (a) 400 K and (b)
power spectrum,
1,000 K, respectively. (c) and (d) show the corresponding
∞ power spectra. In (c), the shaded area between 535 cm−1
G = V(0) V(t) exp(iωt)dt, (2) and 580 cm−1 covers two major peaks, indicating the break-
q,s Z h · iq,s downofthephononquasiparticlepicture. In(d),thevertical
0
dashed line at 556 cm−1 indicates the frequency of the well-
should have a Lorentzian-type line shape [30, 31]. The definedphonon quasiparticle.
renormalized phonon frequency ω , and the linewidth,
q,s
Γ canthen be obtainedas discussedin more details in
q,s
e
the Supplemental Materials. and normal modes were calculated using density func-
The concept of phonon quasiparticle reduces the com- tional perturbation theory (DFPT) [35].
plex problem of interacting anharmonic phonons to an Before proceeding, we should clarify the general un-
effective non-interacting system [32], such that the con- derstanding of the hcp bcc transition. The low T
ventional kinetic gas model and, to a great extent, the and low P hcp structur→e relates to the bcc structure
theory of harmonic phonons are still applicable. More- through the zone center transverse optical (TO) mode
over,sincestructuralphasetransitionistriggeredbylat- andamacroscopicstrain. Thismodeconsistsofopposite
tice vibrationsformanycases,insightintothe transition displacements of neighboring (0001) planes along 1010
mechanism can be also obtained by monitoring the vari- and softens with increasing T. This is not necehssarilyi
ation of frequencies and line widths of phonon quasipar- a soft mode transition, but generally a first order tran-
ticles. sition with negative Clapeyron slope [34]. The (0001)
We used the Vienna ab initio simulation package plane transforms into the (110) plane of the bcc phase.
(VASP) [37, 38] to carry out first principles MD simu- Thispicturewasvalidatedbyanearlyvariablecellshape
lations on 4 4 4 supercells (128 atoms) of Be. We molecular dynamics study [45]. The opposite bcc hcp
used the gene×rali×zed gradient approximation of Perdew, transition involves the lowest transverse acoustic→mode
Burke, and Ernzerhof [39] and the projector-augmented (TA2) at q = [1/2,1/2,0], the N point of the Brillouin
wave method [36] with an associated plane-wave basis zone, marked by an open circle in Fig. 2(a). With this
set energy cutoff of 350 eV. For metallic Be, the finite in mind we monitor closely the behavior of these modes
temperature Mermin functional [47, 48] was used. Sim- with changing T.
ulations were carried out at a series of volumes (V):
Wefirstinvestigatethebehaviorofphononquasiparti-
6.21 < V < 8.71˚A3/atom for bcc Be and 6.35 < V < clesatdifferenttemperatures. ForbccBe,phononquasi-
8.83˚A3/atomforhcpBe. ForhcpBe,properaspectratio particles are not well defined at low T as found in other
(c/a)isadoptedtoobtaingoodhydrostaticconditionsfor systems [32], but recovered at high T for unstable (soft)
specific volume and temperature. Temperatures ranging modes e.g., TA2. The analysis of this mode is shown in
from 300to 2,800K are controlledthroughthe Nos´e dy- Fig. S2 of the Supplementary material. It is more inter-
namics [41]. The considered volumes and temperatures esting to notice that at low T, phononquasiparticles are
result in a pressure range of 0 < P < 30 GPa. For not well defined evenfor certain stable modes with posi-
eachvolumeandtemperature,multipleindependentMD tive harmonic frequencies. Fig. 1 shows V(0) V(t)
q,s
h · i
runs(5parallelreplica)wereperformedtoimprovephase andthe correspondingpowerspectra ofthe TA1 phonon
space sampling quality that also allow for evaluation of mode at N calculated at 400 K (Fig. 1(a)) and 1,000
statistical uncertainties. Each MD run lasted 50 ps and K (Fig. 1(b)). This mode is marked by an open square
used a time step of 1 fs. Harmonic phonon frequencies in Fig. 2(a). Although V(0) V(t) at 400 K dis-
q,s
h · i
3
playsanoscillatorybehavior,theamplitudedecayisnon-
monotonic (Fig. 1(a)). Consequently, the power spec- ((aa))
600 bcc
trum has two major peaks within the shaded area as
shown in Fig. 1(c). This indicates that the frequency of
TA1
thismodecannotbewellconstrained,orequivalently,the 1)400
-m
corresponding phonon quasiparticle is not well-defined. c
(
In contrast, at 1,000K V(0) V(t) exhibits a nicely 200
q,s
decayingoscillatorybehahvior,F·ig.1(ib). Thecorrespond- TA2
ing power spectrum now has a well-defined Lorentzian 0
line shape with a single andwell defined peak, Fig.1(d). 1,000K
It is thus straightforward to identify the renormalized DFPT (0K)
-200
frequency of this mode as 556 cm−1. Similarly, all other H P N

quasiparticlemode frequenciessampledbythe 4x4x4su- (b)
percell are equally well defined. As previously indicated 600
hcp
[30, 31], these renormalized phonon frequencies plus the 1)
normalmodesenable the calculationofthe renormalized -m
force constant matrix and complete phonon dispersions. (c400 TO
This quantitative characterizationof phonon quasiparti-
cles and renormalized phonon dispersion provide a solid 200
foundation for studying thermal properties. 1,000K
DFPT (0K)
Figure 2(a) compares the anharmonic phonon disper- 0
sion of bcc Be at 1,000 K with the harmonic phonon dis- K M A

persion calculated using DFPT. Results are obtained at V = Veq (0K) P = 0 GPa
the static equilibrium volume of bcc Be, 7.81˚A3/atom. (c)480 (d)480
450 450
Therearenoticeabledifferencesbetweenthe anharmonic
astnadblhear(msooftn)icTpAh2onbornandcishpearlsoinongs.thIen pΓartiNculalirn,ethsetaubni-- -1cm)349200 hcp, TO, q = -1(cm) 34 9200 hcp, TO, q =
lizes when high temperature anharmoni−c effects are ac- (q,s210 q,s210
counted for. This indicates that bcc Be is stabilized by
180 180
anharmonic effects. To gain further insight into anhar-
moniceffects,weanalyzeT-dependentphononfrequency 150 bcc, TA2, q = N 150 bcc, TA2, q = N
shifts. Fig. 2(c) shows that the frequency of the bcc
0 400 800 1200 1600 0 400 800 1200 1600
zone edge phonon mode at q = N associated with the T(K) T(K)

TA2branchcalculatedatfixedvolumevariesnon-linearly
FIG.2. (a)Anharmonicphonondispersion at1,000 K(blue
with T. Lowest order many-body perturbation theory
solid curves) and harmonic phonon dispersion calculated us-
(MBPT) [42] predicts a linear frequency shift with T. ing DFPT (grey dashed curves) both at V = 7.81˚A3/atom,
Therefore, as expected, higher order anharmonic effects the static bcc Be equilibrium volume. The two transverse
ignoredin the perturbative treatmentplay animportant branches are labeled TA1 and TA2, respectively. (b) Anhar-
role here. monic phonon dispersion calculated at 1,000 K (blue solid
curves) and harmonic phonon dispersion calculated using
Thecalculatedanharmonicphonondispersionoverthe DFPT (grey dashed curves) both at V = 7.89˚A3/atom, the
whole Brillouin zone makes it possible to calculate the
static hcp Be equilibrium volume. (c) Temperature depen-
free energyinthe thermodynamiclimit (N ). Since dent frequency shifts of the TA2 q= N mode [open circle in
→∞
hcp Be is stable at low T for the entire pressure rangeof (a)] and of the TO q=Γ mode [open circle in (b)] calculate
interest, the lattice thermal properties have been stud- at constant volume and (d) at constant zero pressure. The
ied within the QHA [23–25] without further examina- procedureto convert from constant volumeto constant pres-
sure was described in [30]. The vertical dashed line in (d)
tion of the validity of the approximation. This naturally
indicates the hcp/bcc transition temperature.
brings up a question: how important are anharmonic ef-
fectsinthefreeenergyinthisseeminglystablestructure?
Fig.2(b)comparestheanharmonicphonondispersionat
T = 1,000 K and the harmonic phonon dispersion of quency shifts with increasingT canbe positive,negative
hcp Be calculatedata fixedvolume of 7.89˚A3/atomcor- or nearly zero, demonstrating the complexity of lattice
responding to zero static pressure. The differences, al- anharmonic effects. What is important here is that the
though not alarming, are still significant in most of the large frequency shifts in hcp Be should be incorporated
Brillouin zone. A detailed analysis of individual phonon intothefreeenergycalculationformoreaccurateevalua-
modesintheSupplementarymaterialrevealsthatthefre- tionsofthermodynamicpropertiesandphaseboundaries.
4
FIG.3. (a)FreeenergyF(V,T)versusvolumesofBeforbcc
(solid symbols)andhcp (openedsymbols) phasesat different
temperatures. (b) Free energy G(P,T) versus temperature
of Be for bcc (solid lines) and hcp (dashed lines) phases at FIG. 4. Phase diagram of Be. The melting line is adopted
different pressures. In (a) and (b), the error bars are too from Robert [24]. The experimental results for the hcp/bcc
small to bevisible. phase boundary (solid line [5] and dashed line [6]) and the
experimental point for bcc Be (opened square) [4] are shown
for comparison.
The large frequency shifts not only reveal pronounced
anharmoniceffectsbutalsoshedlightonthemicroscopic
mechanism of this phase transition. As mentioned ear- Supplementary material for more details). Our analy-
lier,the hcpandbccstructuresarerelatedbyacombina- sis of phonon quasiparticles demonstrates that they are
tion of phonon displacements and a macroscopic strain well defined for both phases for T 1,000 K, there-
≥
[43]. Together they provide a path for the hcp bcc fore the choice of T . F(V,T ) = E(V,T )+T S(V,T ),
→ 0 0 0 0 0
transition. The frequency of the zone center TO mode where E(V,T ) is the internal energy obtained from the
0
drops significantly from 467 to 405 cm−1 when T in- MDsimulation. Fig.3(a)displaysthecalculatedfreeen-
creasesfrom0to1,600K(seeFig.2(c)). Thisobservation ergies for both bcc and hcp phases. It is seen that at
is consistent with expectations based on the anticipated T 1,200K, V >V , and consistently, the common
bcc hcp
≥
transformation mechanism [34, 44, 45]. We note that tangent to these curves starts to have negative slope, in-
although the frequency shift is very large at 1,600 K, dicatingatransitionfromhcptobccatpositivepressure.
the picture of phonon quasiparticle is still valid for the We note that the volumetric variations of both phases
hcp phase (see the Supplementary material for detailed also provide clues to understand the predicted hcp/bcc
analysis). As mentioned earlier,fromthe Brillouinzone- transition at very high P (e.g., 400 GPa) [23–25]. More
folding relation,the correspondingmode in bcc Be is the detailsofthevariationofV andV isshowninFig.S6
bcc hcp
zoneedgeTA2 mode atq=N,whosepropertyisshown the Supplementary material.
in Fig. 2(c) and (d).
ItismoreconvenienttoconvertF(V,T)intoG(P,T)=
We now demonstrate that anharmonic effects are crit-
F(V,T)+P(V,T)V to obtain the phase boundary (see
ical for obtaining this hcp bcc phase boundary. When
→ the Supplementary material for details). Fig. 3(b) dis-
using anharmonic T-dependent phonon dispersions, the
playsG(P,T)forbothbccandhcpphases. AteachP,the
QHA free energy formula is no longer valid whereas the
intersectionofG andG givesthehcp/bcctransition
bcc hcp
entropyformulaisstillapplicable[46]. Therefore,wefirst
temperature. The resulting hcp bcc phase boundary
calculate the vibrational entropy [30, 32], →
shown in Fig. 4 together with the predicted/measured
melting line, reveal several important aspects: (i) bcc
S =k [(n +1)ln(n +1) n lnn ], (3)
vib B qs qs − qs qs Be is stable only when T approaches the melting point,
X
qs
in good agreement with direct experimental measure-
withn =[exp(¯hω˜ /k T) 1]−1,andobtainthe total ments. Forexample,atambientpressurewherethemelt-
qs q,s B
− ing temperature is 1550 K, bcc Be was observed at
free energy as:
∼
T > 1500 K [4, 5]. (ii) The hcp bcc phase boundary
→
T has a positive Clapeyron slope of 41 4 K/GPa, close
F(V,T)=F(V,T ) S(T′)dT′. (4) ±
0 −Z to the experimental value (43 7 K/GPa) reported by
T0 Abey [5]. Francois et al. rep±orted a negative Clapey-
where T is 1,000 K, S is the total entropy including ron slope through an indirect measurement of the T-
0
both vibrational and electronic contributions (See the dependent resistivity [6]. Our results call for future ex-
5
periments to clarify this issue. (iii) The hcp bcc 14, L1 (1984).
→
phase transition occurs only in a narrow pressure range [9] V.Vijayakumar,B.K.Godwal,Y.K.Vohra,S.K.Sikka,
of 0 < P < 11 2 GPa. Our results are consistent and R. Chidambaram, J. Phys. F: Met. Phys. 14, L65
± (1984).
with recent experiments that do not observe bcc Be at
[10] A. R.Marder, Science 142, 664 (1963).
high pressures [8–12]. One report did not observe the
[11] K. Nakano, Y. Akahama, and H. Kawamura, J. Phys.:
bcc phase at pressures as low as 8 GPa [13] whereas our
Condens. Matter. 14, 10569 (2002).
results suggest that the bcc phase should exist up to 11 [12] N.Velisavljevic,G.N.Chesnut,Y.K.Vohra,S.T.Weir,
GPa. However, the predicted temperature range over V.Malba,andJ.Akella,Phys.Rev.B65,172107(2002).
which the bcc phase is favorable is vanishingly narrow. [13] A. Lazicki, A. Dewaele, P. Loubeyre, and M. Mezouar,
Thus, it is possible that this narrow temperature range Phys. Rev. B 86, 174118 (2012).
[14] S. H.Glenzer, G. Gregori, R.W.Lee, F. J. Rogers, S.W.
was missed in experiments or that our calculations over-
Pollaine, and O. L. Landen,Phys. Rev. Lett. 90, 175002
estimate the stability field of the bcc phase by few GPa.
(2003).
In summary, using the concept of phonon quasiparti- [15] S. H. Glenzer, O. L. Landen, P. Neumayer et al., Phys.
cle, we have investigated the controversial hcp bcc Rev. Lett. 98, 065002 (2007).
phase boundary of Be. We find that bcc Be is sta→bilized [16] H.J.Lee,P.Neumayer,J.Castoretal.,Phys.Rev. Lett.
102, 115001 (2009).
at low pressures and high temperatures by anharmonic
[17] I. Vobornik, J. Fujii, M. Hochstrasser et al., Phys. Rev.
effects. For hcp Be, anharmonic effects on phonon prop-
Lett. 99, 166403 (2007)
ertiesarealsosignificant. Usinganharmonicphonondis- [18] S. W. Haan et al.,Phys. Plasmas 18, 051001 (2011).
persions, we evaluated the free energies of both phases [19] D. Swift, D. Paisley, and M. Knudon, AIP Conf. Proc.
and showed that the bcc phase emerges as a pre-melting 706, 119 (2003).
phenomenon at relatively low pressures. Our results for [20] K. L. Wilson, R. A Causey, W. L. Hsu, B. E. Mills, M.
the hcp bcc phase boundary are consistentwith most F. Smith, and J. B. Whitley, J. Vac. Sci. Technol. A 8,
→ 1750 (1990).
experimentalobservations[4,5,8–13]inallimportantas-
[21] D. S. Clark, S.W. Haan, and J. D. Salmonson, Phys.
pects. The temperature range over which the bcc phase
Plasmas 15, 056305 (2008).
exists shrinks with increasing pressure and vanishes at a
[22] D.C.Wilson,P.A.Bradley,N.M.Hoffmanetal.,Phys.
theoretical pressure of about 11 GPa. This narrow tem- Plasmas 5, 1953 (1998); J. L. Kline, S. A. Yi, A. N.
peraturerangeexplainsthedifficultyinobservingthebcc Simakov et al.,Phys. Plasmas 23, 056310 (2015).
phase experimentally. [23] G. Robert and A. Sollier, J. Phys. IV France 134, 257
(2006).
This work was supported by NSFC under Grant Nos.
[24] G.Robert,P.Legrand,andS.Bernard,Phys.Rev.B82,
U1530401, 41474069 and 11328401, and NSAF under
104118 (2010).
Grant No. U1530258. RW was supported by the [25] F. Luo, L.-C. Cai, X.-R. Chen, F.-Q. Jing, and D. Alf´e,
Hareaus visiting professorship award from the Univer- J. Appl. Phys. 111, 053503 (2012).
sity of Frankfurt and by NSF grants EAR-1319368 and [26] L. X. Benedict, T. Ogitsu, A. Trave, C. J. Wu, P.
-1348066. PZ was supported by the NSF grant DMR- A. Sterne, and E. Schwegler, Phys. Rev. B 79, 064106
(2009).
1506669. ComputationswereperformedatBeijingCom-
[27] G.V.Sin′koandN.A.Smirnov,Phys.Rev.B71,214108
putational Science Research Center and Minnesota Su-
(2005).
percomputing Institute.
[28] K. K´adas, L. Vitos, R. Ahuja, B. Johansson, and J.
Koll´ar, Phys. Rev. B 76, 235109 (2007).
[29] K. K´adas, L. Vitos, B. Johansson, and J. Koll´ar, Phys.
Rev. B 75, 035132 (2007).
[30] D.-B. Zhang, T. Sun, and R. M. Wentzcovitch, Phys.
∗ Corresponding author; [email protected] Rev. Lett. 112, 058501 (2014).
[1] J. F. Cannon, J. Phys. Chem. Ref. Data, 3, 781 (1974). [31] T. Sun, D.-B. Zhang, and R. M. Wentzcovitch, Phys.
[2] K. Persson, M. Ekman, V. Ozolinsˇ, Phys. Rev. B, 61, Rev. B 89, 094109 (2014).
11221 (2000). [32] T. Sun, X. Shen, and P. B. Allen, Phys. Rev. B 82,
[3] D. A. Young, Phase Diagrams of the Elements, Univer- 224304 (2010).
sity of California, 1991. [33] A.J.C.Ladd,B.Moran, andW.G.Hoover,Phys. Rev.
[4] A.J.Martin and A.Moore, J. Less-Common Met.1, 85 B 34, 5058 (1986).
(1959). [34] J.D.Althof,P.B.Allen,R..M.Wentzcovitch,andJ.A.
[5] A. Abey, LLNL Report No. UCRL-53567, 1984 (unpub- Moriarty, Phys. Rev. B 48, 13253 (1993).
lished). [35] S.Baroni,S.D.Gironcoli,A.D.Corso,andP.Giannozzi,
[6] M. Francois and M. Contre, in Proceedings of the Con- Rev. Mod. Phys. 73, 515 (2001); Giannozzi, et al., J.
ference internationale sur la metallurgie du beryllium, Phys. Condens. Matter 21, 395502 (2009).
Grenoble,(PressesUniversitairesdeFrance,Paris,1965). [36] P. E. Bl¨ochl, Phys. Rev.B 50, 17953 (1994).
[7] W. J. Evans, M. J. Lipp, H. Cynn, C. S. Yoo, M. So- [37] G.Kresse,J.Furthm¨oller,Phys.Rev.B54,11169(1999).
mayazulu,D.H¨ausermann,G.ShenandV.Prakapenka, [38] G. KresseandJ.Hafner,Phys. Rev. B49,14251 (1994).
¯
Phys. Rev. B 72, 094113 (2005). [39] J. P.Perdew, K.Burke,M. Ernzerhof, Phys. Rev. B 77,
[8] L.C.MingandM.H.Manghnani,J.Phys.F:Met.Phys. 3865 (1996).
6
[40] A.Togo and I. Tanaka, Scr. Mater., 108, 1 (2015). 5571 (1988).
[41] S. Nos´e, J. Chem. Phys. 81, 511 (1984); W.G. Hoover, [45] R. M. Wentzcovitch,Phys. Rev. B 50, 10358 (1994).
¯
Phys. Rev. A 31, 1695 (1985). [46] D.C.Wallace,ThermodynamicsofCrystals(Wiley,New
[42] A. A.Maradudin and A. E. Fein, Phys. Rev. 128, 2589 York, 1972).
(1962). [47] N. D. Mermin, Phys. Rev. 137, A1441 (1965).
[43] G.Grimvall,B.Magyari-Kop¨e, VidvudsOzolin¸ˇs, andK. [48] R.M.Wentzcovitch,J.L.Martins,andP.B.Allen,Phys.
A.Persson, Rev. Mod. Phys. 84, 945 (2012). Rev. B 45, 11372 (1992).
[44] R.M.WentzcovitchandM.L.Cohen,Phys. Rev. B,37,