1# coding: utf-8
2# Copyright (c) Pymatgen Development Team.
3# Distributed under the terms of the MIT License.
4
5
6import unittest
7
8from pymatgen.core.lattice import Lattice
9from pymatgen.core.structure import Structure
10from pymatgen.symmetry.kpath import KPathSeek
11from pymatgen.util.testing import PymatgenTest
12
13try:
14    from seekpath import get_path  # type: ignore
15except ImportError:
16    get_path = None
17
18
19class KPathSeekTest(PymatgenTest):
20    @unittest.skipIf(get_path is None, "No seek path present.")
21    def test_kpath_generation(self):
22        triclinic = [1, 2]
23        monoclinic = range(3, 16)
24        orthorhombic = range(16, 75)
25        tetragonal = range(75, 143)
26        rhombohedral = range(143, 168)
27        hexagonal = range(168, 195)
28        cubic = range(195, 231)
29
30        species = ["K", "La", "Ti"]
31        coords = [[0.345, 5, 0.77298], [0.1345, 5.1, 0.77298], [0.7, 0.8, 0.9]]
32        for i in range(230):
33            sg_num = i + 1
34            if sg_num in triclinic:
35                lattice = Lattice(
36                    [
37                        [3.0233057319441246, 1, 0],
38                        [0, 7.9850357844548681, 1],
39                        [0, 1.2, 8.1136762279561818],
40                    ]
41                )
42            elif sg_num in monoclinic:
43                lattice = Lattice.monoclinic(2, 9, 1, 99)
44            elif sg_num in orthorhombic:
45                lattice = Lattice.orthorhombic(2, 9, 1)
46            elif sg_num in tetragonal:
47                lattice = Lattice.tetragonal(2, 9)
48            elif sg_num in rhombohedral:
49                lattice = Lattice.hexagonal(2, 95)
50            elif sg_num in hexagonal:
51                lattice = Lattice.hexagonal(2, 9)
52            elif sg_num in cubic:
53                lattice = Lattice.cubic(2)
54
55            struct = Structure.from_spacegroup(sg_num, lattice, species, coords)
56            kpath = KPathSeek(struct)  # Throws error if something doesn't work, causing test to fail.
57            kpoints = kpath.get_kpoints()  # noqa: F841
58
59    @unittest.skipIf(get_path is None, "No seek path present.")
60    def test_kpath_acentered(self):
61        species = ["K", "La", "Ti"]
62        coords = [[0.345, 5, 0.77298], [0.1345, 5.1, 0.77298], [0.7, 0.8, 0.9]]
63        lattice = Lattice.orthorhombic(2, 9, 1)
64        struct = Structure.from_spacegroup(38, lattice, species, coords)
65        kpath = KPathSeek(struct)
66
67        kpoints = kpath._kpath["kpoints"]
68        labels = list(kpoints.keys())
69
70        self.assertEqual(
71            sorted(labels),
72            sorted(
73                [
74                    "B_0",
75                    "B_2",
76                    "DELTA_0",
77                    "F_0",
78                    "GAMMA",
79                    "G_0",
80                    "G_2",
81                    "R",
82                    "R_2",
83                    "S",
84                    "T",
85                    "T_2",
86                    "Y",
87                    "Z",
88                    "Z_2",
89                ]
90            ),
91        )
92
93        self.assertAlmostEqual(kpoints["GAMMA"][0], 0.0)
94        self.assertAlmostEqual(kpoints["GAMMA"][1], 0.0)
95        self.assertAlmostEqual(kpoints["GAMMA"][2], 0.0)
96
97        self.assertAlmostEqual(kpoints["Y"][0], 0.5)
98        self.assertAlmostEqual(kpoints["Y"][1], 0.5)
99        self.assertAlmostEqual(kpoints["Y"][2], 0.0)
100
101        self.assertAlmostEqual(kpoints["T"][0], 0.5)
102        self.assertAlmostEqual(kpoints["T"][1], 0.5)
103        self.assertAlmostEqual(kpoints["T"][2], 0.5)
104
105        self.assertAlmostEqual(kpoints["T_2"][0], 0.5)
106        self.assertAlmostEqual(kpoints["T_2"][1], 0.5)
107        self.assertAlmostEqual(kpoints["T_2"][2], -0.5)
108
109        self.assertAlmostEqual(kpoints["Z"][0], 0.0)
110        self.assertAlmostEqual(kpoints["Z"][1], 0.0)
111        self.assertAlmostEqual(kpoints["Z"][2], 0.5)
112
113        self.assertAlmostEqual(kpoints["Z_2"][0], 0.0)
114        self.assertAlmostEqual(kpoints["Z_2"][1], 0.0)
115        self.assertAlmostEqual(kpoints["Z_2"][2], -0.5)
116
117        self.assertAlmostEqual(kpoints["S"][0], 0.0)
118        self.assertAlmostEqual(kpoints["S"][1], 0.5)
119        self.assertAlmostEqual(kpoints["S"][2], 0.0)
120
121        self.assertAlmostEqual(kpoints["R"][0], 0.0)
122        self.assertAlmostEqual(kpoints["R"][1], 0.5)
123        self.assertAlmostEqual(kpoints["R"][2], 0.5)
124
125        self.assertAlmostEqual(kpoints["R_2"][0], 0.0)
126        self.assertAlmostEqual(kpoints["R_2"][1], 0.5)
127        self.assertAlmostEqual(kpoints["R_2"][2], -0.5)
128
129        self.assertAlmostEqual(kpoints["DELTA_0"][0], -0.25308641975308643)
130        self.assertAlmostEqual(kpoints["DELTA_0"][1], 0.25308641975308643)
131        self.assertAlmostEqual(kpoints["DELTA_0"][2], 0.0)
132
133        self.assertAlmostEqual(kpoints["F_0"][0], 0.25308641975308643)
134        self.assertAlmostEqual(kpoints["F_0"][1], 0.7469135802469136)
135        self.assertAlmostEqual(kpoints["F_0"][2], 0.0)
136
137        self.assertAlmostEqual(kpoints["B_0"][0], -0.25308641975308643)
138        self.assertAlmostEqual(kpoints["B_0"][1], 0.25308641975308643)
139        self.assertAlmostEqual(kpoints["B_0"][2], 0.5)
140
141        self.assertAlmostEqual(kpoints["B_2"][0], -0.25308641975308643)
142        self.assertAlmostEqual(kpoints["B_2"][1], 0.25308641975308643)
143        self.assertAlmostEqual(kpoints["B_2"][2], -0.5)
144
145        self.assertAlmostEqual(kpoints["G_0"][0], 0.25308641975308643)
146        self.assertAlmostEqual(kpoints["G_0"][1], 0.7469135802469136)
147        self.assertAlmostEqual(kpoints["G_0"][2], 0.5)
148
149        self.assertAlmostEqual(kpoints["G_2"][0], 0.25308641975308643)
150        self.assertAlmostEqual(kpoints["G_2"][1], 0.7469135802469136)
151        self.assertAlmostEqual(kpoints["G_2"][2], -0.5)
152
153
154if __name__ == "__main__":
155    unittest.main()
156