1# Copyright (C) 2015 Atsushi Togo 2# All rights reserved. 3# 4# This file is part of spglib. 5# 6# Redistribution and use in source and binary forms, with or without 7# modification, are permitted provided that the following conditions 8# are met: 9# 10# * Redistributions of source code must retain the above copyright 11# notice, this list of conditions and the following disclaimer. 12# 13# * Redistributions in binary form must reproduce the above copyright 14# notice, this list of conditions and the following disclaimer in 15# the documentation and/or other materials provided with the 16# distribution. 17# 18# * Neither the name of the spglib project nor the names of its 19# contributors may be used to endorse or promote products derived 20# from this software without specific prior written permission. 21# 22# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 23# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 24# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS 25# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE 26# COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, 27# INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, 28# BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 29# LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER 30# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT 31# LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN 32# ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE 33# POSSIBILITY OF SUCH DAMAGE. 34 35from . import _spglib as spg 36import numpy as np 37 38class SpglibError(object): 39 message = "no error" 40 41spglib_error = SpglibError() 42 43def get_version(): 44 _set_no_error() 45 return tuple(spg.version()) 46 47def get_symmetry(cell, symprec=1e-5, angle_tolerance=-1.0): 48 """This gives crystal symmetry operations from a crystal structure. 49 50 Args: 51 cell: Crystal structrue given either in Atoms object or tuple. 52 In the case given by a tuple, it has to follow the form below, 53 (Lattice parameters in a 3x3 array (see the detail below), 54 Fractional atomic positions in an Nx3 array, 55 Integer numbers to distinguish species in a length N array, 56 (optional) Collinear magnetic moments in a length N array), 57 where N is the number of atoms. 58 Lattice parameters are given in the form: 59 [[a_x, a_y, a_z], 60 [b_x, b_y, b_z], 61 [c_x, c_y, c_z]] 62 symprec: 63 float: Symmetry search tolerance in the unit of length. 64 angle_tolerance: 65 float: Symmetry search tolerance in the unit of angle deg. 66 If the value is negative, an internally optimized routine 67 is used to judge symmetry. 68 69 Return: 70 dictionary: Rotation parts and translation parts. 71 72 'rotations': Gives the numpy 'intc' array of the rotation matrices. 73 'translations': Gives the numpy 'double' array of fractional 74 translations with respect to a, b, c axes. 75 """ 76 _set_no_error() 77 78 lattice, positions, numbers, magmoms = _expand_cell(cell) 79 if lattice is None: 80 return None 81 82 multi = 48 * len(positions) 83 rotation = np.zeros((multi, 3, 3), dtype='intc') 84 translation = np.zeros((multi, 3), dtype='double') 85 86 # Get symmetry operations 87 if magmoms is None: 88 dataset = get_symmetry_dataset(cell, 89 symprec=symprec, 90 angle_tolerance=angle_tolerance) 91 if dataset is None: 92 return None 93 else: 94 return {'rotations': dataset['rotations'], 95 'translations': dataset['translations'], 96 'equivalent_atoms': dataset['equivalent_atoms']} 97 else: 98 equivalent_atoms = np.zeros(len(magmoms), dtype='intc') 99 num_sym = spg.symmetry_with_collinear_spin(rotation, 100 translation, 101 equivalent_atoms, 102 lattice, 103 positions, 104 numbers, 105 magmoms, 106 symprec, 107 angle_tolerance) 108 _set_error_message() 109 if num_sym == 0: 110 return None 111 else: 112 return {'rotations': np.array(rotation[:num_sym], 113 dtype='intc', order='C'), 114 'translations': np.array(translation[:num_sym], 115 dtype='double', order='C'), 116 'equivalent_atoms': equivalent_atoms} 117 118def get_symmetry_dataset(cell, 119 symprec=1e-5, 120 angle_tolerance=-1.0, 121 hall_number=0): 122 """Search symmetry dataset from an input cell. 123 124 Args: 125 cell, symprec, angle_tolerance: 126 See the docstring of get_symmetry. 127 hall_number: If a serial number of Hall symbol (>0) is given, 128 the database corresponding to the Hall symbol is made. 129 130 Return: 131 A dictionary is returned. 132 133 number: 134 int: International space group number 135 international: 136 str: International symbol 137 hall: 138 str: Hall symbol 139 choice: 140 str: Centring, origin, basis vector setting 141 transformation_matrix: 142 3x3 float matrix: 143 Transformation matrix from input lattice to standardized lattice 144 L^original = L^standardized * Tmat 145 origin shift: 146 float vecotr: Origin shift from standardized to input origin 147 rotations, translations: 148 3x3 int matrix, float vector: 149 Rotation matrices and translation vectors. Space group 150 operations are obtained by 151 [(r,t) for r, t in zip(rotations, translations)] 152 wyckoffs: 153 List of characters: Wyckoff letters 154 std_lattice, std_positions, std_types: 155 3x3 float matrix, Nx3 float vectors, list of int: 156 Standardized unit cell 157 pointgroup: 158 str: Pointgroup symbol 159 160 If it fails, None is returned. 161 """ 162 _set_no_error() 163 164 lattice, positions, numbers, _ = _expand_cell(cell) 165 if lattice is None: 166 return None 167 168 spg_ds = spg.dataset(lattice, positions, numbers, hall_number, 169 symprec, angle_tolerance) 170 if spg_ds is None: 171 _set_error_message() 172 return None 173 174 keys = ('number', 175 'hall_number', 176 'international', 177 'hall', 178 'choice', 179 'transformation_matrix', 180 'origin_shift', 181 'rotations', 182 'translations', 183 'wyckoffs', 184 'equivalent_atoms', 185 'std_lattice', 186 'std_types', 187 'std_positions', 188 # 'pointgroup_number', 189 'pointgroup') 190 dataset = {} 191 for key, data in zip(keys, spg_ds): 192 dataset[key] = data 193 194 dataset['international'] = dataset['international'].strip() 195 dataset['hall'] = dataset['hall'].strip() 196 dataset['choice'] = dataset['choice'].strip() 197 dataset['transformation_matrix'] = np.array( 198 dataset['transformation_matrix'], dtype='double', order='C') 199 dataset['origin_shift'] = np.array(dataset['origin_shift'], dtype='double') 200 dataset['rotations'] = np.array(dataset['rotations'], 201 dtype='intc', order='C') 202 dataset['translations'] = np.array(dataset['translations'], 203 dtype='double', order='C') 204 letters = "abcdefghijklmnopqrstuvwxyzABCDEFGHIJKLMNOPQRSTUVWXYZ" 205 dataset['wyckoffs'] = [letters[x] for x in dataset['wyckoffs']] 206 dataset['equivalent_atoms'] = np.array(dataset['equivalent_atoms'], 207 dtype='intc') 208 dataset['std_lattice'] = np.array(np.transpose(dataset['std_lattice']), 209 dtype='double', order='C') 210 dataset['std_types'] = np.array(dataset['std_types'], dtype='intc') 211 dataset['std_positions'] = np.array(dataset['std_positions'], 212 dtype='double', order='C') 213 dataset['pointgroup'] = dataset['pointgroup'].strip() 214 215 _set_error_message() 216 return dataset 217 218def get_spacegroup(cell, symprec=1e-5, angle_tolerance=-1.0, symbol_type=0): 219 """Return space group in international table symbol and number as a string. 220 221 If it fails, None is returned. 222 """ 223 _set_no_error() 224 225 dataset = get_symmetry_dataset(cell, 226 symprec=symprec, 227 angle_tolerance=angle_tolerance) 228 if dataset is None: 229 return None 230 231 spg_type = get_spacegroup_type(dataset['hall_number']) 232 if symbol_type == 1: 233 return "%s (%d)" % (spg_type['schoenflies'], dataset['number']) 234 else: 235 return "%s (%d)" % (spg_type['international_short'], dataset['number']) 236 237def get_hall_number_from_symmetry(rotations, translations, symprec=1e-5): 238 """Hall number is obtained from a set of symmetry operations. 239 240 If it fails, None is returned. 241 """ 242 243 hall_number = spg.hall_number_from_symmetry(rotations, translations, symprec) 244 return hall_number 245 246def get_spacegroup_type(hall_number): 247 """Translate Hall number to space group type information. 248 249 If it fails, None is returned. 250 """ 251 _set_no_error() 252 253 keys = ('number', 254 'international_short', 255 'international_full', 256 'international', 257 'schoenflies', 258 'hall_symbol', 259 'choice', 260 'pointgroup_schoenflies', 261 'pointgroup_international', 262 'arithmetic_crystal_class_number', 263 'arithmetic_crystal_class_symbol') 264 spg_type_list = spg.spacegroup_type(hall_number) 265 _set_error_message() 266 267 if spg_type_list is not None: 268 spg_type = dict(zip(keys, spg_type_list)) 269 for key in spg_type: 270 if key != 'number' and key != 'arithmetic_crystal_class_number': 271 spg_type[key] = spg_type[key].strip() 272 return spg_type 273 else: 274 return None 275 276def get_pointgroup(rotations): 277 """Return point group in international table symbol and number. 278 279 The symbols are mapped to the numbers as follows: 280 1 "1 " 281 2 "-1 " 282 3 "2 " 283 4 "m " 284 5 "2/m " 285 6 "222 " 286 7 "mm2 " 287 8 "mmm " 288 9 "4 " 289 10 "-4 " 290 11 "4/m " 291 12 "422 " 292 13 "4mm " 293 14 "-42m " 294 15 "4/mmm" 295 16 "3 " 296 17 "-3 " 297 18 "32 " 298 19 "3m " 299 20 "-3m " 300 21 "6 " 301 22 "-6 " 302 23 "6/m " 303 24 "622 " 304 25 "6mm " 305 26 "-62m " 306 27 "6/mmm" 307 28 "23 " 308 29 "m-3 " 309 30 "432 " 310 31 "-43m " 311 32 "m-3m " 312 """ 313 _set_no_error() 314 315 # (symbol, pointgroup_number, transformation_matrix) 316 pointgroup = spg.pointgroup(np.array(rotations, dtype='intc', order='C')) 317 _set_error_message() 318 return pointgroup 319 320def standardize_cell(cell, 321 to_primitive=False, 322 no_idealize=False, 323 symprec=1e-5, 324 angle_tolerance=-1.0): 325 """Return standardized cell. 326 327 Args: 328 cell, symprec, angle_tolerance: 329 See the docstring of get_symmetry. 330 to_primitive: 331 bool: If True, the standardized primitive cell is created. 332 no_idealize: 333 bool: If True, it is disabled to idealize lengths and angles of 334 basis vectors and positions of atoms according to crystal 335 symmetry. 336 Return: 337 The standardized unit cell or primitive cell is returned by a tuple of 338 (lattice, positions, numbers). 339 If it fails, None is returned. 340 """ 341 _set_no_error() 342 343 lattice, _positions, _numbers, _ = _expand_cell(cell) 344 if lattice is None: 345 return None 346 347 # Atomic positions have to be specified by scaled positions for spglib. 348 num_atom = len(_positions) 349 positions = np.zeros((num_atom * 4, 3), dtype='double', order='C') 350 positions[:num_atom] = _positions 351 numbers = np.zeros(num_atom * 4, dtype='intc') 352 numbers[:num_atom] = _numbers 353 num_atom_std = spg.standardize_cell(lattice, 354 positions, 355 numbers, 356 num_atom, 357 to_primitive * 1, 358 no_idealize * 1, 359 symprec, 360 angle_tolerance) 361 _set_error_message() 362 363 if num_atom_std > 0: 364 return (np.array(lattice.T, dtype='double', order='C'), 365 np.array(positions[:num_atom_std], dtype='double', order='C'), 366 np.array(numbers[:num_atom_std], dtype='intc')) 367 else: 368 return None 369 370def refine_cell(cell, symprec=1e-5, angle_tolerance=-1.0): 371 """Return refined cell. 372 373 The standardized unit cell is returned by a tuple of 374 (lattice, positions, numbers). 375 If it fails, None is returned. 376 """ 377 _set_no_error() 378 379 lattice, _positions, _numbers, _ = _expand_cell(cell) 380 if lattice is None: 381 return None 382 383 # Atomic positions have to be specified by scaled positions for spglib. 384 num_atom = len(_positions) 385 positions = np.zeros((num_atom * 4, 3), dtype='double', order='C') 386 positions[:num_atom] = _positions 387 numbers = np.zeros(num_atom * 4, dtype='intc') 388 numbers[:num_atom] = _numbers 389 num_atom_std = spg.refine_cell(lattice, 390 positions, 391 numbers, 392 num_atom, 393 symprec, 394 angle_tolerance) 395 _set_error_message() 396 397 if num_atom_std > 0: 398 return (np.array(lattice.T, dtype='double', order='C'), 399 np.array(positions[:num_atom_std], dtype='double', order='C'), 400 np.array(numbers[:num_atom_std], dtype='intc')) 401 else: 402 return None 403 404def find_primitive(cell, symprec=1e-5, angle_tolerance=-1.0): 405 """Primitive cell is searched in the input cell. 406 407 The primitive cell is returned by a tuple of (lattice, positions, numbers). 408 If it fails, None is returned. 409 """ 410 _set_no_error() 411 412 lattice, positions, numbers, _ = _expand_cell(cell) 413 if lattice is None: 414 return None 415 416 num_atom_prim = spg.primitive(lattice, 417 positions, 418 numbers, 419 symprec, 420 angle_tolerance) 421 _set_error_message() 422 423 if num_atom_prim > 0: 424 return (np.array(lattice.T, dtype='double', order='C'), 425 np.array(positions[:num_atom_prim], dtype='double', order='C'), 426 np.array(numbers[:num_atom_prim], dtype='intc')) 427 else: 428 return None 429 430def get_symmetry_from_database(hall_number): 431 """Return symmetry operations corresponding to a Hall symbol. 432 433 The Hall symbol is given by the serial number in between 1 and 530. 434 The symmetry operations are given by a dictionary whose keys are 435 'rotations' and 'translations'. 436 If it fails, None is returned. 437 """ 438 _set_no_error() 439 440 rotations = np.zeros((192, 3, 3), dtype='intc') 441 translations = np.zeros((192, 3), dtype='double') 442 num_sym = spg.symmetry_from_database(rotations, translations, hall_number) 443 _set_error_message() 444 445 if num_sym is None: 446 return None 447 else: 448 return {'rotations': 449 np.array(rotations[:num_sym], dtype='intc', order='C'), 450 'translations': 451 np.array(translations[:num_sym], dtype='double', order='C')} 452 453############ 454# k-points # 455############ 456def get_grid_point_from_address(grid_address, mesh): 457 """Return grid point index by tranlating grid address""" 458 _set_no_error() 459 460 return spg.grid_point_from_address(np.array(grid_address, dtype='intc'), 461 np.array(mesh, dtype='intc')) 462 463 464def get_ir_reciprocal_mesh(mesh, 465 cell, 466 is_shift=None, 467 is_time_reversal=True, 468 symprec=1e-5): 469 """Return k-points mesh and k-point map to the irreducible k-points. 470 471 The symmetry is serched from the input cell. 472 473 Args: 474 mesh: 475 int array (3,): Uniform sampling mesh numbers 476 cell, symprec: 477 See the docstring of get_symmetry. 478 is_shift: 479 int array (3,): [0, 0, 0] gives Gamma center mesh and value 1 gives 480 half mesh shift. 481 is_time_reversal: 482 bool: Time reversal symmetry is included or not. 483 484 Return: 485 mapping_table: 486 int array (N,): Grid point mapping table to ir-gird-points 487 grid_address: 488 int array (N, 3): Address of all grid points 489 """ 490 _set_no_error() 491 492 lattice, positions, numbers, _ = _expand_cell(cell) 493 if lattice is None: 494 return None 495 496 mapping = np.zeros(np.prod(mesh), dtype='intc') 497 grid_address = np.zeros((np.prod(mesh), 3), dtype='intc') 498 if is_shift is None: 499 is_shift = [0, 0, 0] 500 if spg.ir_reciprocal_mesh( 501 grid_address, 502 mapping, 503 np.array(mesh, dtype='intc'), 504 np.array(is_shift, dtype='intc'), 505 is_time_reversal * 1, 506 lattice, 507 positions, 508 numbers, 509 symprec) > 0: 510 return mapping, grid_address 511 else: 512 return None 513 514def get_stabilized_reciprocal_mesh(mesh, 515 rotations, 516 is_shift=None, 517 is_time_reversal=True, 518 qpoints=None): 519 """Return k-point map to the irreducible k-points and k-point grid points . 520 521 The symmetry is searched from the input rotation matrices in real space. 522 523 Args: 524 mesh: 525 int array (3,): Uniform sampling mesh numbers 526 is_shift: 527 int array (3,): [0, 0, 0] gives Gamma center mesh and value 1 gives 528 half mesh shift. 529 is_time_reversal: 530 bool: Time reversal symmetry is included or not. 531 qpoints: 532 float array (N ,3) or (3,): 533 Stabilizer(s) in the fractional coordinates. 534 535 Return: 536 mapping_table: 537 int array (N,): Grid point mapping table to ir-gird-points 538 grid_address: 539 int array (N, 3): Address of all grid points 540 """ 541 _set_no_error() 542 543 mapping_table = np.zeros(np.prod(mesh), dtype='intc') 544 grid_address = np.zeros((np.prod(mesh), 3), dtype='intc') 545 if is_shift is None: 546 is_shift = [0, 0, 0] 547 if qpoints is None: 548 qpoints = np.array([[0, 0, 0]], dtype='double', order='C') 549 else: 550 qpoints = np.array(qpoints, dtype='double', order='C') 551 if qpoints.shape == (3,): 552 qpoints = np.array([qpoints], dtype='double', order='C') 553 554 if spg.stabilized_reciprocal_mesh( 555 grid_address, 556 mapping_table, 557 np.array(mesh, dtype='intc'), 558 np.array(is_shift, dtype='intc'), 559 is_time_reversal * 1, 560 np.array(rotations, dtype='intc', order='C'), 561 qpoints) > 0: 562 return mapping_table, grid_address 563 else: 564 return None 565 566def get_grid_points_by_rotations(address_orig, 567 reciprocal_rotations, 568 mesh, 569 is_shift=np.zeros(3, dtype='intc')): 570 """Rotation operations in reciprocal space ``reciprocal_rotations`` are applied 571 to a grid point ``grid_point`` and resulting grid points are returned. 572 """ 573 _set_no_error() 574 575 rot_grid_points = np.zeros(len(reciprocal_rotations), dtype='intc') 576 spg.grid_points_by_rotations( 577 rot_grid_points, 578 np.array(address_orig, dtype='intc'), 579 np.array(reciprocal_rotations, dtype='intc', order='C'), 580 np.array(mesh, dtype='intc'), 581 np.array(is_shift, dtype='intc')) 582 583 return rot_grid_points 584 585def get_BZ_grid_points_by_rotations(address_orig, 586 reciprocal_rotations, 587 mesh, 588 bz_map, 589 is_shift=np.zeros(3, dtype='intc')): 590 """Rotation operations in reciprocal space ``reciprocal_rotations`` are applied 591 to a grid point ``grid_point`` and resulting grid points are returned. 592 """ 593 _set_no_error() 594 595 rot_grid_points = np.zeros(len(reciprocal_rotations), dtype='intc') 596 spg.BZ_grid_points_by_rotations( 597 rot_grid_points, 598 np.array(address_orig, dtype='intc'), 599 np.array(reciprocal_rotations, dtype='intc', order='C'), 600 np.array(mesh, dtype='intc'), 601 np.array(is_shift, dtype='intc'), 602 bz_map) 603 604 return rot_grid_points 605 606def relocate_BZ_grid_address(grid_address, 607 mesh, 608 reciprocal_lattice, # column vectors 609 is_shift=np.zeros(3, dtype='intc')): 610 """Grid addresses are relocated inside Brillouin zone. 611 Number of ir-grid-points inside Brillouin zone is returned. 612 It is assumed that the following arrays have the shapes of 613 bz_grid_address[prod(mesh + 1)][3] 614 bz_map[prod(mesh * 2)] 615 where grid_address[prod(mesh)][3]. 616 Each element of grid_address is mapped to each element of 617 bz_grid_address with keeping element order. bz_grid_address has 618 larger memory space to represent BZ surface even if some points 619 on a surface are translationally equivalent to the other points 620 on the other surface. Those equivalent points are added successively 621 as grid point numbers to bz_grid_address. Those added grid points 622 are stored after the address of end point of grid_address, i.e. 623 624 |-----------------array size of bz_grid_address---------------------| 625 |--grid addresses similar to grid_address--|--newly added ones--|xxx| 626 627 where xxx means the memory space that may not be used. Number of grid 628 points stored in bz_grid_address is returned. 629 bz_map is used to recover grid point index expanded to include BZ 630 surface from grid address. The grid point indices are mapped to 631 (mesh[0] * 2) x (mesh[1] * 2) x (mesh[2] * 2) space (bz_map). 632 """ 633 _set_no_error() 634 635 bz_grid_address = np.zeros( 636 ((mesh[0] + 1) * (mesh[1] + 1) * (mesh[2] + 1), 3), dtype='intc') 637 bz_map = np.zeros( 638 (2 * mesh[0]) * (2 * mesh[1]) * (2 * mesh[2]), dtype='intc') 639 num_bz_ir = spg.BZ_grid_address( 640 bz_grid_address, 641 bz_map, 642 grid_address, 643 np.array(mesh, dtype='intc'), 644 np.array(reciprocal_lattice, dtype='double', order='C'), 645 np.array(is_shift, dtype='intc')) 646 647 return bz_grid_address[:num_bz_ir], bz_map 648 649def delaunay_reduce(lattice, eps=1e-5): 650 """Run Delaunay reduction 651 652 Args: 653 lattice: Lattice parameters in the form of 654 [[a_x, a_y, a_z], 655 [b_x, b_y, b_z], 656 [c_x, c_y, c_z]] 657 symprec: 658 float: Tolerance to check if volume is close to zero or not and 659 if two basis vectors are orthogonal by the value of dot 660 product being close to zero or not. 661 662 Returns: 663 if the Delaunay reduction succeeded: 664 Reduced lattice parameters are given as a numpy 'double' array: 665 [[a_x, a_y, a_z], 666 [b_x, b_y, b_z], 667 [c_x, c_y, c_z]] 668 otherwise None is returned. 669 """ 670 _set_no_error() 671 672 delaunay_lattice = np.array(np.transpose(lattice), 673 dtype='double', order='C') 674 result = spg.delaunay_reduce(delaunay_lattice, float(eps)) 675 _set_error_message() 676 677 if result == 0: 678 return None 679 else: 680 return np.array(np.transpose(delaunay_lattice), 681 dtype='double', order='C') 682 683def niggli_reduce(lattice, eps=1e-5): 684 """Run Niggli reduction 685 686 Args: 687 lattice: Lattice parameters in the form of 688 [[a_x, a_y, a_z], 689 [b_x, b_y, b_z], 690 [c_x, c_y, c_z]] 691 eps: 692 float: Tolerance to check if difference of norms of two basis 693 vectors is close to zero or not and if two basis vectors are 694 orthogonal by the value of dot product being close to zero or 695 not. The detail is shown at 696 https://atztogo.github.io/niggli/. 697 698 Returns: 699 if the Niggli reduction succeeded: 700 Reduced lattice parameters are given as a numpy 'double' array: 701 [[a_x, a_y, a_z], 702 [b_x, b_y, b_z], 703 [c_x, c_y, c_z]] 704 otherwise None is returned. 705 """ 706 _set_no_error() 707 708 niggli_lattice = np.array(np.transpose(lattice), dtype='double', order='C') 709 result = spg.niggli_reduce(niggli_lattice, float(eps)) 710 _set_error_message() 711 712 if result == 0: 713 return None 714 else: 715 return np.array(np.transpose(niggli_lattice), dtype='double', order='C') 716 717def get_error_message(): 718 return spglib_error.message 719 720def _expand_cell(cell): 721 if isinstance(cell, tuple): 722 lattice = np.array(np.transpose(cell[0]), dtype='double', order='C') 723 positions = np.array(cell[1], dtype='double', order='C') 724 numbers = np.array(cell[2], dtype='intc') 725 if len(cell) > 3: 726 magmoms = np.array(cell[3], dtype='double') 727 else: 728 magmoms = None 729 else: 730 lattice = np.array(cell.get_cell().T, dtype='double', order='C') 731 positions = np.array(cell.get_scaled_positions(), 732 dtype='double', order='C') 733 numbers = np.array(cell.get_atomic_numbers(), dtype='intc') 734 magmoms = None 735 736 if _check(lattice, positions, numbers, magmoms): 737 return (lattice, positions, numbers, magmoms) 738 else: 739 return (None, None, None, None) 740 741def _check(lattice, positions, numbers, magmoms): 742 if lattice.shape != (3, 3): 743 return False 744 if positions.ndim != 2: 745 return False 746 if positions.shape[1] != 3: 747 return False 748 if numbers.ndim != 1: 749 return False 750 if len(numbers) != positions.shape[0]: 751 return False 752 if magmoms is not None: 753 if magmoms.ndim != 1: 754 return False 755 if len(magmoms) != len(numbers): 756 return False 757 return True 758 759def _set_error_message(): 760 spglib_error.message = spg.error_message() 761 762def _set_no_error(): 763 spglib_error.message = "no error" 764