Inconsistent forces between rhombohedral and hexagonal unit cells

Top Page > Browsing

Inconsistent forces between rhombohedral and hexagonal unit cells

Date: 2023/08/22 09:25
Name: Hiroaki Tanaka <hiroaki-tanaka@issp.u-tokyo.ac.jp>: Dear developers,

I am trying to obtain a relaxed crystal structure of Bi2Te3 (R-3m) but encountered a problem in the force calculations.

For rhombohedral Bi2Te3, primitive (rhombohedral) and conventional (hexagonal) unit cells are defined.
The following are the unit vectors and atom positions used in my calculations.
The unit cell is from experimental data and the atom positions are optimized for the hexagonal unit cell.

[Rhombohedral]
Atoms.SpeciesAndCoordinates.Unit FRAC
<Atoms.SpeciesAndCoordinates
1 Bi 0.0998 0.0998 0.0998 7.5 7.5 45 0 45 0 0 off
2 Bi 0.9002 0.9002 0.9002 7.5 7.5 45 0 45 0 0 off
3 Te 0.5000 0.5000 0.5000 8.0 8.0 45 0 45 0 0 off
4 Te 0.7082 0.7082 0.7082 8.0 8.0 45 0 45 0 0 off
5 Te 0.2918 0.2918 0.2918 8.0 8.0 45 0 45 0 0 off
Atoms.SpeciesAndCoordinates>
Atoms.UnitVectors.Unit Ang
<Atoms.UnitVectors
0.0000 -2.5334 10.1533
2.1940 1.2667 10.1533
-2.1940 1.2667 10.1533
Atoms.UnitVectors>

[Hexagonal]
Atoms.SpeciesAndCoordinates.Unit FRAC
<Atoms.SpeciesAndCoordinates
1 Te 0.6667 0.3333 0.0415 8.0 8.0 45 0 45 0 0 off
2 Bi 0.0000 0.0000 0.0998 7.5 7.5 45 0 45 0 0 off
3 Te 0.3333 0.6667 0.1667 8.0 8.0 45 0 45 0 0 off
4 Bi 0.6667 0.3333 0.2335 7.5 7.5 45 0 45 0 0 off
5 Te 0.0000 0.0000 0.2918 8.0 8.0 45 0 45 0 0 off
6 Te 0.3333 0.6667 0.3749 8.0 8.0 45 0 45 0 0 off
7 Bi 0.6667 0.3333 0.4331 7.5 7.5 45 0 45 0 0 off
8 Te 0.0000 0.0000 0.5000 8.0 8.0 45 0 45 0 0 off
9 Bi 0.3333 0.6667 0.5669 7.5 7.5 45 0 45 0 0 off
10 Te 0.6667 0.3333 0.6251 8.0 8.0 45 0 45 0 0 off
11 Te 0.0000 0.0000 0.7082 8.0 8.0 45 0 45 0 0 off
12 Bi 0.3333 0.6667 0.7665 7.5 7.5 45 0 45 0 0 off
13 Te 0.6667 0.3333 0.8333 8.0 8.0 45 0 45 0 0 off
14 Bi 0.0000 0.0000 0.9002 7.5 7.5 45 0 45 0 0 off
15 Te 0.3333 0.6667 0.9585 8.0 8.0 45 0 45 0 0 off
Atoms.SpeciesAndCoordinates>
Atoms.UnitVectors.Unit Ang
<Atoms.UnitVectors
4.3880 0.0000 0.0000
-2.1940 3.8001 0.0000
0.0000 0.0000 30.4600
Atoms.UnitVectors>

After the SCF calculations, I got the following force tables from the output file.

[Rhombohedral]
***********************************************************
***********************************************************
xyz-coordinates (Ang) and forces (Hartree/Bohr)
***********************************************************
***********************************************************

<coordinates.forces
5
1 Bi 0.00000 0.00000 3.03990 -0.000000005821 -0.000152357335 0.024998932515
2 Bi -0.00000 -0.00000 27.42000 0.000000005800 0.000152357326 -0.024998698278
3 Te 0.00000 0.00000 15.22995 0.000000000031 0.000000000007 0.000000039453
4 Te -0.00000 -0.00000 21.57170 0.000000008205 0.000060426458 0.044143587947
5 Te 0.00000 -0.00000 8.88820 -0.000000008207 -0.000060426459 -0.044143787574
coordinates.forces>

[Hexagonal]
***********************************************************
***********************************************************
xyz-coordinates (Ang) and forces (Hartree/Bohr)
***********************************************************
***********************************************************

<coordinates.forces
15
1 Te 2.19422 1.26657 1.26409 -0.000444575515 0.000316880728 0.000796635259
2 Bi 0.00000 0.00000 3.03991 0.000092425753 -0.000206437842 -0.000302517604
3 Te -0.00022 2.53353 5.07768 0.000200511268 -0.000113816020 -0.000698647098
4 Bi 2.19422 1.26657 7.11241 -0.000071626983 0.000192938792 0.000356141423
5 Te 0.00000 0.00000 8.88823 0.000154993103 -0.000149936755 0.000287501055
6 Te -0.00022 2.53353 11.41945 0.000082117010 0.000012512273 -0.000297319214
7 Bi 2.19422 1.26657 13.19223 0.000076283439 -0.000196720739 0.000351584429
8 Te 0.00000 0.00000 15.23000 -0.000000000014 0.000000000002 -0.000000002536
9 Bi -0.00022 2.53353 17.26777 -0.000076283410 0.000196720734 -0.000351613081
10 Te 2.19422 1.26657 19.04055 -0.000082117044 -0.000012512263 0.000297318003
11 Te 0.00000 0.00000 21.57177 -0.000154993133 0.000149936761 -0.000287479267
12 Bi -0.00022 2.53353 23.34759 0.000071626979 -0.000192938776 -0.000356047133
13 Te 2.19422 1.26657 25.38232 -0.000200511282 0.000113816014 0.000698649600
14 Bi 0.00000 0.00000 27.42009 -0.000092425732 0.000206437833 0.000302451483
15 Te -0.00022 2.53353 29.19591 0.000444575549 -0.000316880751 -0.000796655455
coordinates.forces>

Forces obtained for the hexagonal unit cell seem to have no problem.
They are small enough because atom positions are optimized beforehand.
However, forces for the rhombohedral unit cell exhibit strangely large values (0.02 or 0.04 Hartree/Bohr).

I could not understand why such a difference occurs only by changing the unit cell.
I would appreciate it if you identify the reason.

Best regards,
Hiroaki Tanaka

Page: [1]

Re: Inconsistent forces between rhombohedral and hexagonal unit cells ( No.1 )

Date: 2023/08/22 14:23
Name: Naoya Yamaguchi

Hi,

Although the atomic positions are found to be shifted by comparing one structure with the other (e.g. `3.03990` and `3.03991`), issues like this often occur when the given unit cell is slender, and it may be caused by numerical errors coming from description of the density gradient.

In my experiences, there are several workarounds.

1. Avoid `GGA-PBE` for `scf.XcType`

2. Use a denser real space grid through `scf.Ngrid` instead of `scf.energycutoff`

As an example of `scf.Ngrid` in a slender cell, you can refer to No. 8 of https://www.openmx-square.org/forum/patio.cgi?mode=view&no=2649 .

3. Not to use the slender cell

Regards,
Naoya Yamaguchi

Re: Inconsistent forces between rhombohedral and hexagonal unit cells ( No.2 )

Date: 2023/08/22 15:05
Name: Hiroaki Tanaka

Dear Naoya Yamaguchi,

Thank you for your prompt reply.
I follow your advices and find that denser Ngrid works well.

Best regards,
Hiroaki Tanaka

Page: [1]