# -*- coding: utf-8 -*-
#
# The high-throughput toolkit (httk)
# Copyright (C) 2012-2015 Rickard Armiento
# Some parts imported from cif2cell, (C) Torbjörn Björkman
#
# This program is free software: you can redistribute it and/or modify
# it under the terms of the GNU Affero General Public License as
# published by the Free Software Foundation, either version 3 of the
# License, or (at your option) any later version.
#
# This program is distributed in the hope that it will be useful,
# but WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
# GNU Affero General Public License for more details.
#
# You should have received a copy of the GNU Affero General Public License
# along with this program. If not, see <http://www.gnu.org/licenses/>.
import sys
from httk.core import is_sequence, breath_first_idxs, FracVector, FracScalar, MutableFracVector
from httk.atomistic.cell import Cell
from httk.atomistic.unitcellsites import UnitcellSites
from httk.atomistic import spacegrouputils
from math import sqrt, acos, cos, sin, pi
from httk.atomistic.data import periodictable
from fractions import Fraction
from httk.atomistic.spacegrouputils import crystal_system_from_hall
[docs]def sort_coordgroups(coordgroups, individual_data):
counts = [len(x) for x in coordgroups]
newcoordgroups = []
newindividual_data = []
for group in range(len(counts)):
order = sorted(range(counts[group]), key=lambda x: (coordgroups[group][x][0], coordgroups[group][x][1], coordgroups[group][x][2]))
newcoordgroups.append(coordgroups[group][order])
if individual_data is not None:
newindividual_data.append([individual_data[group][i] for i in order])
if individual_data is None:
return coordgroups.stack(newcoordgroups), newindividual_data
else:
return coordgroups.stack(newcoordgroups), None
[docs]def niggli_to_metric(niggli):
m = niggli.noms
# Since the niggli matrix contains 2*the product of the off diagonal elements, we increase the denominator by cell.denom*2
return FracVector(((2*m[0][0], m[1][2], m[1][1]), (m[1][2], 2*m[0][1], m[1][0]), (m[1][1], m[1][0], 2*m[0][2])), niggli.denom*2).simplify()
[docs]def metric_to_niggli(cell):
m = cell.noms
return FracVector(((m[0][0], m[1][1], m[2][2]), (2*m[1][2], 2*m[0][2], 2*m[0][1])), cell.denom).simplify()
[docs]def coords_and_occupancies_to_coordgroups_and_assignments(coords, occupancies):
if len(occupancies) != len(coords):
raise Exception("Occupations and coords needs to be of the same length.")
if len(occupancies) == 0:
return [], FracVector((), 1)
coordgroups = []
group_occupancies = []
for i in range(len(occupancies)):
try:
idx = group_occupancies.index(occupancies[i])
except ValueError:
idx = len(group_occupancies)
group_occupancies.append(occupancies[i])
coordgroups.append([])
coordgroups[idx].append(coords[i])
return coordgroups, group_occupancies
[docs]def coords_to_coordgroups(coords, counts):
coordgroups = []
idx = 0
for count in counts:
coordgroups.append(coords[idx:count+idx])
idx += count
return coordgroups
[docs]def coordgroups_to_coords(coordgroups):
coords = []
counts = [len(x) for x in coordgroups]
for group in coordgroups:
for i in range(len(group)):
coords.append(group[i])
coords = FracVector.stack_vecs(coords)
return coords, counts
[docs]def niggli_to_cell_old(niggli_matrix, orientation=1):
cell = FracVector.use(niggli_matrix)
niggli_matrix = niggli_matrix.to_floats()
s11, s22, s33 = niggli_matrix[0][0], niggli_matrix[0][1], niggli_matrix[0][2]
s23, s13, s12 = niggli_matrix[1][0]/2.0, niggli_matrix[1][1]/2.0, niggli_matrix[1][2]/2.0
a, b, c = sqrt(s11), sqrt(s22), sqrt(s33)
alpha_rad, beta_rad, gamma_rad = acos(s23/(b*c)), acos(s13/(c*a)), acos(s12/(a*b))
iv = 1-cos(alpha_rad)**2-cos(beta_rad)**2-cos(gamma_rad)**2+2*cos(alpha_rad)*cos(beta_rad)*cos(gamma_rad)
# Handle that iv may be very, very slightly < 0 by the floating point accuracy limit
if iv > 0:
v = sqrt(iv)
else:
v = 0.0
if c*v < 1e-14:
raise Exception("niggli_to_cell: Physically unreasonable cell, cell vectors degenerate or very close to degenerate.")
if orientation < 0:
cell = [[-a, 0.0, 0.0],
[-b*cos(gamma_rad), -b*sin(gamma_rad), 0.0],
[-c*cos(beta_rad), -c*(cos(alpha_rad)-cos(beta_rad)*cos(gamma_rad))/sin(gamma_rad), -c*v/sin(gamma_rad)]]
else:
cell = [[a, 0.0, 0.0],
[b*cos(gamma_rad), b*sin(gamma_rad), 0.0],
[c*cos(beta_rad), c*(cos(alpha_rad)-cos(beta_rad)*cos(gamma_rad))/sin(gamma_rad), c*v/sin(gamma_rad)]]
for i in range(3):
for j in range(3):
cell[i][j] = round(cell[i][j], 14)
return cell
[docs]def niggli_to_basis(niggli_matrix, orientation=1):
#cell = FracVector.use(niggli_matrix)
niggli_matrix = niggli_matrix.to_floats()
s11, s22, s33 = niggli_matrix[0][0], niggli_matrix[0][1], niggli_matrix[0][2]
s23, s13, s12 = niggli_matrix[1][0]/2.0, niggli_matrix[1][1]/2.0, niggli_matrix[1][2]/2.0
a, b, c = sqrt(s11), sqrt(s22), sqrt(s33)
alpha_rad, beta_rad, gamma_rad = acos(s23/(b*c)), acos(s13/(c*a)), acos(s12/(a*b))
iv = 1-cos(alpha_rad)**2-cos(beta_rad)**2-cos(gamma_rad)**2+2*cos(alpha_rad)*cos(beta_rad)*cos(gamma_rad)
# Handle that iv may be very, very slightly < 0 by the floating point accuracy limit
if iv > 0:
v = sqrt(iv)
else:
v = 0.0
if c*v < 1e-14:
raise Exception("niggli_to_cell: Physically unreasonable cell, cell vectors degenerate or very close to degenerate.")
if orientation < 0:
basis = [[-a, 0.0, 0.0],
[-b*cos(gamma_rad), -b*sin(gamma_rad), 0.0],
[-c*cos(beta_rad), -c*(cos(alpha_rad)-cos(beta_rad)*cos(gamma_rad))/sin(gamma_rad), -c*v/sin(gamma_rad)]]
else:
basis = [[a, 0.0, 0.0],
[b*cos(gamma_rad), b*sin(gamma_rad), 0.0],
[c*cos(beta_rad), c*(cos(alpha_rad)-cos(beta_rad)*cos(gamma_rad))/sin(gamma_rad), c*v/sin(gamma_rad)]]
for i in range(3):
for j in range(3):
basis[i][j] = round(basis[i][j], 14)
return basis
[docs]def basis_to_niggli(basis):
basis = FracVector.use(basis)
A = basis.noms
det = basis.det()
if det == 0:
raise Exception("basis_to_niggli: singular cell matrix.")
if det > 0:
orientation = 1
else:
orientation = -1
s11 = A[0][0]*A[0][0]+A[0][1]*A[0][1]+A[0][2]*A[0][2]
s22 = A[1][0]*A[1][0]+A[1][1]*A[1][1]+A[1][2]*A[1][2]
s33 = A[2][0]*A[2][0]+A[2][1]*A[2][1]+A[2][2]*A[2][2]
s23 = A[1][0]*A[2][0]+A[1][1]*A[2][1]+A[1][2]*A[2][2]
s13 = A[0][0]*A[2][0]+A[0][1]*A[2][1]+A[0][2]*A[2][2]
s12 = A[0][0]*A[1][0]+A[0][1]*A[1][1]+A[0][2]*A[1][2]
new = FracVector.create(((s11, s22, s33), (2*s23, 2*s13, 2*s12)), denom=basis.denom**2).simplify()
return new, orientation
[docs]def basis_determinant(basis):
try:
return basis.det()
except AttributeError:
C = basis
det = C[0][0]*C[1][1]*C[2][2] + C[0][1]*C[1][2]*C[2][0] + C[0][2]*C[1][0]*C[2][1] - C[0][2]*C[1][1]*C[2][0] - C[0][1]*C[1][0]*C[2][2] - C[0][0]*C[1][2]*C[2][1]
return det
[docs]def basis_vol_to_scale(basis, vol):
det = float(basis_determinant(basis))
if abs(det) < 1e-12:
raise Exception("basis_to_niggli: Too singular cell matrix.")
return (float(vol)/(abs(det)))**(1.0/3.0)
[docs]def basis_scale_to_vol(basis, scale):
det = float(basis_determinant(basis))
if abs(det) < 1e-12:
raise Exception("basis_to_niggli: Too singular cell matrix.")
return (((scale)**3)*(abs(det)))
[docs]def niggli_vol_to_scale(niggli_matrix, vol):
niggli_matrix = FracVector.use(niggli_matrix)
metric = niggli_to_metric(niggli_matrix)
volsqr = metric.det()
if volsqr == 0:
raise Exception("niggli_vol_to_scale: singular cell matrix.")
det = sqrt(float(volsqr))
return (float(vol)/det)**(1.0/3.0)
[docs]def niggli_scale_to_vol(niggli_matrix, scale):
niggli_matrix = FracVector.use(niggli_matrix)
metric = niggli_to_metric(niggli_matrix)
volsqr = metric.det()
if volsqr == 0:
raise Exception("niggli_vol_to_scale: singular cell matrix.")
det = sqrt(float(volsqr))
if abs(det) < 1e-12:
raise Exception("niggli_scale_to_vol: singular cell matrix.")
return (((scale)**3)*det)
[docs]def niggli_to_lengths_angles(niggli_matrix):
s11, s22, s33 = niggli_matrix[0][0], niggli_matrix[0][1], niggli_matrix[0][2]
s23, s13, s12 = niggli_matrix[1][0]/2.0, niggli_matrix[1][1]/2.0, niggli_matrix[1][2]/2.0
a, b, c = sqrt(s11), sqrt(s22), sqrt(s33)
alpha, beta, gamma = acos(s23/(b*c))*180/pi, acos(s13/(c*a))*180/pi, acos(s12/(a*b))*180/pi
return [a, b, c], [alpha, beta, gamma]
[docs]def lengths_angles_to_niggli(lengths, angles):
niggli_matrix = [[None, None, None], [None, None, None]]
niggli_matrix[0][0] = lengths[0]*lengths[0]
niggli_matrix[0][1] = lengths[1]*lengths[1]
niggli_matrix[0][2] = lengths[2]*lengths[2]
niggli_matrix[1][0] = 2.0*lengths[1]*lengths[2]*cos(angles[0]*pi/180)
niggli_matrix[1][1] = 2.0*lengths[0]*lengths[2]*cos(angles[1]*pi/180)
niggli_matrix[1][2] = 2.0*lengths[0]*lengths[1]*cos(angles[2]*pi/180)
return niggli_matrix
[docs]def cartesian_to_reduced(cell, coordgroups):
cell = FracVector.use(cell)
cellinv = cell.inv()
newcoordgroups = []
for coordgroup in coordgroups:
newcoordgroup = coordgroup*cellinv
newcoordgroups.append(newcoordgroup)
return newcoordgroups
[docs]def structure_to_p1structure(struct, backends=['ase']):
for backend in backends:
if backend == 'ase':
try:
from httk.external import ase_ext
return ase_ext.structure_to_p1structure(struct)
except ImportError:
raise
pass
raise Exception("structure_to_p1structure: None of the available backends available.")
[docs]def structure_to_sgstructure(struct, backends=['platon']):
for backend in backends:
if backend == 'platon':
try:
from httk.external import platon
return platon.structure_to_sgstructure(struct)
except ImportError:
pass
raise Exception("structure_to_sgstructure: None of the available backends available.")
[docs]def reduced_to_cartesian(cell, coordgroups):
cell = FracVector.use(cell)
newcoordgroups = []
for coordgroup in coordgroups:
newcoordgroup = coordgroup*cell
newcoordgroups.append(newcoordgroup)
return newcoordgroups
[docs]def abstract_symbol(count):
if count < 27:
return chr(64+count)
my = (count-1) % 26
rec = (count-1) // 26
return abstract_symbol(rec)+chr(97+my)
# def parse_assignment(assignment):
# if isinstance(assignment,(tuple,list)):
# if len(assignment)>=3:
# return FracVector.use(assignment[2])
# return FracVector.use(assignment[1])
# return FracVector.create(1)
#def assignments_tidy(assignments):
# newassignments=[]
# for assignment in assignments:
# newassignments.append((periodictable.atomic_symbol(assignment),periodictable.atomic_number(assignment),parse_assignment(assignment),))
# return tuple(newassignments)
[docs]def coordgroups_cartesian_to_reduced(coordgroups, basis):
basis = FracVector.use(basis)
cellinv = basis.inv()
newcoordgroups = []
for coordgroup in coordgroups:
newcoordgroup = coordgroup*cellinv
newcoordgroups.append(newcoordgroup)
return newcoordgroups
[docs]def coords_and_counts_to_coordgroups(coords, counts):
coordgroups = []
idx = 0
for count in counts:
coordgroups.append(coords[idx:count+idx])
idx += count
return coordgroups
[docs]def coordgroups_and_assignments_to_coords_and_occupancies(coordgroups, assignments):
coords = coordgroups_to_coords(coordgroups)
occupancies = []
for i in range(len(coordgroups)):
occupancies += [assignments[i]]*len(coordgroups[i])
return coords, occupancies
[docs]def structure_reduced_uc_to_representative(struct, backends=['isotropy', 'fake']):
for backend in backends:
if backend == 'isotropy':
try:
from httk.external import isotropy_ext
sys.stderr.write("Warning: need to run symmetry finder. This may take a while.\n")
struct = isotropy_ext.struct_process_with_isotropy(struct)
return struct
except ImportError:
pass
if backend == 'fake':
try:
from httk.atomistic import Structure
sys.stderr.write("Warning: using 'fake' symmetry finder that never finds any symmetry.\n")
newstruct = Structure.create(
assignments=struct.assignments,
rc_cell=struct.uc_cell,
rc_reduced_coords=struct.uc_reduced_coords,
rc_counts=struct.uc_counts,
uc_cell=struct.uc_cell,
uc_reduced_coords=struct.uc_reduced_coords,
uc_counts=struct.uc_counts,
spacegroup='P 1',
tags=struct.get_tags(), refs=struct.get_refs(), periodicity=struct.uc_sites.pbc)
return newstruct
except ImportError:
raise
pass
raise Exception("structure_to_sgstructure: None of the available backends available.")
[docs]def coordgroups_reduced_uc_to_representative(coordgroups, basis, backends=['isotropy']):
sys.stderr.write("WARNING: coordgroups_reduced_uc_to_representative: running untested code...\n")
for backend in backends:
if backend == 'isotropy':
try:
from httk.external import isotropy_ext
struct = isotropy_ext.uc_reduced_coordgroups_process_with_isotropy(coordgroups, basis)
return struct
except ImportError:
raise
pass
#if backend == 'platon':
# try:
# from httk.external import platon_ext
# return platon_ext.coordgroups_reduced_uc_to_representative(coordgroups, basis, hall_symbol)
# except ImportError:
# raise
# pass
raise Exception("structure_to_sgstructure: None of the available backends available.")
[docs]def internal_coordgroups_reduced_rc_to_unitcellsites(coordgroups, basis, hall_symbol, eps=0.001):
symops = spacegrouputils.get_symops(hall_symbol)
newcoordgroups = []
for coordgroup in coordgroups:
newcoordgroup = []
for symop in symops:
rotcoords = coordgroup*(symop[0].T())
for coord in rotcoords:
finalcoord = (coord+symop[1]).normalize()
if finalcoord not in newcoordgroup:
for checkcoord in newcoordgroup:
if (checkcoord-finalcoord).normalize_half().lengthsqr() < eps:
break
else:
newcoordgroup += [finalcoord]
newcoordgroup = sorted(newcoordgroup, key=lambda x: (x[0], x[1], x[2]))
newcoordgroups += [newcoordgroup]
return newcoordgroups, basis
[docs]def coordgroups_reduced_rc_to_unitcellsites(coordgroups, basis, hall_symbol, backends=['cif2cell', 'internal', 'ase']):
# TODO: Make own, or use cif2cell to generate reduced unitcell
if hall_symbol == 'P 1':
return UnitcellSites.create(reduced_coordgroups=coordgroups), Cell.create(basis)
for backend in backends:
if backend == 'internal':
newcoordgroups, newcell = internal_coordgroups_reduced_rc_to_unitcellsites(coordgroups, basis, hall_symbol)
return UnitcellSites.create(reduced_coordgroups=newcoordgroups), Cell.create(basis)
if backend == 'cif2cell':
try:
from httk.external import cif2cell_ext
return cif2cell_ext.coordgroups_reduced_rc_to_unitcellsites(coordgroups, basis, hall_symbol)
except ImportError:
pass
if backend == 'ase':
try:
from httk.external import ase_glue
return ase_glue.coordgroups_reduced_rc_to_unitcellsites(coordgroups, basis, hall_symbol)
except ImportError:
raise
pass
raise Exception("structure_to_p1structure: None of the available backends available.")
[docs]def coordswap(fromidx, toidx, cell, coordgroups):
new_coordgroups = []
for group in coordgroups:
coords = MutableFracVector.from_FracVector(group)
rows = coords[:, toidx]
coords[:, toidx] = coords[:, fromidx]
coords[:, fromidx] = rows
new_coordgroups.append(coords.to_FracVector())
coordgroups = FracVector.create(new_coordgroups)
cell = MutableFracVector.from_FracVector(cell)
row = cell[toidx]
cell[toidx] = cell[fromidx]
cell[fromidx] = row
cell = cell.to_FracVector()
return (cell, coordgroups)
[docs]def clean_coordgroups_and_assignments(coordgroups, assignments):
if len(assignments) != len(coordgroups):
raise Exception("Occupations and coords needs to be of the same length.")
if len(assignments) == 0:
return [], FracVector((), 1)
new_coordgroups = []
new_assignments = []
for i in range(len(assignments)):
for j in range(len(new_assignments)):
if assignments[i] == new_assignments[j]:
idx = j
new_coordgroups[idx] = FracVector.chain_vecs([new_coordgroups[idx], coordgroups[i]])
break
else:
new_coordgroups.append(coordgroups[i])
new_assignments.append(assignments[i])
idx = len(coordgroups)-1
return new_coordgroups, new_assignments
[docs]def occupations_and_coords_to_assignments_and_coordgroups(occupationscoords, occupations):
if len(occupationscoords) == 0:
return [], FracVector((), 1)
occupationscoords = FracVector.use(occupationscoords)
new_coordgroups = []
new_assignments = []
for i in range(len(occupations)):
for j in range(len(new_assignments)):
if occupations[i] == new_assignments[j]:
new_coordgroups[j].append(occupationscoords[i])
break
else:
new_coordgroups.append([occupationscoords[i]])
new_assignments.append(occupations[i])
new_coordgroups = FracVector.create(new_coordgroups)
return new_assignments, new_coordgroups
[docs]def coordgroups_and_assignments_to_symbols(coordgroups, assignmentobj):
"""
Return a list of atomic symbols, repeated as needed
"""
symbols = []
for i in range(len(coordgroups)):
name = assignmentobj.symbols[i]
symbols += [name]*len(coordgroups[i])
return symbols
[docs]def structure_tidy(struct, backends=['platon']):
for backend in backends:
if backend == 'platon':
try:
from httk.external import platon_ext
return platon_ext.structure_tidy(struct)
except ImportError:
raise
pass
raise Exception("structure_tidy: None of the available backends available.")
# def get_primitive_basis_transform(hall_symbol, niggli_matrix):
# half = Fraction(1, 2)
# lattice_symbol = hall_symbol.lstrip("-")[0][0]
# crystal_system = crystal_system_from_hall(hall_symbol)
#
# lattrans = None
# unit = FracVector.create([[1, 0, 0],
# [0, 1, 0],
# [0, 0, 1]])
#
# if crystal_system == 'cubic':
# if lattice_symbol == 'P':
# lattrans = unit
# elif lattice_symbol == 'F':
# lattrans = FracVector.create([[half, half, 0],
# [half, 0, half],
# [0, half, half]])
# elif lattice_symbol == 'I':
# lattrans = FracVector.create([[half, half, half],
# [-half, half, half],
# [-half, -half, half]])
# elif (crystal_system == 'hexagonal' or crystal_system == 'trigonal'):
# if lattice_symbol == 'P':
# # Conventional cell should already be primitive hexagonal one
# lattrans = unit
# elif lattice_symbol == 'R':
# # Conventional cell should already be primitive rhombohedral one
# lattrans = unit
#
# elif crystal_system == 'tetragonal':
# if lattice_symbol == 'P':
# lattrans = unit
# elif lattice_symbol == 'I':
# lattrans = FracVector.create([[half, -half, half],
# [half, half, half],
# [-half, -half, half]])
# elif crystal_system == 'orthorhombic':
# if lattice_symbol == 'P':
# lattrans = unit
# elif lattice_symbol == 'B':
# lattrans = FracVector.create([[half, half, 0],
# [-half, half, 0],
# [0, 0, 1]])
# elif lattice_symbol == 'F' or lattice_symbol == 'A' or lattice_symbol == 'B' or lattice_symbol == 'C':
# lattrans = FracVector.create([[half, 0, half],
# [half, half, 0],
# [0, half, half]])
# elif lattice_symbol == 'I':
# lattrans = FracVector.create([[half, half, half],
# [-half, half, half],
# [-half, -half, half]])
#
# elif crystal_system == 'monoclinic':
# if lattice_symbol == 'P':
# lattrans = unit
# elif lattice_symbol == 'B':
# lattrans = FracVector.create([[half, 0, -half],
# [0, 1, 0],
# [half, 0, half]])
#
# elif crystal_system == 'triclinic':
# lattrans = unit
#
# else:
# raise Exception("structureutils.get_primitive_basis_transform: unknown crystal system, "+str(crystal_system))
#
# if lattrans is None:
# raise Exception("structureutils.get_primitive_basis_transform: no match for lattice transform.")
#
# return lattrans
# Imported from cif2cell by Torbjörn Björkman, uctools.py and heavily modified
# def get_primitive_basis_transform(basis, hall_symbol, eps=0.001):
# basis = FracVector.use(basis)
# niggli, orientation = basis_to_niggli(basis)
# length, angles = niggli_to_lengths_angles(niggli)
# a, b, c = length
# alpha, beta, gamma = angles
#
# zero = Fraction(0)
# half = Fraction(1, 2)
# one = Fraction(1)
# third = Fraction(1, 3)
# lattice_symbol = hall_symbol.lstrip("-")[0][0]
# system = crystal_system_from_hall(hall_symbol)
#
# if lattice_symbol == 'I':
# # Body centered
# #transvecs = [FracVector.create([zero, zero, zero]),
# # FracVector.create([half, half, half])]
# if system == 'cubic':
# lattrans = FracVector.create([[-half, half, half],
# [half, -half, half],
# [half, half, -half]])
# else:
# lattrans = FracVector.create([[one, zero, zero],
# [zero, one, zero],
# [half, half, half]])
# elif lattice_symbol == 'F':
# # face centered
# #transvecs = [FracVector.create([zero, zero, zero]),
# # FracVector.create([half, half, zero]),
# # FracVector.create([half, zero, half]),
# # FracVector.create([zero, half, half])]
# lattrans = FracVector.create([[half, half, zero],
# [half, zero, half],
# [zero, half, half]])
# elif lattice_symbol == 'A':
# # A-centered
# #transvecs = [FracVector.create([zero, zero, zero]),
# # FracVector.create([zero, half, half])]
# lattrans = FracVector.create([[one, zero, zero],
# [zero, half, -half],
# [zero, half, half]])
# elif lattice_symbol == 'B':
# # B-centered
# #transvecs = [FracVector.create([zero, zero, zero]),
# # FracVector.create([half, zero, half])]
# lattrans = FracVector.create([[half, zero, -half],
# [zero, one, zero],
# [half, zero, half]])
# elif lattice_symbol == 'C':
# # C-centered
# #transvecs = [FracVector.create([zero, zero, zero]),
# # FracVector.create([half, half, zero])]
# lattrans = FracVector.create([[half, -half, zero],
# [half, half, zero],
# [zero, zero, one]])
# #elif hall_symbol in Hex2RhombHall or hall_symbol in Rhomb2HexHall:
# #elif hall_symbol in Hex2RhombHall or hall_symbol in Rhomb2HexHall:
# #if hall_symbol in Rhomb2HexHall and abs(gamma - 120) > eps:
# # rhomb2hextrans = FracVector.create([[2*third, third, third],
# # [-third, third, third],
# # [-third, -2*third, third]])
# #else:
# # rhomb2hextrans = 1
# #if abs(gamma - 120) < eps:
# # # rhombohedral from hexagonal setting
# # # rhombohedral = True
# # # transvecs = [FracVector.create([zero, zero, zero]),
# # # FracVector.create([third, 2*third, 2*third]),
# # # FracVector.create([2*third, third, third])]
# # lattrans = FracVector.create([[2*third, third, third],
# # [-third, third, third],
# # [-third, -2*third, third]])
# #else:
# # #transvecs = [FracVector.create([zero, zero, zero])]
# # lattrans = FracVector.create([[1, 0, 0],
# # [0, 1, 0],
# # [0, 0, 1]])
# #lattrans = lattrans*rhomb2hextrans
# else:
# #transvecs = [FracVector.create([zero, zero, zero])]
# lattrans = FracVector.create([[1, 0, 0],
# [0, 1, 0],
# [0, 0, 1]])
# # Transform to primitive cell
# return lattrans
[docs]def main():
cell = FracVector.create([[1, 1, 0], [1, 0, 1], [0, 1, 1]])
coordgroups = FracVector.create([[[2, 3, 5], [3, 5, 4]], [[4, 6, 7]]])
assignments = [2, 5]
print(cell, coordgroups)
cell, coordgroups = coordswap(0, 2, cell, coordgroups)
print(cell, coordgroups)
pass
if __name__ == "__main__":
main()