giovannipizzi / seekpath

A module to obtain and visualize k-vector coefficients and obtain band paths in the Brillouin zone of crystal structures
Other
114 stars 45 forks source link

Test failure with spglib 1.12.2.post0 #57

Closed risicle closed 5 years ago

risicle commented 5 years ago

Hi, we're trying to maintain seekpath packages for Nix/NixOS. When we upgraded spglib from 1.10.4.11 to 1.12.2.post0 we noticed that a single test in seekpath 1.8.4, test_oI2Y started failing. The error being:

======================================================================
FAIL: test_oI2Y (seekpath.test_paths_hpkot.TestPaths3D_HPKOT)
----------------------------------------------------------------------
Traceback (most recent call last):
  File "/build/seekpath-1.8.4/seekpath/test_paths_hpkot.py", line 706, in test_oI2Y
    self.base_test(ext_bravais="oI2", with_inv=True)
  File "/build/seekpath-1.8.4/seekpath/test_paths_hpkot.py", line 366, in base_test
    self.assertEqual(res['bravais_lattice_extended'], ext_bravais)
AssertionError: 'oI3' != 'oI2'
- oI3
?   ^
+ oI2
?   ^

Same error on python 2.7 & 3.7, on macos & linux.

giovannipizzi commented 5 years ago

Thanks for the report! Indeed I can reproduce the error.

Last version that works: spglib==1.11.2.post1 First version that fails: spglib==1.12.0.post5

Minimal example to reproduce:

import numpy as np
import json

class NumpyEncoder(json.JSONEncoder):
    def default(self, obj):
        if isinstance(obj, np.ndarray):
            return obj.tolist()
        return json.JSONEncoder.default(self, obj)

def test_oI2Y():
    """
    Obtain the k-path for a test system with inversion symmetry and
    Bravais lattice (extended) oI2.
    """
    import spglib
    from seekpath import hpkot

    system = (
        [
            [8.186532843027958, 0.0, 0.0], 
            [0.0, 5.827100673574356, 0.0], 
            [0.0, 0.0, 5.8365707664294835]
        ], [
            [0.2905378000000001, 0.0, 0.0], 
            [0.2094621999999999, 0.0, 0.5], 
            [0.2094621999999999, 0.5, 0.5], 
            [0.7905378000000001, 0.0, 0.5], 
            [0.0, 0.75, 0.1751481400000001], 
            [0.5, 0.75, 0.32485186], 
            [0.7905378000000001, 0.5, 0.5], 
            [0.7094621999999999, 0.5, 0.0], 
            [0.7094621999999999, 0.0, 0.0], 
            [0.2905378, 0.5, 0.0], 
            [0.5, 0.25, 0.6751481400000001], 
            [0.0, 0.25, 0.8248518599999999], 
            [0.0, 0.75, 0.7411616600000002], 
            [0.5, 0.75, 0.7588383399999998], 
            [0.5, 0.25, 0.2411616600000002], 
            [0.0, 0.25, 0.2588383399999998], 
            [0.25, 0.75, 0.25], 
            [0.25, 0.25, 0.75], 
            [0.75, 0.25, 0.75], 
            [0.75, 0.75, 0.25]
        ], 
        [8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 8, 
        8, 38, 38, 38, 38, 72, 72, 72, 72])

    sym_ds = spglib.get_symmetry_dataset(system)
    print(json.dumps(sym_ds, indent=2, sort_keys=True, cls=NumpyEncoder))

    res = hpkot.get_path(system, with_time_reversal=False)
    assert res['bravais_lattice_extended'] == 'oI2'

if __name__ == "__main__":
    test_oI2Y()

I also report the output of the code (in particular, the json-serialized get_symmetry_dataset for the two versions):

In particular, the main issue is that in the "working" version, we have from spglib:

"std_lattice": [
    [8.186532843027958, 0.0, 0.0], 
    [0.0, 5.827100673574356, 0.0], 
    [0.0, 0.0, 5.8365707664294835]
  ]

while in the "failing" version:

"std_lattice": [
    [5.827100673574355, 0.0, 0.0], 
    [0.0, 8.186532843027958, 0.0], 
    [0.0, 0.0, 5.8365707664294835]
  ]

and since the order is different, the categorisation (2 vs. 3) by seek path is different.

I will investigate with the author of spglib.

atztogo commented 5 years ago

I am the author of spglib. This is probably due to spglib. oI2 or oI3 may be undetermined by the current implementation. So a detailed change introduced recently changed the choice. But this choice may/should be able to be determined. I will look at it deeper and get back here when I can resolve this issue.

giovannipizzi commented 5 years ago

Thanks a lot @atztogo!

atztogo commented 5 years ago

In orthorhombic system, I found four space group types, Cmma (67), Ccca(68), Ibca(73), Imma(74), were wrongly treated. They were exceptions than my previous definitions, i.e., my understanding on crystallography was not enough. This fix is in the develop branch (https://github.com/atztogo/spglib/commit/c61d10914c1260571605a1fcef0a62fbb1fdedcd). Previously degrees of freedom of a,b,c choice relied on number of hall symbols for each space group type. But now, degrees of freedom of a, b, c choice rely on the index of L(g) in N_\epsilon(g) of highest symmetry or equivalently the index of L(g) in N_A(g) (ITA 2016, 3.5. Normalizers of space groups and their use in crystallography).

The similar problem may exist in monoclinic system. I will work on it and once it will be solved, I will make a spglib release.

atztogo commented 5 years ago

I have released spglib-1.13.0 with the fix above. For monoclinic system, I made small update, but not like the orthorhombic case above. This version is better for the standardization.

giovannipizzi commented 5 years ago

Hi @atztogo thanks! I'm starting to prepare a new version of seekpath pointing to the new version of spglib.

A question: do you agree, then, that there is no oI space group with inversion symmetry that can have a being the largest vector (a>b and a>c)?

The four space groups with oI and inversion symmetry are:

For sure in 71 and 74 we always apply a<b. This is also the case for 72 and 73, right?

(And in the case of oI without inversion symmetry, instead, this is possible because there is space group Ima2 so the order of the three vectors is fixed by the symmetry, and we do not have freedom to choose their ordering, right?)

atztogo commented 5 years ago

Yes, most strictly, they are implemented as follows.

71: Immm, a < b + 1e-10, b < c + 1e-10 72: Ibam, a < b + 1e-10 73: Ibca, a < b + 1e-10, b < c + 1e-10 74: Imma, a < b + 1e-10

For Ima2, yes, it is fixed by symmetry. The fact that you wrote in the parenthesis is correct for 3D crystal. But I can't explain the reason in short works. It is not the matter of centring or inversion symmetry. Non-symmorphic operations make thing tricky. I just relied on the data base of affine normalizer to see how many degrees of freedom exist. Pbca is tricky. We have choices abc, bca, cab, but flip of two axies is not allowed. For this, I implement

61: Pbca, a < b + 1e-10, a < c + 1e-10.

giovannipizzi commented 5 years ago

Fixed in release v1.9.0 (in particular, in commit a13047801b119bc8859b6ae124e9fef77b5037af)