[Pw_forum] Question: inconsistency for d3 between symmetry on and off?

Jiawei Zhou (周嘉炜) zhoujw20 at gmail.com
Tue Apr 22 19:41:58 CEST 2014


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://lists.quantum-espresso.org/pipermail/users/attachments/20140422/0924e91f/attachment.html>


More information about the users mailing list