Dear all,
I was planning to use d3 to do Gamma point linewidth calculation. I start with
silicon, which is an example provided in the phonon package. I successively run
pw, ph and d3 and then I can get the d3 dynamical matrix. The input files and
the final output from d3 are shown below: (the version I use is 5.0.3)
pw.x:
&control
calculation= 'scf'
wf_collect=.true.
restart_mode= 'from_scratch'
pseudo_dir= '/home/jiawei/pseudo/'
outdir= '/scratch/jiawei/silicon/'
prefix= 'Si-ph-k7-e20-a10.21-pp:pz-vbc'
tprnfor= .true.
tstress= .true.
/
&system
ibrav= 2, celldm(1)= 10.21, nat= 2, ntyp= 1
ecutwfc=20
/
&electrons
conv_thr= 1.0d-12
diagonalization= 'david'
mixing_beta= 0.7
/
ATOMIC_SPECIES
Si 28.0855 Si.pz-vbc.UPF
ATOMIC_POSITIONS
Si 0.00 0.00 0.00
Si 0.25 0.25 0.25
K_POINTS {automatic}
7 7 7 0 0 0
ph.x:
phonons of Si at Gamma
&inputph
recover=.false.
tr2_ph=1.0d-14
prefix= 'Si-ph-k7-e20-a10.21-pp:pz-vbc'
epsil=.false.
ldisp=.false.
amass(1)=28.0855
outdir='/scratch/jiawei/silicon/'
fildyn='Si.dyn'
fildrho='Si.drhoG'
/
0.0 0.0 0.0
d3.x:
anharm at BZ
&inputph
prefix= 'Si-ph-k7-e20-a10.21-pp:pz-vbc'
amass(1)=28.0855
outdir='/scratch/jiawei/silicon/'
fildyn='Si.anhG'
fildrho='Si.drhoG'
fild0rho='Si.drhoG'
wraux=.true.
/
0.000 0.000 0.000
Si.anhG:
Derivative of the force constants
1 2 2 10.2100000 0.0000000 0.0000000 0.0000000 0.0000000 0.0000000
1 'Si ' 28.085500000000000
1 1 0.0000000 0.0000000 0.0000000
2 1 0.2500000 0.2500000 0.2500000
Third derivative in cartesian axes
q = ( 0.000000000 0.000000000 0.000000000 )
modo: 1
1 1
-0.416333634234E-16 0.000000000000E+00 -0.192592994439E-33
0.000000000000E+00 0.192592994439E-33 0.000000000000E+00
-0.416333634234E-16 0.000000000000E+00 0.192592994439E-33
0.000000000000E+00 0.385847337244E+00 0.000000000000E+00
0.277555756156E-16 0.000000000000E+00 0.385847337244E+00
0.000000000000E+00 0.192592994439E-33 0.000000000000E+00
......
(here I just show the first several numbers)
Up to now there is no problem. But I also wanted to know the wavefunction on
the whole Brillouin zone. Because pw.x will automatically find the IBZ and only
calculate wavefunction on that reduced k-mesh, I want to turn off the symmetry
and see what I can get. There are two switches I found that I can use:
nosym_evc and noinv. So I set them both to be true: (inside the &system control
flag)
nosym_evc=.true.
noinv=.true.
So the pw.x will not use any symmetry and do the calculation in the FBZ. The
result from pw.x clearly shows no symmetry found and 343(7*7*7) k-points are
used. The results from phonon code are also the same.(dynamical matrix at Gamma
point as well as the frequency). However, when it comes to d3, the results are
different. It runs without any problem, but the final d3.anhG file I get is
below:
d3.anhG (symmetry turnned off)
modo: 1
1 1
0.000000000000E+00 0.000000000000E+00 0.832667268469E-16
0.000000000000E+00 0.832667268469E-16 0.000000000000E+00
0.277555756156E-16 0.000000000000E+00 -0.277555756156E-16
0.000000000000E+00 -0.747798812221E+00 0.000000000000E+00
0.111022302463E-15 0.000000000000E+00 -0.747798812221E+00
0.000000000000E+00 0.832667268469E-16 0.000000000000E+00
1 2
0.000000000000E+00 0.000000000000E+00 0.111022302463E-15
0.000000000000E+00 -0.111022302463E-15 0.000000000000E+00
0.111022302463E-15 0.000000000000E+00 0.111022302463E-15
0.000000000000E+00 -0.177369122767E+01 0.000000000000E+00
-0.111022302463E-15 0.000000000000E+00 -0.177369122767E+01
0.000000000000E+00 -0.111022302463E-15 0.000000000000E+00
.......
The results are different from the previous one. I wonder why this happens?
Does the symmetry helps to constrain the d3 dynamical matrix elements, so that
it cannot be simply neglected like that is done in pw and ph?
The other question goes into more details. I tried to see what effects the
symmetry will have on wavefunctions. Originally I thought as long as we get the
ground-state wavefunction for different bands on a k-mesh from pw.x, we can
know any of its value on any k, including k+q.(these are used by phonon I
think). At gamma point, it should be no different. I noticed that phq_init
subroutine read wavefunctions (from files to evc and evq, each corresponding to
point k and k+q if q is not 0), so I output them to an extra file to see how
they look like. For this I use 4.2.1 version (because I cannot change the
already installed QE on the cluster). When symmetry is active, I get:
q = ( 0.000000000000000E+000 , 0.000000000000000E+000 ,
0.000000000000000E+000 )
nksq = 20
************************************************
npol = 1
npwx = 411
nbnd = 4
ik # = 1 -----------------------------------
0.69783903837624816635E+00 -0.64649301435265482674E+00
0.69783903837624816635E+00 -0.64649301435265482674E+00
-0.19216558783363027322E-16 0.87586028340881804830E-16
-0.19216558783363027322E-16 0.87586028340881804830E-16
0.55275566596732210915E-16 0.32076646606791784316E-17
0.55275566596732210915E-16 0.32076646606791784316E-17
0.23060842986941295573E-17 -0.42753486774548910442E-16
0.23060842986941295573E-17 -0.42753486774548910442E-16
0.99051076923549141728E-01 0.37832014530605620820E-02
0.99051076923549141728E-01 0.37832014530605620820E-02
......(these shows the first several numbers for G at k=0,q=0,ibnd=1,ig=1..ngm)
When symmetry is turnned off, I get:
q = ( 0.000000000000000E+000 , 0.000000000000000E+000 ,
0.000000000000000E+000 )
nksq = 343
************************************************
npol = 1
npwx = 411
nbnd = 4
ik # = 1 -----------------------------------
-0.49664369043450279362E+00 0.81134307529730043118E+00
-0.49664369043450279362E+00 0.81134307529730043118E+00
-0.27818613173485040276E-10 0.54661997423613342026E-10
-0.27818613173485040276E-10 0.54661997423613342026E-10
0.43853898353156132494E-10 -0.42976651079895819979E-10
0.43853898353156132494E-10 -0.42976651079895819979E-10
0.63231205525296367916E-11 -0.10330070651313073490E-10
0.63231205525296367916E-11 -0.10330070651313073490E-10
-0.96373137459236146718E-01 0.23187212493507927680E-01
-0.96373137459236146718E-01 0.23187212493507927680E-01
.......(at k=0,q=0,ibnd=1,ig=1..ngm)
Which confuses me is that they look quite different, because I orginally
thought they should be the same as long as we are looking at the same point. I
thought maybe this has something to do with the d3 code. But I have no idea.
Besides, although phonon code output different wavefunctions at the same
k-point, they give the same final results (dynamical matrix and frequencies).
So I am wondering how phonon do the sum over k-mesh and why the wavefunctions
would be different for the same point and band.
I will appreciate a lot for any of your suggestions. Thanks!
Regards,
Jiawei
Massachusetts Institute of Technology
---------------------------------
Jiawei Zhou, zhoujw20 at gmail.com
2014/4/22
-------------- next part --------------
An HTML attachment was scrubbed...
URL:
http://pwscf.org/pipermail/pw_forum/attachments/20140422/0924e91f/attachment.html
