1# coding: utf-8 2# Copyright (c) Pymatgen Development Team. 3# Distributed under the terms of the MIT License. 4 5""" 6Magnetic space groups. 7""" 8 9import os 10import sqlite3 11import textwrap 12from array import array 13from fractions import Fraction 14 15import numpy as np 16from monty.design_patterns import cached_class 17 18from pymatgen.core.operations import MagSymmOp 19from pymatgen.electronic_structure.core import Magmom 20from pymatgen.symmetry.groups import SymmetryGroup, in_array_list 21from pymatgen.symmetry.settings import JonesFaithfulTransformation 22from pymatgen.util.string import transformation_to_string 23 24__author__ = "Matthew Horton, Shyue Ping Ong" 25 26 27MAGSYMM_DATA = os.path.join(os.path.dirname(__file__), "symm_data_magnetic.sqlite") 28 29 30@cached_class 31class MagneticSpaceGroup(SymmetryGroup): 32 """ 33 Representation of a magnetic space group. 34 """ 35 36 def __init__(self, id, setting_transformation="a,b,c;0,0,0"): 37 """ 38 Initializes a MagneticSpaceGroup from its Belov, Neronova and 39 Smirnova (BNS) number supplied as a list or its label supplied 40 as a string. To create a magnetic structure in pymatgen, the 41 Structure.from_magnetic_spacegroup() method can be used, which 42 relies on this class. 43 44 The main difference between magnetic space groups and normal 45 crystallographic space groups is the inclusion of a time reversal 46 operator that acts on an atom's magnetic moment. This is 47 indicated by a prime symbol (') next to the respective symmetry 48 operation in its label, e.g. the standard crystallographic 49 space group Pnma has magnetic subgroups Pn'ma, Pnm'a, Pnma', 50 Pn'm'a, Pnm'a', Pn'ma', Pn'm'a'. 51 52 The magnetic space groups are classified as one of 4 types 53 where G = magnetic space group, and F = parent crystallographic 54 space group: 55 56 1. G=F no time reversal, i.e. the same as corresponding 57 crystallographic group 58 2. G=F+F1', "grey" groups, where avg. magnetic moment is zero, 59 e.g. a paramagnet in zero ext. mag. field 60 3. G=D+(F-D)1', where D is an equi-translation subgroup of F of 61 index 2, lattice translations do not include time reversal 62 4. G=D+(F-D)1', where D is an equi-class subgroup of F of index 2 63 64 There are two common settings for magnetic space groups, BNS 65 and OG. In case 4, the BNS setting != OG setting, and so a 66 transformation to go between the two settings is required: 67 specifically, the BNS setting is derived from D, and the OG 68 setting is derived from F. 69 70 This means that the OG setting refers to the unit cell if magnetic 71 order is neglected, and requires multiple unit cells to reproduce 72 the full crystal periodicity when magnetic moments are present. 73 This does not make the OG setting, in general, useful for 74 electronic structure calculations and the BNS setting is preferred. 75 However, this class does contain information on the OG setting and 76 can be initialized from OG labels or numbers if required. 77 78 Conventions: ITC monoclinic unique axis b, monoclinic cell choice 1, 79 hexagonal axis for trigonal groups, origin choice 2 for groups with 80 more than one origin choice (ISO-MAG). 81 82 Raw data comes from ISO-MAG, ISOTROPY Software Suite, iso.byu.edu 83 http://stokes.byu.edu/iso/magnetic_data.txt 84 with kind permission from Professor Branton Campbell, BYU 85 86 Data originally compiled from: 87 (1) Daniel B. Litvin, Magnetic Group Tables (International Union 88 of Crystallography, 2013) www.iucr.org/publ/978-0-9553602-2-0. 89 (2) C. J. Bradley and A. P. Cracknell, The Mathematical Theory of 90 Symmetry in Solids (Clarendon Press, Oxford, 1972). 91 92 See http://stokes.byu.edu/iso/magneticspacegroupshelp.php for more 93 information on magnetic symmetry. 94 95 :param id: BNS number supplied as list of 2 ints or BNS label as 96 str or index as int (1-1651) to iterate over all space groups""" 97 98 self._data = {} 99 100 # Datafile is stored as sqlite3 database since (a) it can be easily 101 # queried for various different indexes (BNS/OG number/labels) and (b) 102 # allows binary data to be stored in a compact form similar to that in 103 # the source data file, significantly reducing file size. 104 # Note that a human-readable JSON format was tested first but was 20x 105 # larger and required *much* longer initial loading times. 106 107 # retrieve raw data 108 db = sqlite3.connect(MAGSYMM_DATA) 109 c = db.cursor() 110 if isinstance(id, str): 111 id = "".join(id.split()) # remove any white space 112 c.execute("SELECT * FROM space_groups WHERE BNS_label=?;", (id,)) 113 elif isinstance(id, list): 114 c.execute("SELECT * FROM space_groups WHERE BNS1=? AND BNS2=?;", (id[0], id[1])) 115 elif isinstance(id, int): 116 # OG3 index is a 'master' index, going from 1 to 1651 117 c.execute("SELECT * FROM space_groups WHERE OG3=?;", (id,)) 118 raw_data = list(c.fetchone()) 119 120 # Jones Faithful transformation 121 self.jf = JonesFaithfulTransformation.from_transformation_string("a,b,c;0,0,0") 122 if isinstance(setting_transformation, str): 123 if setting_transformation != "a,b,c;0,0,0": 124 self.jf = JonesFaithfulTransformation.from_transformation_string(setting_transformation) 125 elif isinstance(setting_transformation, JonesFaithfulTransformation): 126 if setting_transformation != self.jf: 127 self.jf = setting_transformation 128 129 self._data["magtype"] = raw_data[0] # int from 1 to 4 130 self._data["bns_number"] = [raw_data[1], raw_data[2]] 131 self._data["bns_label"] = raw_data[3] 132 self._data["og_number"] = [raw_data[4], raw_data[5], raw_data[6]] 133 self._data["og_label"] = raw_data[7] # can differ from BNS_label 134 135 def _get_point_operator(idx): 136 """Retrieve information on point operator (rotation matrix and Seitz label).""" 137 hex = self._data["bns_number"][0] >= 143 and self._data["bns_number"][0] <= 194 138 c.execute( 139 "SELECT symbol, matrix FROM point_operators WHERE idx=? AND hex=?;", 140 (idx - 1, hex), 141 ) 142 op = c.fetchone() 143 op = { 144 "symbol": op[0], 145 "matrix": np.array(op[1].split(","), dtype="f").reshape(3, 3), 146 } 147 return op 148 149 def _parse_operators(b): 150 """Parses compact binary representation into list of MagSymmOps.""" 151 if len(b) == 0: # e.g. if magtype != 4, OG setting == BNS setting, and b == [] for OG symmops 152 return None 153 raw_symops = [b[i : i + 6] for i in range(0, len(b), 6)] 154 155 symops = [] 156 157 for r in raw_symops: 158 point_operator = _get_point_operator(r[0]) 159 translation_vec = [r[1] / r[4], r[2] / r[4], r[3] / r[4]] 160 time_reversal = r[5] 161 op = MagSymmOp.from_rotation_and_translation_and_time_reversal( 162 rotation_matrix=point_operator["matrix"], 163 translation_vec=translation_vec, 164 time_reversal=time_reversal, 165 ) 166 # store string representation, e.g. (2x|1/2,1/2,1/2)' 167 seitz = "({0}|{1},{2},{3})".format( 168 point_operator["symbol"], 169 Fraction(translation_vec[0]), 170 Fraction(translation_vec[1]), 171 Fraction(translation_vec[2]), 172 ) 173 if time_reversal == -1: 174 seitz += "'" 175 symops.append({"op": op, "str": seitz}) 176 177 return symops 178 179 def _parse_wyckoff(b): 180 """Parses compact binary representation into list of Wyckoff sites.""" 181 if len(b) == 0: 182 return None 183 184 wyckoff_sites = [] 185 186 def get_label(idx): 187 if idx <= 25: 188 return chr(97 + idx) # returns a-z when idx 0-25 189 return "alpha" # when a-z labels exhausted, use alpha, only relevant for a few space groups 190 191 o = 0 # offset 192 n = 1 # nth Wyckoff site 193 num_wyckoff = b[0] 194 while len(wyckoff_sites) < num_wyckoff: 195 m = b[1 + o] # multiplicity 196 label = str(b[2 + o] * m) + get_label(num_wyckoff - n) 197 sites = [] 198 for j in range(m): 199 s = b[3 + o + (j * 22) : 3 + o + (j * 22) + 22] # data corresponding to specific Wyckoff position 200 translation_vec = [s[0] / s[3], s[1] / s[3], s[2] / s[3]] 201 matrix = [ 202 [s[4], s[7], s[10]], 203 [s[5], s[8], s[11]], 204 [s[6], s[9], s[12]], 205 ] 206 matrix_magmom = [ 207 [s[13], s[16], s[19]], 208 [s[14], s[17], s[20]], 209 [s[15], s[18], s[21]], 210 ] 211 # store string representation, e.g. (x,y,z;mx,my,mz) 212 wyckoff_str = "({};{})".format( 213 transformation_to_string(matrix, translation_vec), 214 transformation_to_string(matrix_magmom, c="m"), 215 ) 216 sites.append( 217 { 218 "translation_vec": translation_vec, 219 "matrix": matrix, 220 "matrix_magnetic": matrix_magmom, 221 "str": wyckoff_str, 222 } 223 ) 224 225 # only keeping string representation of Wyckoff sites for now 226 # could do something else with these in future 227 wyckoff_sites.append({"label": label, "str": " ".join([s["str"] for s in sites])}) 228 n += 1 229 o += m * 22 + 2 230 231 return wyckoff_sites 232 233 def _parse_lattice(b): 234 """Parses compact binary representation into list of lattice vectors/centerings.""" 235 if len(b) == 0: 236 return None 237 raw_lattice = [b[i : i + 4] for i in range(0, len(b), 4)] 238 239 lattice = [] 240 241 for r in raw_lattice: 242 lattice.append( 243 { 244 "vector": [r[0] / r[3], r[1] / r[3], r[2] / r[3]], 245 "str": "({0},{1},{2})+".format( 246 Fraction(r[0] / r[3]).limit_denominator(), 247 Fraction(r[1] / r[3]).limit_denominator(), 248 Fraction(r[2] / r[3]).limit_denominator(), 249 ), 250 } 251 ) 252 253 return lattice 254 255 def _parse_transformation(b): 256 """Parses compact binary representation into transformation between OG and BNS settings.""" 257 if len(b) == 0: 258 return None 259 # capital letters used here by convention, 260 # IUCr defines P and p specifically 261 P = [[b[0], b[3], b[6]], [b[1], b[4], b[7]], [b[2], b[5], b[8]]] 262 p = [b[9] / b[12], b[10] / b[12], b[11] / b[12]] 263 P = np.array(P).transpose() 264 P_string = transformation_to_string(P, components=("a", "b", "c")) 265 p_string = "{},{},{}".format( 266 Fraction(p[0]).limit_denominator(), 267 Fraction(p[1]).limit_denominator(), 268 Fraction(p[2]).limit_denominator(), 269 ) 270 return P_string + ";" + p_string 271 272 for i in range(8, 15): 273 try: 274 raw_data[i] = array("b", raw_data[i]) # construct array from sql binary blobs 275 except Exception: 276 # array() behavior changed, need to explicitly convert buffer to str in earlier Python 277 raw_data[i] = array("b", str(raw_data[i])) 278 279 self._data["og_bns_transform"] = _parse_transformation(raw_data[8]) 280 self._data["bns_operators"] = _parse_operators(raw_data[9]) 281 self._data["bns_lattice"] = _parse_lattice(raw_data[10]) 282 self._data["bns_wyckoff"] = _parse_wyckoff(raw_data[11]) 283 self._data["og_operators"] = _parse_operators(raw_data[12]) 284 self._data["og_lattice"] = _parse_lattice(raw_data[13]) 285 self._data["og_wyckoff"] = _parse_wyckoff(raw_data[14]) 286 287 db.close() 288 289 @classmethod 290 def from_og(cls, id): 291 """ 292 Initialize from Opechowski and Guccione (OG) label or number. 293 294 :param id: OG number supplied as list of 3 ints or 295 or OG label as str 296 :return: 297 """ 298 299 db = sqlite3.connect(MAGSYMM_DATA) 300 c = db.cursor() 301 if isinstance(id, str): 302 c.execute("SELECT BNS_label FROM space_groups WHERE OG_label=?", (id,)) 303 elif isinstance(id, list): 304 c.execute( 305 "SELECT BNS_label FROM space_groups WHERE OG1=? and OG2=? and OG3=?", 306 (id[0], id[1], id[2]), 307 ) 308 bns_label = c.fetchone()[0] 309 db.close() 310 311 return cls(bns_label) 312 313 def __eq__(self, other): 314 return self._data == other._data 315 316 @property 317 def crystal_system(self): 318 """ 319 :return: Crystal system, e.g., cubic, hexagonal, etc. 320 """ 321 i = self._data["bns_number"][0] 322 if i <= 2: 323 return "triclinic" 324 if i <= 15: 325 return "monoclinic" 326 if i <= 74: 327 return "orthorhombic" 328 if i <= 142: 329 return "tetragonal" 330 if i <= 167: 331 return "trigonal" 332 if i <= 194: 333 return "hexagonal" 334 return "cubic" 335 336 @property 337 def sg_symbol(self): 338 """ 339 :return: Space group symbol 340 """ 341 return self._data["bns_label"] 342 343 @property 344 def symmetry_ops(self): 345 """ 346 Retrieve magnetic symmetry operations of the space group. 347 :return: List of :class:`pymatgen.core.operations.MagSymmOp` 348 """ 349 ops = [op_data["op"] for op_data in self._data["bns_operators"]] 350 351 # add lattice centerings 352 centered_ops = [] 353 lattice_vectors = [l["vector"] for l in self._data["bns_lattice"]] 354 355 for vec in lattice_vectors: 356 if not (np.array_equal(vec, [1, 0, 0]) or np.array_equal(vec, [0, 1, 0]) or np.array_equal(vec, [0, 0, 1])): 357 for op in ops: 358 new_vec = op.translation_vector + vec 359 new_op = MagSymmOp.from_rotation_and_translation_and_time_reversal( 360 op.rotation_matrix, 361 translation_vec=new_vec, 362 time_reversal=op.time_reversal, 363 ) 364 centered_ops.append(new_op) 365 366 ops = ops + centered_ops 367 368 # apply jones faithful transformation 369 ops = [self.jf.transform_symmop(op) for op in ops] 370 371 return ops 372 373 def get_orbit(self, p, m, tol=1e-5): 374 """ 375 Returns the orbit for a point and its associated magnetic moment. 376 377 Args: 378 p: Point as a 3x1 array. 379 m: A magnetic moment, compatible with 380 :class:`pymatgen.electronic_structure.core.Magmom` 381 tol: Tolerance for determining if sites are the same. 1e-5 should 382 be sufficient for most purposes. Set to 0 for exact matching 383 (and also needed for symbolic orbits). 384 385 Returns: 386 (([array], [array])) Tuple of orbit for point and magnetic moments for orbit. 387 """ 388 orbit = [] 389 orbit_magmoms = [] 390 m = Magmom(m) 391 for o in self.symmetry_ops: 392 pp = o.operate(p) 393 pp = np.mod(np.round(pp, decimals=10), 1) 394 mm = o.operate_magmom(m) 395 if not in_array_list(orbit, pp, tol=tol): 396 orbit.append(pp) 397 orbit_magmoms.append(mm) 398 return orbit, orbit_magmoms 399 400 def is_compatible(self, lattice, tol=1e-5, angle_tol=5): 401 """ 402 Checks whether a particular lattice is compatible with the 403 *conventional* unit cell. 404 405 Args: 406 lattice (Lattice): A Lattice. 407 tol (float): The tolerance to check for equality of lengths. 408 angle_tol (float): The tolerance to check for equality of angles 409 in degrees. 410 """ 411 # function from pymatgen.symmetry.groups.SpaceGroup 412 abc = lattice.lengths 413 angles = lattice.angles 414 crys_system = self.crystal_system 415 416 def check(param, ref, tolerance): 417 return all(abs(i - j) < tolerance for i, j in zip(param, ref) if j is not None) 418 419 if crys_system == "cubic": 420 a = abc[0] 421 return check(abc, [a, a, a], tol) and check(angles, [90, 90, 90], angle_tol) 422 if crys_system == "hexagonal" or (crys_system == "trigonal" and self.sg_symbol.endswith("H")): 423 a = abc[0] 424 return check(abc, [a, a, None], tol) and check(angles, [90, 90, 120], angle_tol) 425 if crys_system == "trigonal": 426 a = abc[0] 427 return check(abc, [a, a, a], tol) 428 if crys_system == "tetragonal": 429 a = abc[0] 430 return check(abc, [a, a, None], tol) and check(angles, [90, 90, 90], angle_tol) 431 if crys_system == "orthorhombic": 432 return check(angles, [90, 90, 90], angle_tol) 433 if crys_system == "monoclinic": 434 return check(angles, [90, None, 90], angle_tol) 435 return True 436 437 def data_str(self, include_og=True): 438 """ 439 Get description of all data, including information for OG setting. 440 :return: str 441 """ 442 443 # __str__() omits information on OG setting to reduce confusion 444 # as to which set of symops are active, this property gives 445 # all stored data including OG setting 446 447 desc = {} # dictionary to hold description strings 448 description = "" 449 450 # parse data into strings 451 452 # indicate if non-standard setting specified 453 if self.jf != JonesFaithfulTransformation.from_transformation_string("a,b,c;0,0,0"): 454 description += "Non-standard setting: .....\n" 455 description += self.jf.__repr__() 456 description += "\n\nStandard setting information: \n" 457 458 desc["magtype"] = self._data["magtype"] 459 desc["bns_number"] = ".".join(map(str, self._data["bns_number"])) 460 desc["bns_label"] = self._data["bns_label"] 461 desc["og_id"] = ( 462 "\t\tOG: " + ".".join(map(str, self._data["og_number"])) + " " + self._data["og_label"] 463 if include_og 464 else "" 465 ) 466 desc["bns_operators"] = " ".join([op_data["str"] for op_data in self._data["bns_operators"]]) 467 468 desc["bns_lattice"] = ( 469 " ".join([lattice_data["str"] for lattice_data in self._data["bns_lattice"][3:]]) 470 if len(self._data["bns_lattice"]) > 3 471 else "" 472 ) # don't show (1,0,0)+ (0,1,0)+ (0,0,1)+ 473 474 desc["bns_wyckoff"] = "\n".join( 475 [ 476 textwrap.fill( 477 wyckoff_data["str"], 478 initial_indent=wyckoff_data["label"] + " ", 479 subsequent_indent=" " * len(wyckoff_data["label"] + " "), 480 break_long_words=False, 481 break_on_hyphens=False, 482 ) 483 for wyckoff_data in self._data["bns_wyckoff"] 484 ] 485 ) 486 487 desc["og_bns_transformation"] = ( 488 "OG-BNS Transform: ({})\n".format(self._data["og_bns_transform"]) 489 if desc["magtype"] == 4 and include_og 490 else "" 491 ) 492 493 bns_operators_prefix = "Operators{}: ".format(" (BNS)" if desc["magtype"] == 4 and include_og else "") 494 bns_wyckoff_prefix = "Wyckoff Positions{}: ".format(" (BNS)" if desc["magtype"] == 4 and include_og else "") 495 496 # apply textwrap on long lines 497 desc["bns_operators"] = textwrap.fill( 498 desc["bns_operators"], 499 initial_indent=bns_operators_prefix, 500 subsequent_indent=" " * len(bns_operators_prefix), 501 break_long_words=False, 502 break_on_hyphens=False, 503 ) 504 505 description += ( 506 "BNS: {d[bns_number]} {d[bns_label]}{d[og_id]}\n" 507 "{d[og_bns_transformation]}" 508 "{d[bns_operators]}\n" 509 "{bns_wyckoff_prefix}{d[bns_lattice]}\n" 510 "{d[bns_wyckoff]}" 511 ).format(d=desc, bns_wyckoff_prefix=bns_wyckoff_prefix) 512 513 if desc["magtype"] == 4 and include_og: 514 515 desc["og_operators"] = " ".join([op_data["str"] for op_data in self._data["og_operators"]]) 516 517 # include all lattice vectors because (1,0,0)+ (0,1,0)+ (0,0,1)+ 518 # not always present in OG setting 519 desc["og_lattice"] = " ".join([lattice_data["str"] for lattice_data in self._data["og_lattice"]]) 520 521 desc["og_wyckoff"] = "\n".join( 522 [ 523 textwrap.fill( 524 wyckoff_data["str"], 525 initial_indent=wyckoff_data["label"] + " ", 526 subsequent_indent=" " * len(wyckoff_data["label"] + " "), 527 break_long_words=False, 528 break_on_hyphens=False, 529 ) 530 for wyckoff_data in self._data["og_wyckoff"] 531 ] 532 ) 533 534 og_operators_prefix = "Operators (OG): " 535 536 # apply textwrap on long lines 537 desc["og_operators"] = textwrap.fill( 538 desc["og_operators"], 539 initial_indent=og_operators_prefix, 540 subsequent_indent=" " * len(og_operators_prefix), 541 break_long_words=False, 542 break_on_hyphens=False, 543 ) 544 545 description += ( 546 "\n{d[og_operators]}\n" "Wyckoff Positions (OG): {d[og_lattice]}\n" "{d[og_wyckoff]}" 547 ).format(d=desc) 548 elif desc["magtype"] == 4: 549 description += "\nAlternative OG setting exists for this space group." 550 551 return description 552 553 def __str__(self): 554 """ 555 String representation of the space group, specifying the setting 556 of the space group, its magnetic symmetry operators and Wyckoff 557 positions. 558 :return: str 559 """ 560 return self.data_str(include_og=False) 561 562 563def _write_all_magnetic_space_groups_to_file(filename): 564 """ 565 Write all magnetic space groups to a human-readable text file. 566 Should contain same information as text files provided by ISO-MAG. 567 :param filename: 568 :return: 569 """ 570 s = ( 571 "Data parsed from raw data from:\n" 572 "ISO-MAG, ISOTROPY Software Suite, iso.byu.edu\n" 573 "http://stokes.byu.edu/iso/magnetic_data.txt\n" 574 "Used with kind permission from Professor Branton Campbell, BYU\n\n" 575 ) 576 all_msgs = [] 577 for i in range(1, 1652): 578 all_msgs.append(MagneticSpaceGroup(i)) 579 for msg in all_msgs: 580 s += "\n{}\n\n--------\n".format(msg.data_str()) 581 with open(filename, "w") as f: 582 f.write(s) 583