Stress-strain curve by VASP

[jun-jie_shi1]

Hi,
we recently simulated the uniaxial stress-strain curve of metal copper by
changing the lattice constant in one direction while relaxing the lattice
constant in the other two directions. It turns out that the curve displays
some similarity to that of the copper extracted from the experimental
results in the reference. However, the magnitude differs by a factor of
10^4 (MPa compared to 10 GPa calculated by VASP).

VASP-Copper.png

[/size]

Could you kindly give us any advice on the discrepancies of the results, or
is VASP suitable for calculating stress-strain curve in this way?

[henrique_miranda]

I find it very, very suspicious to have a factor of 1e4 between experiment and calculation.
This type of elastic properties (in particular for a simple metal like Cu) should be fairly well reproduced using DFT.

Check if the deformation you are applying to the POSCAR corresponds to what you have on the x-axis and if you are reading the stress tensor from the correct place in the output and in the correct units.
I might be able to help you further if you share some or all of the POSCAR, INCAR, KPOINTS and OUTCAR files you use to perform the VASP calculations.

[jun-jie_shi1]

Dear Miranda,
Thanks for your reply. Here are the calculation details of simulating uniaxial stress-strain curve of copper.

1) After structural relaxation, the crystal constant of copper is 3.63 Å, which is nearly consistent with that of copper experimentally (3.61 Å).
POSCAR:
Cu4
1.00000000000000
3.6341213118432316 0.0000000000000000 -0.0000000000000000
0.0000000000000000 3.6341213118432316 -0.0000000000000000
-0.0000000000000000 -0.0000000000000000 3.6341213118432316
Cu
4
Direct
-0.0000000000000000 -0.0000000000000000 0.0000000000000000
0.0000000000000000 0.5000000000000000 0.5000000000000000
0.5000000000000000 0.0000000000000000 0.5000000000000000
0.5000000000000000 0.5000000000000000 0.0000000000000000

2) The calculation details of 50% strain of copper:
① INCAR:
SYSTEM = VASP
ISTART = 0
NWRITE = 2
PREC = Accurate
ENCUT = 500
GGA = PE
ADDGRID = .TRUE.
NSW = 200
ISIF = 3
ISYM = 2
IBRION = 2
NELM = 80
EDIFF = 1E-06
EDIFFG = -0.001
ALGO = Normal
LDIAG = .TRUE.
LREAL = Auto
ISMEAR = 0
SIGMA = 0.03
LORBIT = 11
ICHARG = 2
LWAVE = .FALSE.
LCHARG = .FALSE.
ICORELEVEL = 1
NCORE = 4

② OPTCELL：To fix x axial
000
010
001

③ KPOINTS：
Monkhorst Pack
0
Gamma
20 20 20
.0 .0 .0

④ After structural relaxation, the CONTCAR (x:+50 %):
Cu4
1.00000000000000
5.4511819677648470 0.0000000000000000 0.0000000000000000
0.0000000000000000 3.3884547989523495 0.0000000000000000
0.0000000000000000 0.0000000000000000 3.3884547989523495
Cu
4
Direct
0.0000000000000000 0.0000000000000000 -0.0000000000000000
-0.0000000000000000 0.5000000000000000 0.5000000000000000
0.5000000000000000 0.0000000000000000 0.5000000000000000
0.5000000000000000 0.5000000000000000 0.0000000000000000

⑤ OUTCAR：The stress -212.53 kBar was read from it.
FORCE on cell =-STRESS in cart. coord. units (eV):
Direction XX YY ZZ XY YZ ZX
--------------------------------------------------------------------------------------
Alpha Z 186.95545 186.95545 186.95545
Ewald -1794.47009 -1077.17259 -1077.17259 0.00000 0.00000 -0.00000
Hartree 263.88562 697.02291 697.02291 -0.00000 -0.00000 -0.00000
E(xc) -256.79638 -256.20166 -256.20166 0.00000 0.00000 -0.00000
Local 949.13637 -189.24371 -189.24371 0.00000 -0.00000 -0.00000
n-local -98.42952 -91.38961 -91.87284 -1.34843 5.26052 0.06349
augment 381.04408 374.97363 374.97364 0.00000 0.00000 0.00000
Kinetic 360.37217 361.07951 349.50803 -2.43231 1.34407 -0.57365
Fock 0.00000 0.00000 0.00000 0.00000 0.00000 0.00000
-------------------------------------------------------------------------------------
Total -8.30230 -0.00343 -0.00343 -0.00000 0.00000 0.00000
in kB -212.52747 -0.08772 -0.08772 -0.00000 0.00000 0.00000
external pressure = -70.90 kB Pullay stress = 0.00 kB

Could you kindly give us some advice on this puzzlement? Thank you sincerely for your time and assistance again.

[henrique_miranda]

Ok, I don't see anything fundamentally wrong with your calculation.

Did you try to compute the Young's modulus
https://en.wikipedia.org/wiki/Young%27s_modulus
from the stress-strain curve that you get from VASP and then with the experimental data?
They will also be off by 2 orders of magnitude, but then you way to double-check whether the data has the correct units by comparing with the value from a table (110 GPa from the Wikipedia page above)