a problem in calculation elastic constants

Phonons, DFPT, electron-phonon, electric-field response, mechanical response…

Moderators: mverstra, joaocarloscabreu

Locked
nannan
Posts: 10
Joined: Thu Mar 21, 2013 10:09 am

a problem in calculation elastic constants

Post by nannan » Sat Sep 21, 2013 5:47 am

Hello there,
I have been confused by a question regarding elastic calculation for a long time. The case is that my input file was as following,
#Response function calculation for:
# * rigid-atom elastic tensor
# * rigid-atom piezoelectric tensor
# * interatomic force constants at gamma
# * Born effective charges
ndtset 3
# Set 1 : Initial self-consistent run

iscf1 7
kptopt1 1
tolvrs1 1.0d-15 #need excellent convergence of GS quantities for RF runs

# Set 2 : Calculate the ddk wf's - needed for piezoelectric tensor and
# Born effective charges in dataset 3

getwfk2 -1
iscf2 -3 #this option is needed for ddk
kptopt2 2 #use time-reversal symmetry only for k points
nqpt2 1 #one wave vector will be specified
qpt2 0 0 0 #need to specify gamma point
rfelfd2 2 #set for ddk wf's only
rfdir2 1 1 1 #full set of directions needed
tolwfr2 1.0d-17 #only wf convergence can be monitored here

# Set 3 : response-function calculations for all needed perturbations

getddk3 -1
getwfk3 -2
iscf3 7
kptopt3 2 #use time-reversal symmetry only for k points
nqpt3 1
qpt3 0 0 0
rfphon3 1 #do atomic displacement perturbation
rfatpol3 1 11 #do for all atoms
rfstrs3 3 #do strain perturbation
rfdir3 1 1 1 #the full set of directions is needed
tolvrs3 1.0d-10 #need reasonable convergence of 1st-order quantities

#Common input data


# acell COPY RELAXED RESULT FROM PREVIOUS CALCULATION
# Here is a default value, for automatic testing : suppress it and fill the previous line
acell 9.1919708857E+00 9.1919708857E+00 9.1919708857E+00
rprim 7.7815978695E-01 4.4927076247E-01 4.3888851428E-01
-7.7815978695E-01 4.4927076247E-01 4.3888851428E-01
-2.2013175300E-16 -8.9854152494E-01 4.3888851428E-01

#Definition of the atom types and atoms
ntypat 3
znucl 8 15 31
natom 11
typat 1 1 1 1 1 1 1 2 3 3 3

#Starting approximation for atomic positions in REDUCED coordinates
#based on ideal tetrahedral bond angles

# xred COPY RELAXED RESULT FROM PREVIOUS CALCULATION
# Here is a set of default values, for automatic testing : suppress it and fill the previous line
xred 5.3817283486E-01 1.5849785681E-01 5.3817283486E-01
-4.4075828733E-03 -4.4075828733E-03 -4.4075828733E-03
2.4341909448E-01 2.4341909448E-01 2.4341909448E-01
3.8914047774E-01 7.7294502161E-01 3.8914047774E-01
8.2564542231E-01 1.5908814016E-01 8.2564542231E-01
1.5849785681E-01 5.3817283486E-01 5.3817283486E-01
7.7294502161E-01 3.8914047774E-01 3.8914047774E-01
1.5908814016E-01 8.2564542231E-01 8.2564542231E-01
5.3817283486E-01 5.3817283486E-01 1.5849785681E-01
3.8914047774E-01 3.8914047774E-01 7.7294502161E-01
8.2564542231E-01 8.2564542231E-01 1.5908814016E-01

#Gives the number of bands, explicitely (do not take the default)
nband 55 # For an insulator (if described correctly as an
# insulator by DFT), conduction bands should not
# be included in response-function calculations

#Definition of the plane wave basis set
ecut 35.0 # Maximum kinetic energy cutoff (Hartree)
ecutsm 25.0 # Smoothing energy needed for lattice paramete
# optimization. This will be retained for
# consistency throughout.

#Definition of the k-point grid
kptopt 1 # Use symmetry and treat only inequivalent points
ngkpt 4 4 4 # 4x4x4 Monkhorst-Pack grid
nshiftk 1 # Use one copy of grid only (default)
shiftk 0.5 0.5 0.5 # This choice of origin for the k point grid
# preserves the hexagonal symmetry of the grid,
# which would be broken by the default choice.

#Definition of the self-consistency procedure
diemac 9.0 # Model dielectric preconditioner
iscf 7 # Use conjugate-gradient SCF cycle
nstep 1000 # Maxiumum number of SCF iterations

the value of ecut is about 35.0 Hatree, and the kptopt is convergent at 6×6×6. When I optimize the geometry structure of crystal, the less value of ecut, the more exact to the crystal parameter. But if the value of ecut less than 28.0 Hatree, the elastic matrix is as simple cubic not that of trigallium structure. If ecut is more than 30.0 Hatree, the elastic matrix is as the trigallium structure which has 6 independent elastic constants. But they all higher than the experiment value, even is twice to that.
For example,C11 is 167, but this work is 329.07266, 335.058, 329.07761. this trend is same as others.
The warning in ananddboutput file is as following,
Acoustic sum rule violation met : the eigenvalues of acoustic mode are too large at Gamma point.
Increase cutoff energy or k-point sampling.
I have tried to increase ecut or kptopt, but elastic values nearly never changed. I cannot find the which varies I should change, so please give a hand to guide how can I do ? Thanks!
nannan

ilukacevic
Posts: 271
Joined: Sat Jan 16, 2010 12:05 pm
Location: Dept. of Physics, University J. J. Strossmayer, Osijek, Croatia
Contact:

Re: a problem in calculation elastic constants

Post by ilukacevic » Mon Sep 23, 2013 8:19 am

Hi!

1] Have in mind that all parameters in your relaxation and elastic input files have to be the same (except the ones corresponding to DFPT).

2] Is your relaxation converged? Which pressure do you get in the end of the 1st dataset of your elastic output?

3] If 1] and 2] are ok, it could be your pseudopotential.

4] Your ecutsm is huge. Are you sure that you need that large value for your relaxation? A common value is around 0.5 Ha.

5] I don't understand your sentence:
nannan wrote:When I optimize the geometry structure of crystal, the less value of ecut, the more exact to the crystal parameter.

Please, clarify.

Hope this helps,

Igor

Locked