Recently when I was trying to run a self-consistent calculation of a hcp structure, I first relaxed my hexagonal unit-cell (with 2 atoms), and then determined the optimized density and c/a ratio under certain pressure.
After that I switched to an orthorhombic hcp unit cell (with same density and c/a ratio), containing 4 atoms, and expanded it to a supercell with 144 atoms.
However, when I made a self-consistent calculation of this expanded orthorhombic supercell, there are large forces exerted on y direction for each atom, roughly 0.01 eV/Angstrom. While the forces on x and z direction are all negligible, which was unusual.
By comparison, I expanded my hexagonal unit cell(2 atoms ) to a supercell with 96 atoms, and then made the same self-consistent calculation. To my surprise the forces on y-direction almost disappeared.
These two hcp supercells(orthorhombic and hexagonal) shared exactly same density and c/a ratio.
Could it be possible that something went wrong with VASP compilation so that force calculation is not accurate for orthorhombic structure?
We would really appreciate it if anyone could offer help.
Sincerely,
Jizhou
force calculation of VASP
Moderators: Global Moderator, Moderator
-
martin.schlipf
- Global Moderator
- Posts: 542
- Joined: Fri Nov 08, 2019 7:18 am
Re: force calculation of VASP
The most common mistake here is a slight change in the structure. I don't know if you generated the supercell manually or used a tool for it, but what I can recommend is visualizing the structure with a tool that allows you to measure distances (e.g. VESTA). Then verify that all bonds are exactly the same length.
Alternatively, depending on the method used, there may be cases where forces can appear in a lower symmetry structures. In that case the distorted structure would actually have a lower energy than the undistorted one, but because of symmetrization, you never explore that structure in the primitive cell. You can manually distort the hexagonal supercell and see if the atoms go back to the original position.
Alternatively, depending on the method used, there may be cases where forces can appear in a lower symmetry structures. In that case the distorted structure would actually have a lower energy than the undistorted one, but because of symmetrization, you never explore that structure in the primitive cell. You can manually distort the hexagonal supercell and see if the atoms go back to the original position.
Martin Schlipf
VASP developer
-
jizhou_wu1
- Newbie
- Posts: 5
- Joined: Wed Dec 09, 2020 12:21 am
Re: force calculation of VASP
Dear Martin,
Thank you so much for your quick reply. Indeed I made scf calculation with symmetric supercell with 144 atoms, expanded from orthorhombic unit cell (with 4 atoms) of hcp structure. The structure belongs to symmetry group 194, which has been tested with FINDSYM website.
Here's my orthorhombic hcp unitcell.
However the OUTCAR file shows that each atom is subject to a force around 0.007 eV/Ang along y-direction, which is unusual. I used 9x9x9 Monkhorst-Pack k-point sampling, far enough to converge.
Thanks again for helping me.
Sincerely,
Jizhou
Thank you so much for your quick reply. Indeed I made scf calculation with symmetric supercell with 144 atoms, expanded from orthorhombic unit cell (with 4 atoms) of hcp structure. The structure belongs to symmetry group 194, which has been tested with FINDSYM website.
Here's my orthorhombic hcp unitcell.
Code: Select all
POSCAR_unitcell
0.908509279684186
1.956555212672139 0.000000000000000 0.000000000000000
0.000000000000000 3.388853031405246 0.000000000000000
0.000000000000000 0.000000000000000 3.156314790820479
Be
4
Direct
0.000000000000000 0.000000000000000 0.000000000000000
0.500000000000000 0.500000000000000 0.000000000000000
0.000000000000000 0.333333333333333 0.500000000000000
0.500000000000000 0.833333333333333 0.500000000000000
Code: Select all
POSITION TOTAL-FORCE (eV/Angst)
-----------------------------------------------------------------------------------
0.00000 0.00000 0.00000 0.000000 0.007090 0.000000
0.00000 0.00000 2.86754 0.000000 0.007090 0.000000
0.00000 0.00000 5.73508 -0.000000 0.007090 0.000000
0.00000 3.07880 0.00000 0.000000 0.007090 0.000000
0.00000 3.07880 2.86754 -0.000000 0.007090 -0.000000
0.00000 3.07880 5.73508 0.000000 0.007090 -0.000000
0.00000 6.15761 0.00000 -0.000000 0.007090 -0.000000
0.00000 6.15761 2.86754 -0.000000 0.007090 -0.000000
0.00000 6.15761 5.73508 -0.000000 0.007090 -0.000000
1.77755 0.00000 0.00000 0.000000 0.007090 -0.000000
1.77755 0.00000 2.86754 0.000000 0.007090 -0.000000
1.77755 0.00000 5.73508 0.000000 0.007090 0.000000
1.77755 3.07880 0.00000 -0.000000 0.007090 -0.000000
1.77755 3.07880 2.86754 -0.000000 0.007090 -0.000000
1.77755 3.07880 5.73508 -0.000000 0.007090 -0.000000
1.77755 6.15761 0.00000 0.000000 0.007090 0.000000
1.77755 6.15761 2.86754 0.000000 0.007090 0.000000
1.77755 6.15761 5.73508 0.000000 0.007090 0.000000
3.55510 0.00000 0.00000 0.000000 0.007090 -0.000000
3.55510 0.00000 2.86754 0.000000 0.007090 -0.000000
3.55510 0.00000 5.73508 0.000000 0.007090 -0.000000
3.55510 3.07880 0.00000 0.000000 0.007090 0.000000
3.55510 3.07880 2.86754 -0.000000 0.007090 -0.000000
3.55510 3.07880 5.73508 0.000000 0.007090 -0.000000
3.55510 6.15761 0.00000 0.000000 0.007090 -0.000000
3.55510 6.15761 2.86754 0.000000 0.007090 0.000000
3.55510 6.15761 5.73508 0.000000 0.007090 0.000000
5.33265 0.00000 0.00000 0.000000 0.007090 0.000000
5.33265 0.00000 2.86754 0.000000 0.007090 0.000000
5.33265 0.00000 5.73508 0.000000 0.007090 0.000000
5.33265 3.07880 0.00000 0.000000 0.007090 0.000000
5.33265 3.07880 2.86754 0.000000 0.007090 -0.000000
5.33265 3.07880 5.73508 0.000000 0.007090 0.000000
5.33265 6.15761 0.00000 -0.000000 0.007090 -0.000000
5.33265 6.15761 2.86754 -0.000000 0.007090 -0.000000
5.33265 6.15761 5.73508 -0.000000 0.007090 -0.000000
0.88877 1.53940 0.00000 0.000000 0.007090 0.000000
0.88877 1.53940 2.86754 0.000000 0.007090 0.000000
0.88877 1.53940 5.73508 -0.000000 0.007090 0.000000
0.88877 4.61821 0.00000 -0.000000 0.007090 -0.000000
0.88877 4.61821 2.86754 -0.000000 0.007090 -0.000000
0.88877 4.61821 5.73508 -0.000000 0.007090 -0.000000
0.88877 7.69701 0.00000 0.000000 0.007090 0.000000
0.88877 7.69701 2.86754 0.000000 0.007090 0.000000
0.88877 7.69701 5.73508 0.000000 0.007090 -0.000000
2.66632 1.53940 0.00000 -0.000000 0.007090 0.000000
2.66632 1.53940 2.86754 -0.000000 0.007090 0.000000
2.66632 1.53940 5.73508 -0.000000 0.007090 -0.000000
2.66632 4.61821 0.00000 -0.000000 0.007090 -0.000000
2.66632 4.61821 2.86754 -0.000000 0.007090 -0.000000
2.66632 4.61821 5.73508 -0.000000 0.007090 0.000000
2.66632 7.69701 0.00000 0.000000 0.007090 -0.000000
2.66632 7.69701 2.86754 0.000000 0.007090 -0.000000
2.66632 7.69701 5.73508 0.000000 0.007090 0.000000
4.44387 1.53940 0.00000 0.000000 0.007090 -0.000000
4.44387 1.53940 2.86754 0.000000 0.007090 0.000000
4.44387 1.53940 5.73508 0.000000 0.007090 -0.000000
4.44387 4.61821 0.00000 0.000000 0.007090 -0.000000
4.44387 4.61821 2.86754 0.000000 0.007090 0.000000
4.44387 4.61821 5.73508 0.000000 0.007090 0.000000
4.44387 7.69701 0.00000 0.000000 0.007090 0.000000
4.44387 7.69701 2.86754 0.000000 0.007090 -0.000000
4.44387 7.69701 5.73508 0.000000 0.007090 -0.000000
6.22142 1.53940 0.00000 0.000000 0.007090 0.000000
6.22142 1.53940 2.86754 0.000000 0.007090 0.000000
6.22142 1.53940 5.73508 0.000000 0.007090 0.000000
6.22142 4.61821 0.00000 -0.000000 0.007090 -0.000000
6.22142 4.61821 2.86754 -0.000000 0.007090 0.000000
6.22142 4.61821 5.73508 -0.000000 0.007090 0.000000
6.22142 7.69701 0.00000 0.000000 0.007090 0.000000
6.22142 7.69701 2.86754 0.000000 0.007090 0.000000
6.22142 7.69701 5.73508 0.000000 0.007090 0.000000
0.00000 2.05254 1.43377 -0.000000 -0.007090 -0.000000
0.00000 2.05254 4.30131 0.000000 -0.007090 -0.000000
0.00000 2.05254 7.16885 0.000000 -0.007090 -0.000000
0.00000 5.13134 1.43377 -0.000000 -0.007090 -0.000000
0.00000 5.13134 4.30131 -0.000000 -0.007090 -0.000000
0.00000 5.13134 7.16885 0.000000 -0.007090 -0.000000
0.00000 8.21015 1.43377 -0.000000 -0.007090 0.000000
0.00000 8.21015 4.30131 -0.000000 -0.007090 0.000000
0.00000 8.21015 7.16885 -0.000000 -0.007090 0.000000
1.77755 2.05254 1.43377 -0.000000 -0.007090 -0.000000
1.77755 2.05254 4.30131 -0.000000 -0.007090 0.000000
1.77755 2.05254 7.16885 -0.000000 -0.007090 0.000000
1.77755 5.13134 1.43377 0.000000 -0.007090 -0.000000
1.77755 5.13134 4.30131 0.000000 -0.007090 -0.000000
1.77755 5.13134 7.16885 -0.000000 -0.007090 0.000000
1.77755 8.21015 1.43377 0.000000 -0.007090 0.000000
1.77755 8.21015 4.30131 0.000000 -0.007090 0.000000
1.77755 8.21015 7.16885 0.000000 -0.007090 -0.000000
3.55510 2.05254 1.43377 -0.000000 -0.007090 0.000000
3.55510 2.05254 4.30131 -0.000000 -0.007090 0.000000
3.55510 2.05254 7.16885 -0.000000 -0.007090 0.000000
3.55510 5.13134 1.43377 -0.000000 -0.007090 -0.000000
3.55510 5.13134 4.30131 -0.000000 -0.007090 0.000000
3.55510 5.13134 7.16885 -0.000000 -0.007090 0.000000
3.55510 8.21015 1.43377 -0.000000 -0.007090 -0.000000
3.55510 8.21015 4.30131 -0.000000 -0.007090 0.000000
3.55510 8.21015 7.16885 0.000000 -0.007090 -0.000000
5.33265 2.05254 1.43377 0.000000 -0.007090 -0.000000
5.33265 2.05254 4.30131 0.000000 -0.007090 0.000000
5.33265 2.05254 7.16885 -0.000000 -0.007090 0.000000
5.33265 5.13134 1.43377 -0.000000 -0.007090 -0.000000
5.33265 5.13134 4.30131 -0.000000 -0.007090 0.000000
5.33265 5.13134 7.16885 -0.000000 -0.007090 -0.000000
5.33265 8.21015 1.43377 0.000000 -0.007090 -0.000000
5.33265 8.21015 4.30131 0.000000 -0.007090 0.000000
5.33265 8.21015 7.16885 0.000000 -0.007090 0.000000
0.88877 0.51313 1.43377 -0.000000 -0.007090 -0.000000
0.88877 0.51313 4.30131 -0.000000 -0.007090 0.000000
0.88877 0.51313 7.16885 -0.000000 -0.007090 -0.000000
0.88877 3.59194 1.43377 -0.000000 -0.007090 0.000000
0.88877 3.59194 4.30131 -0.000000 -0.007090 -0.000000
0.88877 3.59194 7.16885 -0.000000 -0.007090 0.000000
0.88877 6.67074 1.43377 -0.000000 -0.007090 -0.000000
0.88877 6.67074 4.30131 -0.000000 -0.007090 -0.000000
0.88877 6.67074 7.16885 -0.000000 -0.007090 0.000000
2.66632 0.51313 1.43377 -0.000000 -0.007090 0.000000
2.66632 0.51313 4.30131 0.000000 -0.007090 0.000000
2.66632 0.51313 7.16885 0.000000 -0.007090 -0.000000
2.66632 3.59194 1.43377 -0.000000 -0.007090 0.000000
2.66632 3.59194 4.30131 -0.000000 -0.007090 0.000000
2.66632 3.59194 7.16885 -0.000000 -0.007090 -0.000000
2.66632 6.67074 1.43377 0.000000 -0.007090 0.000000
2.66632 6.67074 4.30131 0.000000 -0.007090 -0.000000
2.66632 6.67074 7.16885 0.000000 -0.007090 0.000000
4.44387 0.51313 1.43377 -0.000000 -0.007090 0.000000
4.44387 0.51313 4.30131 -0.000000 -0.007090 -0.000000
4.44387 0.51313 7.16885 -0.000000 -0.007090 0.000000
4.44387 3.59194 1.43377 -0.000000 -0.007090 -0.000000
4.44387 3.59194 4.30131 -0.000000 -0.007090 0.000000
4.44387 3.59194 7.16885 -0.000000 -0.007090 -0.000000
4.44387 6.67074 1.43377 0.000000 -0.007090 0.000000
4.44387 6.67074 4.30131 -0.000000 -0.007090 0.000000
4.44387 6.67074 7.16885 -0.000000 -0.007090 0.000000
6.22142 0.51313 1.43377 -0.000000 -0.007090 0.000000
6.22142 0.51313 4.30131 -0.000000 -0.007090 -0.000000
6.22142 0.51313 7.16885 -0.000000 -0.007090 -0.000000
6.22142 3.59194 1.43377 -0.000000 -0.007090 -0.000000
6.22142 3.59194 4.30131 -0.000000 -0.007090 -0.000000
6.22142 3.59194 7.16885 -0.000000 -0.007090 -0.000000
6.22142 6.67074 1.43377 0.000000 -0.007090 -0.000000
6.22142 6.67074 4.30131 0.000000 -0.007090 -0.000000
6.22142 6.67074 7.16885 -0.000000 -0.007090 0.000000
-----------------------------------------------------------------------------------
total drift: 0.000000 0.000136 -0.000000
--------------------------------------------------------------------------------------------------------
Sincerely,
Jizhou
-
martin.schlipf
- Global Moderator
- Posts: 542
- Joined: Fri Nov 08, 2019 7:18 am
Re: force calculation of VASP
Can you send me the full set of input files (INCAR, KPOINTS, POTCAR, POSCAR) for both the primitive and the extended system.
Did you try if this also happens if you create a smaller supercell?
Did you try if this also happens if you create a smaller supercell?
Martin Schlipf
VASP developer
-
jizhou_wu1
- Newbie
- Posts: 5
- Joined: Wed Dec 09, 2020 12:21 am
Re: force calculation of VASP
Dear Martin.
Thanks a lot. Indeed I performed self-consistent calculation on both orthorhombic unitcell and supercell. Large forces exist along y-direction for both supercell and primitive cell.
My INCAR file:
My unitcell POSCAR, an orthorhombic cell:
And my KPOINT made use of 15x15x15 Monkhorst-Pack Algorithm:
For POTCAR, I used PAW_PBE Be_sv 06Sep2000.
However in my output file, I found force around 5e-3 eV/Ang along y-direction:
Is it possible that I missed some high-symmetry points in Brillouin zone? Thanks for help.
Sincerely,
Jizhou
Thanks a lot. Indeed I performed self-consistent calculation on both orthorhombic unitcell and supercell. Large forces exist along y-direction for both supercell and primitive cell.
My INCAR file:
Code: Select all
SYSTEM = Be
NWRITE = 1
#Start parameter for this Run (values of parameters are equal to default):
ISTART = 0
ICHARG = 2
ENCUT = 1000 planewave cutoff
NELM = 500 max number of electronic steps
EDIFF = 1e-4
#EDIFFG = -.001 force stopping-criterion for geometry steps
NELMIN = 6
PREC = High
DOS related values:
ISMEAR = 0 -4-tet -1-fermi 1=Methfessel/Paxton 1.order
SIGMA = 0.1 broadening in eV
IALGO = 38
LREAL = .FALSE.
NELMIN = 6
BMIX = 2.0
MAXMIX = 50
POTIM = 0.3
ISIF = 2 (2:force=y stress=y ions=y shape=n volume=n)
IBRION = -1
ISYM = 0
NSW = 0
Code: Select all
POSCAR_HCP_Pri4atoms
0.908509279684186
1.956555212672139 0.000000000000000 0.000000000000000
0.000000000000000 3.388853031405246 0.000000000000000
0.000000000000000 0.000000000000000 3.156314790820479
Be
4
Direct
0.000000000000000 0.000000000000000 0.000000000000000
0.500000000000000 0.500000000000000 0.000000000000000
0.000000000000000 0.333333333333333 0.500000000000000
0.500000000000000 0.833333333333333 0.500000000000000
Code: Select all
Automatic mesh
0
M
15 15 15
0 0 0
However in my output file, I found force around 5e-3 eV/Ang along y-direction:
Code: Select all
POSITION TOTAL-FORCE (eV/Angst)
-----------------------------------------------------------------------------------
0.00000 0.00000 0.00000 0.000013 0.004999 0.000138
0.88877 1.53940 0.00000 -0.000021 0.004989 0.000152
0.00000 1.02627 1.43377 -0.000000 -0.004992 -0.000148
0.88877 2.56567 1.43377 0.000009 -0.004996 -0.000142
Sincerely,
Jizhou
-
martin.schlipf
- Global Moderator
- Posts: 542
- Joined: Fri Nov 08, 2019 7:18 am
Re: force calculation of VASP
Could you clarify what exactly is puzzling you? In the primitive cell you get forces of 5e-3 and in the supercell 7e-3. But you will have used a different k-point mesh for the supercell, so a bit of change in the forces is to be expected.
Also note that you mesh is not very well chosen. Typically you want roughly the same k-point density along all directions, i.e. for an orthorhombic cell (a, b, c) and a grid (NX, NY, NZ), you want that a * NX ~ b * NY ~ c * NZ
Also note that you mesh is not very well chosen. Typically you want roughly the same k-point density along all directions, i.e. for an orthorhombic cell (a, b, c) and a grid (NX, NY, NZ), you want that a * NX ~ b * NY ~ c * NZ
Martin Schlipf
VASP developer
-
jizhou_wu1
- Newbie
- Posts: 5
- Joined: Wed Dec 09, 2020 12:21 am
Re: force calculation of VASP
Dear Martin,
Thanks for your help.
What really puzzled me was the force along y-direction. Since it's a standard hcp cell (symmetry group 194), from my previous experience I would not expect such a large force 1e-3 along the y-direction even with many k-points.
Is there any setting I can toggle so that the force would be eliminated at these high-symmetry locations?
Thanks a lot.
Sincerely,
Jizhou
Thanks for your help.
What really puzzled me was the force along y-direction. Since it's a standard hcp cell (symmetry group 194), from my previous experience I would not expect such a large force 1e-3 along the y-direction even with many k-points.
Is there any setting I can toggle so that the force would be eliminated at these high-symmetry locations?
Thanks a lot.
Sincerely,
Jizhou
-
martin.schlipf
- Global Moderator
- Posts: 542
- Joined: Fri Nov 08, 2019 7:18 am
Re: force calculation of VASP
Could you clarify what you actually want to achieve? If you increase the precision of your calculation, the forces should probably get smaller. If you relax the structure, I'd expect the atoms to move only a tiny bit from their current position. If you don't want to relax the structure, why is it important whether the forces are zero?
Martin Schlipf
VASP developer
-
jizhou_wu1
- Newbie
- Posts: 5
- Joined: Wed Dec 09, 2020 12:21 am
Re: force calculation of VASP
Dear Martin,
Thank you so much for your instruction. Indeed I relaxed my supercell and then made a self-consistent calculation. The force is smaller than 1e-5 eV/Ang.
However, when I ran phonon calculation with my relaxed supercell using DFPT method(IBRION =
, the OUTCAR showed that the forces on certain atoms turned larger than 1e-2 eV/Ang.
Is it reasonable that force is different when I change IBRION from -1 (scf) to 8(DFPT) ? Because for bcc structure forces are always 0 in OUTCAR, unlike what is shown in hcp calculation.
Thank you again for your help.
Sincerely,
Jizhou
Thank you so much for your instruction. Indeed I relaxed my supercell and then made a self-consistent calculation. The force is smaller than 1e-5 eV/Ang.
However, when I ran phonon calculation with my relaxed supercell using DFPT method(IBRION =
, the OUTCAR showed that the forces on certain atoms turned larger than 1e-2 eV/Ang.Is it reasonable that force is different when I change IBRION from -1 (scf) to 8(DFPT) ? Because for bcc structure forces are always 0 in OUTCAR, unlike what is shown in hcp calculation.
Thank you again for your help.
Sincerely,
Jizhou
-
martin.schlipf
- Global Moderator
- Posts: 542
- Joined: Fri Nov 08, 2019 7:18 am
Re: force calculation of VASP
For phonon calculations it is important that you are in the energetic minimum, so the forces should really vanish. What you can always try is to compare the results from finite differences (IBRION = 5 or 6) with the ones from DFPT. If you find there are significant changes that may point to an error in the setup.
I'm not fully clear on the workflow that you use to calculate the phonon modes. Usually you perform an initial relaxation with IBRION = 1 or 2. Then you copy CONTCAR to POSCAR to calculate the phonon modes. Is this how you did it?
I'm not fully clear on the workflow that you use to calculate the phonon modes. Usually you perform an initial relaxation with IBRION = 1 or 2. Then you copy CONTCAR to POSCAR to calculate the phonon modes. Is this how you did it?
Martin Schlipf
VASP developer