<html><head><meta http-equiv="Content-Type" content="text/html; charset=utf-8"></head><body style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class=""><meta http-equiv="Content-Type" content="text/html; charset=utf-8" class=""><div style="word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class="">Dear all:<div class=""><br class=""></div><div class="">I know it is a question that has been asked so many times and I went over lots of them in the PW forum archive. But I still could not fix the problem.</div><div class=""><br class=""></div><div class="">I am calculating phonon dispersion of alpha-Uranium. First, the unit cell is relaxed to equilibrium followed by a SCF calculation. The input for the SCF calculation is copied below:</div><div class=""><div class=""><i class="">&control</i></div><div class=""><i class=""> calculation='scf',</i></div><div class=""><i class=""> restart_mode='from_scratch',</i></div><div class=""><i class=""> pseudo_dir = '/home/peng276/pseudopotential',</i></div><div class=""><i class=""> prefix='Uranium',</i></div><div class=""><i class=""> outdir='.',</i></div><div class=""><i class=""> forc_conv_thr =1E-5</i></div><div class=""><i class=""> etot_conv_thr = 1E-6</i></div><div class=""><i class=""> tstress=.true.</i></div><div class=""><i class=""> tprnfor=.true.</i></div><div class=""><i class=""> nstep=100</i></div><div class=""><i class=""> !verbosity="high"</i></div><div class=""><i class=""> /</i></div><div class=""><i class=""> &system</i></div><div class=""><i class=""> ibrav = 0,</i></div><div class=""><i class=""> celldm(1) = 5.4051</i></div><div class=""><i class=""> nat = 2,</i></div><div class=""><i class=""> ntyp = 1,</i></div><div class=""><i class=""> ecutwfc = 100,</i></div><div class=""><i class=""> occupations='smearing',</i></div><div class=""><i class=""> degauss=0.005</i></div><div class=""><i class=""> /</i></div><div class=""><i class=""> &electrons</i></div><div class=""><i class=""> conv_thr = 1.0e-8</i></div><div class=""><i class=""> electron_maxstep = 100,</i></div><div class=""><i class=""> !mixing_beta=0.3</i></div><div class=""><i class=""> diagonalization='cg'</i></div><div class=""><i class=""> /</i></div><div class=""><i class=""> &ions</i></div><div class=""><i class=""> ion_dynamics = 'bfgs'</i></div><div class=""><i class=""> /</i></div><div class=""><i class=""> &cell</i></div><div class=""><i class=""> cell_dynamics='bfgs',</i></div><div class=""><i class=""> press=0.0,</i></div><div class=""><i class="">/</i></div><div class=""><i class="">ATOMIC_SPECIES</i></div><div class=""><i class=""> U 238.029 U.pbe-spfn-rrkjus_psl.1.0.0.UPF</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class="">CELL_PARAMETERS (alat)</i></div><div class=""><i class=""> 0.490195165 -1.021315393 0.000000000</i></div><div class=""><i class=""> 0.490195165 1.021315393 0.000000000</i></div><div class=""><i class=""> 0.000000000 0.000000000 1.711961721</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class="">ATOMIC_POSITIONS (crystal)</i></div><div class=""><i class="">U 0.9017856414 0.0982143586 0.2500000000</i></div><div class=""><i class="">U 0.0982143586 0.9017856414 0.7500000000</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class="">K_POINTS {automatic}</i></div><div class=""><i class=""> 12 12 12 0 0 0</i></div></div><div class=""><br class=""></div><div class="">Then I start a phonon calculation on a qpoint mesh 4*4*4. The input ph.in is copied below:</div><div class=""><div class=""><i class="">phonons of alpha Uranium</i></div><div class=""><i class="">&inputph</i></div><div class=""><i class="">prefix = 'Uranium',</i></div><div class=""><i class="">amass(1) = 238.029,</i></div><div class=""><i class="">outdir = './',</i></div><div class=""><i class="">fildvscf='Udv',</i></div><div class=""><i class="">fildyn = 'Uranium.dyn',</i></div><div class=""><i class="">ldisp=.true.</i></div><div class=""><i class="">nq1=4,</i></div><div class=""><i class="">nq2=4,</i></div><div class=""><i class="">nq3=4,</i></div><div class=""><i class="">/</i></div></div><div class=""><br class=""></div>The obtained dispersion is plotted in the figure below, in which I see negative frequency of the first acoustic phonon mode near X. Then I tried to improve computational accuracy by tuning the calculation parameters used in both Pw.x and Ph.x, particularly looking at the phonon frequencies of one of the negative frequency phonons at (0.5,0,0). An example output is copied below:<div class=""><i class=""> Diagonalizing the dynamical matrix</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""> q = ( 0.500000000 0.000000000 0.000000000 ) </i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""> **************************************************************************</i></div><div class=""><i class=""> freq ( 1) = -0.939207 [THz] = -31.328581 [cm-1]</i></div><div class=""><i class=""> freq ( 2) = 1.029847 [THz] = 34.352008 [cm-1]</i></div><div class=""><i class=""> freq ( 3) = 2.289485 [THz] = 76.368998 [cm-1]</i></div><div class=""><i class=""> freq ( 4) = 2.490969 [THz] = 83.089795 [cm-1]</i></div><div class=""><i class=""> freq ( 5) = 2.862652 [THz] = 95.487797 [cm-1]</i></div><div class=""><i class=""> freq ( 6) = 3.345818 [THz] = 111.604479 [cm-1]</i></div><div class=""><i class=""> **************************************************************************.</i></div><div class="">The parameters I tried to tune include:</div><div class=""><div class=""><u class="">In SCF calculations:</u></div><div class="">1) varying Ecutwfc from 80 Ry to 140 Ry</div><div class="">2) varying kpoint grid from 8*8*8 to 16*16*16.</div><div class="">3) varying conv_thr from 1e-6 to 1e-10</div><div class="">4) varying degauss from 0.005 to 0.02 (since alpha-U is a metal, the gauss smearing might affect the energy and therefore force constant calculations)</div><div class=""><u class="">In Ph.x calculations:</u></div><div class="">1) varying tr2_ph from 1e-12 to 1e-14</div><div class="">2) varying alpha_mix(1) from 0.3 to 0.7 (In principle this should not matter too much since it is the proportion of old charge density used in updating to the new charge density, no matter what value is used at the physical quantities should be converged in the end)</div><div class=""><br class=""></div><div class="">However, in all these calculations with different combinations of the parameter sets, I am consistently getting negative frequencies at (0.5,0,0). The possible reasons:</div><div class="">1) Unstable structure. It is unlikely this is the cause, the frequency at Gamma even before I impose the acoustic summation rule is very close to 0. Also the forces obtained from SCF calculation on the equilibrium structure (see below) is quite small.</div><div class=""><div class=""><i class=""> Forces acting on atoms (cartesian axes, Ry/au):</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""> atom 1 type 1 force = 0.00000000 0.00001190 0.00000000</i></div><div class=""><i class=""> atom 2 type 1 force = -0.00000000 -0.00001190 0.00000000</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""> Total force = 0.000017 Total SCF correction = 0.000003</i></div><div class=""><i class=""> SCF correction compared to forces is large: reduce conv_thr to get better values</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""> Computing stress (Cartesian axis) and pressure</i></div><div class=""><i class=""><br class=""></i></div><div class=""><i class=""> total stress (Ry/bohr**3) (kbar) P= 0.05</i></div><div class=""><i class=""> 0.00000066 0.00000000 0.00000000 0.10 0.00 0.00</i></div><div class=""><i class=""> 0.00000000 -0.00000037 0.00000000 0.00 -0.05 0.00</i></div><div class=""><i class=""> 0.00000000 0.00000000 0.00000077 0.00 0.00 0.11</i></div><div class=""><i class=""><br class=""></i></div><div class="">2) Incommensurate qpoint grid with kpoint grid. I have seen some discussions on the PW forum that in metallic systems, the incommensurate qpoint and kpoint grid might lead to issues. However, I have tried using 12*12*12 and 16*16*16 kpoint grid in SCF calculations to do the phonon calculation on a 4*4*4 grid, but still get negative frequency at (0.5,0,0). Moreover, in the case of single qpoint phonon calculation, does it matter what kpoint mesh grid is used in the SCF calculation?</div><div class=""><br class=""></div><div class="">3) Unstable structure of the supercell displaced according to the phonon pattern. I think that something like soft modes could cause this negative frequency. However, the published dispersions of alpha-U all show similar shape to mine, with a sudden drop in the acoustic branch near X along the Gamma-X (see attached figure from a literature), but no negative frequency as in my case</div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">The figure is taken from: Manley, M. E., et al. "Phonon dispersion in uranium measured using inelastic x-ray scattering." </span><i style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">Physical Review B</i><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class=""> 67.5 (2003): 052302.</span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class=""><br class=""></span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">4) Spin-orbit coupling. In 5f electron systems, spin-orbital coupling might be important to include in calculations. However, there is still debate on whether it is necessary to include SOC in DFT calculations on alpha-U or not. Moreover, I have calculated dispersion of alpha-U in VASP without SOC and agrees with other published results well. So I would not go too deep into this, plus SOC calculations are expensive.</span></div><div class=""><br class=""></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class=""><br class=""></span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">With all the attempts and thoughts I have made, I doo not know how to resolve this issue. I strongly believe it is a numerical issue and currently trying another pseudopotential (Ulatrsolft PBE instead of the PAW PBE used in all above calculations). I would appreciate if anyone could give me suggestions or advice.</span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class=""><br class=""></span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">Thank you very much!</span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class=""><br class=""></span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">Best</span></div><div class=""><span style="color: rgb(34, 34, 34); font-family: Arial, sans-serif; font-size: 13px; font-variant-ligatures: normal; orphans: 2; widows: 2; background-color: rgb(255, 255, 255); text-decoration-thickness: initial;" class="">Jie</span></div><div class=""><img apple-inline="yes" id="5A6AEE99-E80F-48A6-B0B9-B8183914C4ED" class="" src="cid:0271E602-B8C2-4189-8207-2B26C6798901@itap.purdue.edu"><img apple-inline="yes" id="08BC97B1-2C8D-4341-8891-EF4EE97A8E3F" class="" src="cid:5421F4D3-5EE7-403D-A553-DFCE18F3CBC4@itap.purdue.edu"><div class=""><br class=""></div><div class=""><br class=""><div class="">
<div style="color: rgb(0, 0, 0); letter-spacing: normal; text-align: start; text-indent: 0px; text-transform: none; white-space: normal; word-spacing: 0px; -webkit-text-stroke-width: 0px; word-wrap: break-word; -webkit-nbsp-mode: space; line-break: after-white-space;" class="">————————————————————————————————————————————————<br class=""><div class="">Jie Peng, Postdoctoral research fellow<br class="">School of Materials Engineering, Purdue University, West Lafayette, Indiana<br class=""><a href="mailto:pengjiesjtu@gmail.com" class="">Email: pengjiesjtu@gmail.com</a><br class="">Cell: (240)-495-9445</div></div>
</div>
<br class=""></div></div></div></div></div></body></html>