[Pw_forum] Projection, convergence and imaginary phonon frequencies

[John Kendrick]

Hi,

I have been doing some calculations recently on molecular crystals. 
Because of our interest in the lowest phonon frequencies I have been 
concerned that translational invariance is maintained as far as 
possible.  To this end I have increased the accuracy of the optimisation 
and the scf convergence, to a point where I cannot see how to converge 
the calculation any tighter.  I noticed that the norm conserving 
pseudopotentials were having greater difficulties than the others.  In 
particular for one molecule I have seen significant differences between 
the different projection schemes that ensure translational invariance. 
Particularly the intensities in the infrared a very different.  Below is 
some output from a representative calculation.

Crystal Projection
#  mode   [cm-1]     [THz]      IR           Raman     depol
     1      0.00    0.0000    0.0000         0.0011    0.7500
     2      0.00    0.0000    0.0000         0.0043    0.7500
     3      0.00    0.0000    0.0000         0.0047    0.7500
     4     22.00    0.6596    0.0000         6.7060    0.7500
     5     56.93    1.7067    0.0033         0.1326    0.7500
     6     56.93    1.7067    0.0033         0.1326    0.7500
     7     73.28    2.1969    0.0000         9.8083    0.6361
     8     82.51    2.4737    0.0000         0.0000    0.7422
     9     83.26    2.4960    0.2348         0.0334    0.7500
    10     83.26    2.4960    0.2348         0.0334    0.7500

NOASR No Projection
#  mode   [cm-1]     [THz]      IR           Raman     depol
     1    -52.77   -1.5820    0.0016         0.7753    0.7500
     2     22.00    0.6596    0.0000         6.7060    0.7500
     3     44.67    1.3392    0.0147         0.4140    0.7500
     4     44.67    1.3392    0.0147         0.4140    0.7500
     5     72.20    2.1646    0.0211         0.0034    0.7500
     6     72.20    2.1646    0.0211         0.0034    0.7500
     7     73.28    2.1969    0.0000         9.8083    0.6361
     8     82.51    2.4737    0.0000         0.0000    0.7425
     9     83.79    2.5120    0.0000         0.1142    0.7500
    10     84.69    2.5390    0.2889         0.6525    0.7500

SIMPLE Projection
#  mode   [cm-1]     [THz]      IR           Raman     depol
     1     -1.31   -0.0394    0.0000         0.0120    0.7500
     2     -1.31   -0.0394    0.0000         0.0120    0.7500
     3     -0.92   -0.0277    0.0000         0.0010    0.7500
     4     39.67    1.1893    0.0001         0.4896    0.7500
     5     39.67    1.1893    0.0001         0.4896    0.7500
     6     51.30    1.5380    0.0000         2.9853    0.7120
     7     60.39    1.8104    0.0000         6.5707    0.7500
     8     61.02    1.8293    0.0040         2.8566    0.7500
     9     61.02    1.8293    0.0040         2.8566    0.7500
    10     61.94    1.8570    0.0000         0.0000    0.7499
    11     68.41    2.0507    0.0000         7.7313    0.7500
    12     76.98    2.3078    0.1318         0.0930    0.7500
    13     76.98    2.3078    0.1318         0.0930    0.7500


Should I be looking at
the energy cutoff - 80Ryd?
the kspace integration - 1 1 2?
The unit cell is approximately 9X9X6 angstrom
Other convergence thresholds?
There maybe some meshes used for integration that could still help 
improve the calculation or perhaps an approach to further tightening the 
optimisation. Or this maybe as good as it gets?  I was wondering if 
anyone had any experience they could offer.

Input to the optimisation step looks like this.
  &CONTROL
       restart_mode = 'restart',
       calculation = 'relax' ,
       outdir = '.' ,
       pseudo_dir = '/home/jkendric/Data/PseudoPotentials' ,
       disk_io = 'minimal' ,
       nstep = 200,
       etot_conv_thr = 1.0D-6,
       forc_conv_thr = 1.0D-5,
  /
  &SYSTEM
       ibrav = 0,
       celldm(1) = 1.8897262,
       nat = 58,
       ntyp = 4,
       ecutwfc = 80.D0 ,
       occupations = 'fixed' ,
  /
  &ELECTRONS
       conv_thr = 1.D-11 ,
       mixing_beta = 0.3D0 ,
  /
  &IONS
       ion_dynamics = 'bfgs'
       pot_extrapolation = 'second_order' ,
       wfc_extrapolation = 'second_order' ,
  /
CELL_PARAMETERS cubic
   9.3027 0.0 0.0
   0.0 9.3027 0.0
   0.0 0.0 6.6403
ATOMIC_SPECIES
    C   12.01070  C.blyp-mt.UPF
    N   14.00674  N.blyp-mt.UPF
    O   15.99940  O.blyp-mt.UPF
    H    1.00794  H.blyp-vbc.UPF
ATOMIC_POSITIONS crystal
C        0.500000000   0.500000000   0.000000000
C        0.623684469   0.560707616   0.129462381
H        0.717493110   0.584367933   0.037673510
H        0.590256782   0.657115895   0.209561631
.....
K_POINTS automatic
  1 1 2 0 0 0
Energies in the output file look like this
!    total energy              = -1055.84030152 Ry
!    total energy              = -1055.83327033 Ry
!    total energy              = -1055.84030145 Ry
!    total energy              = -1055.84030152 Ry
!    total energy              = -1055.84030153 Ry
!    total energy              = -1055.84030153 Ry
!    total energy              = -1055.84030154 Ry
!    total energy              = -1055.84030154 Ry
!    total energy              = -1055.84030154 Ry
      Total force =     0.000129     Total SCF correction =     0.000006



Input to the phonon calculations looks like this
Phonons at Gamma Point with raman
  &inputph
   tr2_ph=1.0d-18,
   epsil=.true.,
   lraman=.true.,
   outdir='.',
   fildyn='petn.dynG',
  /
0.0 0.0 0.0

Thanks.

John Kendrick
University of Bradford
Bradford BD13 1DP
UK