jochym / Elastic

A module for ASE for elastic constants calculation.
GNU General Public License v3.0
35 stars 23 forks source link

Fix usage of `spglib` #70

Closed yuzie007 closed 5 months ago

yuzie007 commented 5 months ago

There is a bug in the present master introduced in b88c51c39b5350ad6492afcead861794908c3afd

For the following POSCAR file with the space group type of I4/m

generated by phonopy
   1.0
     5.7424952469569419    0.0000000000000000    0.0000000000000000
     0.0000000000000000    5.7424952469569419    0.0000000000000000
     0.0000000000000000    0.0000000000000000    3.5536083410560715
Ni W
   8    2
Direct
  0.5996722148177035  0.8006704283160834  0.0000000000000000
  0.3006704283160833  0.9003277851822965  0.5000000000000000
  0.6993295716839166  0.0996722148177036  0.5000000000000000
  0.4003277851822965  0.1993295716839166  0.0000000000000000
  0.0996722148177036  0.3006704283160833  0.5000000000000000
  0.8006704283160833  0.4003277851822964  0.0000000000000000
  0.1993295716839167  0.5996722148177036  0.0000000000000000
  0.9003277851822965  0.6993295716839166  0.5000000000000000
  0.0000000000000000  0.0000000000000000  0.0000000000000000
  0.5000000000000000  0.5000000000000000  0.5000000000000000

the following script

import ase.io
import spglib as spg
from elastic.elastic import get_lattice_type

atoms = ase.io.read("POSCAR")
cell = (atoms.cell, atoms.get_scaled_positions(), atoms.numbers)
print(spg.get_spacegroup(cell))
# print(spg.get_symmetry_dataset(cell))
print(get_lattice_type(atoms))

gives at this moment

I4/m (87)
(1, 'Triclinic', 'P1', 1)

meaning the Elastic does not detect the space group correctly.

The reason of this bug is because in spglib we must give "fractional" coordinates, while presently cartesian coordinates are given. (https://spglib.readthedocs.io/en/stable/python-interface.html#crystal-structure-cell)

The present PR fixes this. With this, the code above correctly gives

I4/m (87)
(4, 'Tetragonal', 'I4/m', 87)

I would also like to suggest using spg.get_symmetry_dataset rather than parsing the string of space group type and number on the Elastic side.

Thank you so much @jochym for your kind review in advance!

jochym commented 5 months ago

My mistake! Thanks, @yuzie007 ! Good idea to use get_symmetry_dataset. Merging.