mirror of https://github.com/abinit/abipy.git
Docs + code cleanup
This commit is contained in:
Matteo Giantomassi
parent
5c614179e7
commit
352a781a0d
|
|
@ -8,214 +8,11 @@ import numpy as np
|
|||
import sympy as sp
|
||||
|
||||
#from monty.functools import lazy_property
|
||||
from collections import Sequence, OrderedDict
|
||||
from collections import OrderedDict
|
||||
from monty.termcolor import cprint
|
||||
from abipy.core.mixins import Has_Structure
|
||||
|
||||
|
||||
def _validate_params(params, wpos):
|
||||
"""
|
||||
Check keys in params. Raises ValueError.
|
||||
"""
|
||||
check = params.copy()
|
||||
for k in params:
|
||||
if k not in wpos.names:
|
||||
raise ValueError("Cannot find key: %s in symbol_names: %s" % (k, wpos.names))
|
||||
check.pop(k)
|
||||
if check:
|
||||
raise ValueError("Found unknown symbol names in %s" % check)
|
||||
|
||||
|
||||
class WyckoffPositions(Sequence):
|
||||
"""
|
||||
A numpy array with sympy expression.
|
||||
"""
|
||||
@classmethod
|
||||
def from_string(cls, s):
|
||||
"""
|
||||
Initialize the object from string `s`.
|
||||
"""
|
||||
# Build a matrix of strings with expressions.
|
||||
s = s.replace("\n", "")
|
||||
s = s.replace("(", "").replace(")", "")
|
||||
s = s.replace("[", "").replace("]", "")
|
||||
|
||||
try:
|
||||
matrix_str = np.reshape(s.split(","), (-1, 3))
|
||||
except ValueError:
|
||||
raise ValueError("Wrong string in input, perhaps missing comma. Input string:\n %s" % str(s))
|
||||
|
||||
return cls(matrix_str)
|
||||
|
||||
def __init__(self, matrix_str, letter=None, site_symmetry=None):
|
||||
"""
|
||||
|
||||
Args:
|
||||
matrix_str:
|
||||
letter:
|
||||
site_symmetry:
|
||||
"""
|
||||
self.letter, self.site_symmetry = None, None
|
||||
|
||||
# Build numpy matrix of sympy expressions from the matrix of strings.
|
||||
from sympy import Symbol
|
||||
from sympy.parsing.sympy_parser import parse_expr
|
||||
exprs = []
|
||||
for vec in np.reshape(matrix_str, (-1, 3)):
|
||||
exprs.extend([parse_expr(e) for e in vec])
|
||||
|
||||
self._expr_mat = np.reshape(exprs, (-1, 3))
|
||||
|
||||
# Get set of sympy symbols.
|
||||
symbols = []
|
||||
for expr in self.ravel():
|
||||
for arg in expr.args:
|
||||
if isinstance(arg, Symbol): symbols.append(arg)
|
||||
|
||||
# SymPy symbols and names.
|
||||
self.symbols = list(set(symbols))
|
||||
self.names = [s.name for s in self.symbols]
|
||||
|
||||
def __len__(self):
|
||||
return len(self._expr_mat)
|
||||
|
||||
def __getitem__(self, key):
|
||||
return self._expr_mat.__getitem__(key)
|
||||
|
||||
def __str__(self):
|
||||
"""String representation."""
|
||||
lines = []
|
||||
for vec in self:
|
||||
lines.append("(" + ",".join(str(e) for e in vec) + ")")
|
||||
return ", ".join(lines)
|
||||
|
||||
def ravel(self):
|
||||
"""Flat list, similar to np.ravel"""
|
||||
return self._expr_mat.ravel()
|
||||
|
||||
@property
|
||||
def mult(self):
|
||||
"""Multiplicity"""
|
||||
return len(self)
|
||||
|
||||
@property
|
||||
def num_params(self):
|
||||
"""The number of parameters."""
|
||||
return len(self.symbols)
|
||||
|
||||
def params_from_frac_coords(self, frac_coords):
|
||||
"""
|
||||
Compute the Wyckoff parameters from the fractional coordinates.
|
||||
"""
|
||||
frac_coords = np.reshape(frac_coords, (-1, 3))
|
||||
assert len(self) == len(frac_coords)
|
||||
equations = []
|
||||
for i in range(len(self)):
|
||||
eqs = self[i, :] - frac_coords[i, :]
|
||||
equations.extend(eqs)
|
||||
|
||||
# Solve the problem.
|
||||
from sympy.solvers import solve
|
||||
d = solve(equations, self.symbols)
|
||||
#print("solution: ", d)
|
||||
|
||||
# Return results in a dict with symbol names.
|
||||
if not d:
|
||||
return None
|
||||
|
||||
return {symbol.name: d[symbol] for symbol in self.symbols}
|
||||
|
||||
def error_of_params(self, params, frac_coords):
|
||||
"""
|
||||
Return a numpy array with the difference between the
|
||||
wyckoff positions computed from the sympy expressions and the input frac_coords
|
||||
"""
|
||||
_validate_params(params, self)
|
||||
frac_coords = np.reshape(frac_coords, (-1, 3))
|
||||
assert len(self) == len(frac_coords)
|
||||
|
||||
frac_errors = np.empty((self.mult, 3))
|
||||
for i in range(len(self)):
|
||||
for j in range(3):
|
||||
frac_errors[i,j] = self[i,j].subs(params) - frac_coords[i,j]
|
||||
|
||||
# Modulo lattice vector (with tolerance?)
|
||||
return frac_errors % 1
|
||||
|
||||
|
||||
class Wyckoff(object):
|
||||
"""
|
||||
"""
|
||||
#@classmethod
|
||||
#def from_site2wpos_ordict(cls, site2wpos):
|
||||
# return cls(site2wpos)
|
||||
|
||||
def __init__(self, site2wpos):
|
||||
# TODO: Ordered dict.
|
||||
d = {}
|
||||
for k, s in site2wpos.items():
|
||||
d[k] = WyckoffPositions.from_string(s)
|
||||
|
||||
self.site2wpos = d
|
||||
|
||||
def __str__(self):
|
||||
"""String representation."""
|
||||
lines = []
|
||||
app = lines.append
|
||||
app("Wyckoff positions:")
|
||||
for site, wpos in self.site2wpos.items():
|
||||
app("%s site:\n\t%s" % (site, wpos))
|
||||
|
||||
return "\n".join(lines)
|
||||
|
||||
def generate_structure(self, lattice, site2params, to_unit_cell=False, struct_cls=None):
|
||||
"""
|
||||
Generate structure object from `lattice` and dictionary `site2params`
|
||||
"""
|
||||
if not isinstance(site2params, OrderedDict):
|
||||
raise ValueError("Please use a OrderedDict for site2params so that\n" +
|
||||
"the algorithm used to build the structure is deterministic.")
|
||||
|
||||
# Build species and coords arrays.
|
||||
species, coords = [], []
|
||||
for elsym, wpos in self.site2wpos.items():
|
||||
params = site2params[elsym]
|
||||
_validate_params(params, wpos)
|
||||
for vec_expr in wpos:
|
||||
species.append(elsym)
|
||||
# Evaluate sympy expression with params.
|
||||
coords.append([float(expr.subs(params)) for expr in vec_expr])
|
||||
|
||||
#for sp, c in zip(species, coords):
|
||||
# print(sp, c, type(c), type(c[0]))
|
||||
|
||||
from abipy.core.structure import Structure
|
||||
cls = struct_cls if struct_cls is not None else Structure
|
||||
return cls(lattice, species, coords, validate_proximity=True,
|
||||
to_unit_cell=to_unit_cell, coords_are_cartesian=False, site_properties=None)
|
||||
|
||||
def find_params(self, structure):
|
||||
"""
|
||||
Compute the value of the parameter given a Structure object.
|
||||
"""
|
||||
frac_coords = structure.get_symbol2coords()
|
||||
|
||||
# Solve the problem.
|
||||
params = {}
|
||||
for elsym, wpos in self.site2wpos.items():
|
||||
params[elsym] = wpos.params_from_frac_coords(frac_coords[elsym])
|
||||
|
||||
# Return results in a dict with symbol names.
|
||||
return params
|
||||
|
||||
def error_of_params(self, params, structure):
|
||||
frac_coords = structure.get_symbol2coords()
|
||||
|
||||
for elsym, wpos in self.site2wpos.items():
|
||||
frac_errors = wpos.error_of_params(params[elsym], frac_coords[elsym])
|
||||
print(frac_errors)
|
||||
|
||||
|
||||
class SiteSymmetries(Has_Structure):
|
||||
|
||||
def __init__(self, structure):
|
||||
|
|
@ -265,13 +62,27 @@ class SiteSymmetries(Has_Structure):
|
|||
if len(rotations) != 0:
|
||||
herm_symbol, ptg_num, trans_mat = spglib.get_pointgroup(rotations)
|
||||
|
||||
self.sitesym_labels.append("%s (#%d) ord: %d" % (herm_symbol.strip(), ptg_num, len(rotations)))
|
||||
self.sitesym_labels.append("%s (#%d,nsym:%d)" % (herm_symbol.strip(), ptg_num, len(rotations)))
|
||||
|
||||
#for irred_isite in self.irred_isites:
|
||||
# for isite_eq, rm1, tau, l0 in self.eq_sites[irred_isite]:
|
||||
# # isite_eq = rm1(irred_site - tau) + l0
|
||||
|
||||
@property
|
||||
def structure(self):
|
||||
"""|Structure| object."""
|
||||
return self._structure
|
||||
|
||||
def __str__(self):
|
||||
return self.to_string()
|
||||
|
||||
def to_string(self, verbose=0):
|
||||
"""String representation with verbosity level verbose."""
|
||||
lines = []; app = lines.append
|
||||
app(self.structure.to_string(verbose=verbose, title="Structure"))
|
||||
|
||||
return "\n".join(lines)
|
||||
|
||||
def get_wyckoff_dataframe(self, view="all", select_symbols=None, decimals=5, verbose=0):
|
||||
"""
|
||||
Find Wyckoff positions.
|
||||
|
|
@ -304,7 +115,13 @@ class SiteSymmetries(Has_Structure):
|
|||
l0 = indsym[iatom, isym, :3]
|
||||
l0 = sp.Matrix(l0)
|
||||
m = rm1 * (vector - tau) - (l0 + vector)
|
||||
if verbose: sp.pprint(m)
|
||||
if verbose:
|
||||
print(92 * "=")
|
||||
print("System of linear equations for iatom %d, %s" % (iatom, repr(site)))
|
||||
sp.pprint(m)
|
||||
print("Rotation:")
|
||||
sp.pprint(rm1)
|
||||
|
||||
system.append(m)
|
||||
|
||||
# Solve system of linear equations.
|
||||
|
|
@ -322,17 +139,17 @@ class SiteSymmetries(Has_Structure):
|
|||
d["wyckoff"] = wlabel
|
||||
d["site_symmetry"] = sitesym_labels[iatom]
|
||||
for s in frac_symbols:
|
||||
d[s] = s
|
||||
d[str(s)] = str(s)
|
||||
if len(solutions) > 1:
|
||||
cprint("Found multiple solutions for iatom %d" % iatom, "red")
|
||||
d.update(**solutions[0])
|
||||
d.update({str(k): str(v) for k, v in solutions[0].items()})
|
||||
rows.append(d)
|
||||
|
||||
import pandas as pd
|
||||
df = pd.DataFrame(rows, index=None, columns=list(rows[0].keys()) if rows else None)
|
||||
return df
|
||||
|
||||
def get_tensor_rank2_dataframe(self, view="all", select_symbols=None, decimal=5, verbose=0):
|
||||
def get_tensor_rank2_dataframe(self, view="all", select_symbols=None, decimals=5, verbose=0):
|
||||
"""
|
||||
|
||||
Args:
|
||||
|
|
@ -349,9 +166,11 @@ class SiteSymmetries(Has_Structure):
|
|||
|
||||
sitesym_labels = self.structure.spget_site_symmetries()
|
||||
|
||||
# Symmetric tensor.
|
||||
axx, ayy, azz, axy, axz, ayz = symbols = sp.symbols('axx, ayy, azz, axy, axz, ayz')
|
||||
tensor = sp.Matrix(([axx, axy, axz], [axy, ayy, ayz], [axz, ayz, azz]))
|
||||
# Symmetric tensor in reduced coordinates (direct lattice)
|
||||
# Operations in reduced coords are given by integer matrices and this facilitates the solution
|
||||
# of the system of equations with sympy.
|
||||
Txx, Tyy, Tzz, Txy, Txz, Tyz = symbols = sp.symbols('Txx, Tyy, Tzz, Txy, Txz, Tyz')
|
||||
tensor = sp.Matrix(([Txx, Txy, Txz], [Txy, Tyy, Tyz], [Txz, Tyz, Tzz]))
|
||||
|
||||
indsym = self.structure.indsym
|
||||
rows = []
|
||||
|
|
@ -361,7 +180,14 @@ class SiteSymmetries(Has_Structure):
|
|||
for isym, rotf in enumerate(self.sp_symrel):
|
||||
if indsym[iatom, isym, 3] != iatom: continue
|
||||
m = rotf * tensor * rotf.T - tensor
|
||||
if verbose: sp.pprint(m)
|
||||
if verbose:
|
||||
print(92 * "=")
|
||||
print("System of linear equations for iatom %d, %s" % (iatom, repr(site)))
|
||||
sp.pprint(m)
|
||||
print("Rotation:")
|
||||
sp.pprint(rotf)
|
||||
print(92 * "=")
|
||||
|
||||
system.append(m)
|
||||
|
||||
# Solve system of linear equations.
|
||||
|
|
@ -379,10 +205,10 @@ class SiteSymmetries(Has_Structure):
|
|||
d["wyckoff"] = wlabel
|
||||
d["site_symmetry"] = sitesym_labels[iatom]
|
||||
for s in symbols:
|
||||
d[s] = s
|
||||
d[str(s)] = str(s)
|
||||
if len(solutions) > 1:
|
||||
cprint("Found multiple solutions for iatom %d" % iatom, "red")
|
||||
d.update(**solutions[0])
|
||||
d.update({str(k): str(v) for k, v in solutions[0].items()})
|
||||
rows.append(d)
|
||||
|
||||
import pandas as pd
|
||||
|
|
@ -604,8 +604,7 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
|
||||
return("\n".join(lines))
|
||||
|
||||
def get_conventional_standard_structure(self, international_monoclinic=True,
|
||||
symprec=1e-3, angle_tolerance=5):
|
||||
def get_conventional_standard_structure(self, international_monoclinic=True, symprec=1e-3, angle_tolerance=5):
|
||||
"""
|
||||
Gives a structure with a conventional cell according to certain
|
||||
standards. The standards are defined in :cite:`Setyawan2010`
|
||||
|
|
@ -638,6 +637,18 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
|
||||
return self.__class__.as_structure(standardized_structure)
|
||||
|
||||
def refine(self, symprec=1e-3, angle_tolerance=5):
|
||||
"""
|
||||
Get the refined structure based on detected symmetry. The refined
|
||||
structure is a *conventional* cell setting with atoms moved to the
|
||||
expected symmetry positions.
|
||||
|
||||
Returns: Refined structure.
|
||||
"""
|
||||
sym_finder = SpacegroupAnalyzer(structure=self, symprec=symprec, angle_tolerance=angle_tolerance)
|
||||
new = sym_finder.get_refined_structure()
|
||||
return self.__class__.as_structure(new)
|
||||
|
||||
def abi_sanitize(self, symprec=1e-3, angle_tolerance=5, primitive=True, primitive_standard=False):
|
||||
"""
|
||||
Returns a new structure in which:
|
||||
|
|
@ -658,8 +669,7 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
|
||||
# Refine structure
|
||||
if symprec is not None and angle_tolerance is not None:
|
||||
sym_finder = SpacegroupAnalyzer(structure=structure, symprec=symprec, angle_tolerance=angle_tolerance)
|
||||
structure = sym_finder.get_refined_structure()
|
||||
structure = structure.refine(symprec=symprec, angle_tolerance=angle_tolerance)
|
||||
|
||||
# Convert to primitive structure.
|
||||
if primitive:
|
||||
|
|
@ -932,13 +942,14 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
"Use `overwrite=True` to allow modification."))
|
||||
|
||||
msg = ("Structure object does not have symmetry operations computed from Abinit.\n"
|
||||
"Will call spglib to get symmetry operations.")
|
||||
"Calling spglib to get symmetry operations.")
|
||||
cprint(msg, "magenta")
|
||||
|
||||
spglib_data = SpacegroupAnalyzer(self).get_symmetry_dataset()
|
||||
spgid = spglib_data["number"]
|
||||
symrel, tnons = spglib_data["rotations"], spglib_data["translations"]
|
||||
symafm = [1] * len(symrel) # TODO: Anti-ferromagnetic symmetries are not supported by spglib
|
||||
# TODO: Anti-ferromagnetic symmetries are not supported by spglib
|
||||
symafm = [1] * len(symrel)
|
||||
|
||||
abispg = AbinitSpaceGroup(spgid, symrel, tnons, symafm, has_timerev, inord="C")
|
||||
self.set_abi_spacegroup(abispg)
|
||||
|
|
@ -972,7 +983,7 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
|
||||
@lazy_property
|
||||
def site_symmetries(self):
|
||||
from abipy.core.wyckoff import SiteSymmetries
|
||||
from abipy.core.site_symmetries import SiteSymmetries
|
||||
return SiteSymmetries(self)
|
||||
|
||||
# TODO: site_symmetry or spget_site_symmetries?
|
||||
|
|
@ -984,13 +995,13 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
sitesym_labels = []
|
||||
for iatom, site in enumerate(self):
|
||||
rotations = [symrel[isym] for isym in range(nsym) if
|
||||
indsym[iatom, isym, 3] == iatom and symafm[isym] == 1]
|
||||
indsym[iatom, isym, 3] == iatom and symafm[isym] == +1]
|
||||
# Passing a 0-length rotations list to spglib can segfault.
|
||||
herm_symbol, ptg_num = "1", 1
|
||||
if len(rotations) != 0:
|
||||
herm_symbol, ptg_num, trans_mat = spglib.get_pointgroup(rotations)
|
||||
|
||||
sitesym_labels.append("%s (#%d) nsym:%d" % (herm_symbol.strip(), ptg_num, len(rotations)))
|
||||
sitesym_labels.append("%s (#%d,nsym:%d)" % (herm_symbol.strip(), ptg_num, len(rotations)))
|
||||
|
||||
return sitesym_labels
|
||||
|
||||
|
|
@ -1857,32 +1868,6 @@ class Structure(pymatgen.Structure, NotebookWriter):
|
|||
|
||||
return "\n".join(lines)
|
||||
|
||||
#def ksampling_from_jhudb(self, **kwargs):
|
||||
# from pymatgen.ext.jhu import get_kpoints
|
||||
# __doc__ = get_kpoints.__doc__
|
||||
# kpoints = get_kpoints(self, **kwargs)
|
||||
# print(kpoints)
|
||||
# print(kpoints.style)
|
||||
# print("num_kpts", kpoints.num_kpts)
|
||||
|
||||
# d = {"kptopt": 0,
|
||||
# "kpt": kpoints.kpts,
|
||||
# "nkpt": kpoints.num_kpts,
|
||||
# #"kptnrm": kptnrm,
|
||||
# "wtk": kpoints.kpts_weights,
|
||||
# "shiftk": kpoints.kpts_shift,
|
||||
# "chksymbreak": 0,
|
||||
# }
|
||||
# print(d)
|
||||
|
||||
# #from pymatgen.io.abinit.abiobjects import KSampling
|
||||
# #return KSampling(mode=KSamplingModes.automatic,
|
||||
# # num_kpts= 0,
|
||||
# # kpts=((1, 1, 1),),
|
||||
# # kpt_shifts=(0.5, 0.5, 0.5),
|
||||
# # kpts_weights=None, use_symmetries=True, use_time_reversal=True, chksymbreak=None,
|
||||
# # comment=None)
|
||||
|
||||
def calc_ksampling(self, nksmall, symprec=0.01, angle_tolerance=5):
|
||||
"""
|
||||
Return the k-point sampling from the number of divisions ``nksmall`` to be used for
|
||||
|
|
|
|||
|
|
@ -0,0 +1,56 @@
|
|||
"""Tests for core.site_symmetries module"""
|
||||
from __future__ import print_function, division, unicode_literals, absolute_import
|
||||
|
||||
import unittest
|
||||
import os
|
||||
import numpy as np
|
||||
|
||||
from abipy.core.testing import AbipyTest
|
||||
from abipy.core.structure import Structure
|
||||
from abipy.core.site_symmetries import SiteSymmetries
|
||||
|
||||
import abipy.data as abidata
|
||||
|
||||
class TestSiteSymmetries(AbipyTest):
|
||||
"""Unit tests for SiteSymmetries."""
|
||||
|
||||
def test_si(self):
|
||||
"""Testing wyckoff positions for Si2"""
|
||||
si = Structure.from_file(abidata.cif_file("si.cif"))
|
||||
ss = si.site_symmetries
|
||||
repr(ss); str(ss)
|
||||
assert ss.to_string(verbose=2)
|
||||
df = ss.get_wyckoff_dataframe(verbose=2)
|
||||
self.assert_equal(np.array(df["xfrac"].values, dtype=float), [0, 0.25])
|
||||
self.assert_equal(np.array(df["yfrac"].values, dtype=float), [0, 0.25])
|
||||
self.assert_equal(np.array(df["zfrac"].values, dtype=float), [0, 0.25])
|
||||
#0 -43m (#31) nsym:24 0 0 0
|
||||
#1 -43m (#31) nsym:24 0.250000000000000 0.250000000000000 0.250000000000000
|
||||
|
||||
df = ss.get_tensor_rank2_dataframe(verbose=2)
|
||||
ref = ["Tzz", "Tzz"]
|
||||
self.assert_equal(df["Txx"].values, ref)
|
||||
self.assert_equal(df["Tyy"].values, ref)
|
||||
ref = ["-Tzz/3", "-Tzz/3"]
|
||||
self.assert_equal(df["Txy"].values, ref)
|
||||
self.assert_equal(df["Txz"].values, ref)
|
||||
self.assert_equal(df["Tyz"].values, ref)
|
||||
|
||||
def test_alpha_sio2(self):
|
||||
"""Testing wyckoff positions for alpha-SiO2"""
|
||||
asi02 = Structure.from_file(os.path.join(abidata.dirpath, "refs", "mp-7000_DDB.bz2"))
|
||||
ss = asi02.site_symmetries
|
||||
df = ss.get_wyckoff_dataframe(verbose=2)
|
||||
df[df["element"] == "Si"]
|
||||
self.assert_equal(df["xfrac"].values, ["xfrac", "yfrac", "0.0"])
|
||||
self.assert_equal(df["yfrac"].values, ["0.0", "yfrac", "yfrac"])
|
||||
self.assert_equal(np.array(df["zfrac"].values, dtype=float), [0.833335, 0.5, 0.166665])
|
||||
|
||||
#df = ss.get_tensor_rank2_dataframe()
|
||||
#ref = ["Tzz", "Tzz"]
|
||||
#self.assert_equal(df["Txx"].values, ref)
|
||||
#self.assert_equal(df["Tyy"].values, ref)
|
||||
#ref = ["-Tzz/3", "-Tzz/3"]
|
||||
#self.assert_equal(df["Txy"].values, ref)
|
||||
#self.assert_equal(df["Txz"].values, ref)
|
||||
#self.assert_equal(df["Tyz"].values, ref)
|
||||
|
|
@ -1,61 +0,0 @@
|
|||
"""Tests for core.density module"""
|
||||
from __future__ import print_function, division, unicode_literals, absolute_import
|
||||
|
||||
import unittest
|
||||
|
||||
from collections import OrderedDict
|
||||
from abipy.core.testing import AbipyTest
|
||||
from abipy.core.structure import Lattice
|
||||
from abipy.core.wyckoff import Wyckoff
|
||||
|
||||
|
||||
class TestWyckoff(AbipyTest):
|
||||
"""Unit tests for Wyckoff."""
|
||||
|
||||
def test_sio2(self):
|
||||
"""Testing Wyckoff positions with SiO2 (requires sympy)"""
|
||||
try:
|
||||
import sympy
|
||||
except ImportError:
|
||||
raise unittest.SkipTest("sympy is not installed")
|
||||
|
||||
# http://www.cryst.ehu.es/cgi-bin/cryst/programs/nph-wp-list?gnum=152
|
||||
site2wpos = {}
|
||||
#3 b .2.
|
||||
#site2wpos["Si"] = "(x,0,5/6) (0,x,1/6) (-x,-x,1/2)"
|
||||
#3 a .2.
|
||||
site2wpos["Si"] = "[x,0,1/3], (0,x,2/3), (-x,-x,0)"
|
||||
#6 c 1
|
||||
site2wpos["O"] = """(x,y,z), (-y,x-y,z+1/3), (-x+y,-x,z+2/3), (y,x,-z),
|
||||
(x-y,-y,-z+2/3), (-x,-x+y,-z+1/3)"""
|
||||
|
||||
wyckoff = Wyckoff(site2wpos)
|
||||
repr(wyckoff); str(wyckoff)
|
||||
|
||||
site2params = OrderedDict([
|
||||
("Si", dict(x=0.4763)),
|
||||
("O", dict(x=0.1588, y=0.7439, z=0.4612)),
|
||||
#("O", dict(x=0.1588, y=0.7439, z=0.4612, yyy=3)),
|
||||
])
|
||||
|
||||
lattice = Lattice.hexagonal(a=4.971, c=5.473)
|
||||
for to_unit_cell in [False]:
|
||||
#for to_unit_cell in [False, True]:
|
||||
structure = wyckoff.generate_structure(lattice, site2params, to_unit_cell=to_unit_cell)
|
||||
print("to_unit_cell %s\n" % to_unit_cell, structure)
|
||||
|
||||
#from pymatgen.symmetry.analyzer import SpacegroupAnalyzer
|
||||
#spga = SpacegroupAnalyzer(structure)
|
||||
#structure = spga.get_refined_structure()
|
||||
#print("Refined:", structure)
|
||||
|
||||
wyck_params = wyckoff.find_params(structure)
|
||||
#structure.perturb(distance=0.01)
|
||||
#wyck_params = wyckoff.find_params(structure)
|
||||
|
||||
print("Error of params")
|
||||
wyckoff.error_of_params(wyck_params, structure)
|
||||
|
||||
for site, params in site2params.items():
|
||||
print(wyck_params[site])
|
||||
assert wyck_params[site] == params
|
||||
|
|
@ -46,11 +46,11 @@ class _Component(object):
|
|||
|
||||
class MsqDos(Has_Structure):
|
||||
"""
|
||||
This object stores the generalized phonon DOS with the mean square displacement tensor in cartesian coords.
|
||||
This object stores the generalized phonon DOS with the mean square displacement tensor in CARTESIAN coords.
|
||||
This DOS-like quantity allows one to calculate Debye Waller factors as a function of T
|
||||
by integration with 1/omega and the Bose-Einstein factor.
|
||||
See :cite:`Lee1995` for the further details about the internal implementation and
|
||||
:cite:`Trueblood1996` for the different conventions used by crystallographers
|
||||
:cite:`Trueblood1996` for the different conventions used by crystallographers.
|
||||
"""
|
||||
C = _Component
|
||||
ALL_COMPS = OrderedDict([
|
||||
|
|
@ -73,7 +73,7 @@ class MsqDos(Has_Structure):
|
|||
Arg:
|
||||
structure: |Structure| object.
|
||||
wmesh: Frequency mesh
|
||||
values: (natom, 3, 3, nomega) arrays with generalized DOS.
|
||||
values: (natom, 3, 3, nomega) arrays with generalized DOS in cart coords.
|
||||
amu_symbol: Dictionary element.symbol -> mass in atomic units.
|
||||
"""
|
||||
self._structure = structure
|
||||
|
|
@ -120,7 +120,7 @@ class MsqDos(Has_Structure):
|
|||
Compute mean square displacement for each atom in `iatom_list` as a function of T (bohr^2).
|
||||
|
||||
Args:
|
||||
tmesh: array-like with temperatures in Kelvin degrees
|
||||
tmesh: array-like with temperatures in Kelvin.
|
||||
iatom_list: List of atom sites to comput. None for all.
|
||||
what_list:
|
||||
"""
|
||||
|
|
@ -252,7 +252,7 @@ class MsqDos(Has_Structure):
|
|||
|
||||
def write_cif_file(self, filepath, temp=300):
|
||||
"""
|
||||
Write CIF file with structure info and anisotropic U tensor in CIF format.
|
||||
Write CIF file with structure and anisotropic U tensor in CIF format.
|
||||
|
||||
Args:
|
||||
filepath: Name of CIF file. If None, a temporary filepath is created.
|
||||
|
|
@ -275,7 +275,7 @@ class MsqDos(Has_Structure):
|
|||
In the Vesta GUI, select: Properties -> Atoms -> Show as displament ellipsoids.
|
||||
"""
|
||||
filepath = self.write_cif_file(filepath=None, temp=temp)
|
||||
print("Writing structure + Debye-Waller tensor in CIF format for T = %s to file: %s" % (temp, filepath))
|
||||
cprint("Writing structure + Debye-Waller tensor in CIF format for T = %s to file: %s" % (temp, filepath), "green")
|
||||
from abipy.iotools import Visualizer
|
||||
visu = Visualizer.from_name("vesta")
|
||||
|
||||
|
|
@ -318,15 +318,16 @@ _atom_site_aniso_U_12""".splitlines()
|
|||
return s + "\n".join(aniso_u)
|
||||
|
||||
def check_site_symmetries(self, temp=300, verbose=0):
|
||||
|
||||
"""
|
||||
Check symmetries of the displacement tensor at temperature `temp`.
|
||||
Return maximum error.
|
||||
"""
|
||||
msq = self.get_msq_tmesh([float(temp)], what_list="displ")
|
||||
values = getattr(msq, "displ")
|
||||
values_cart = getattr(msq, "displ")
|
||||
natom = len(self.structure)
|
||||
values = np.reshape(values, (natom, 3, 3))
|
||||
values_cart = np.reshape(values_cart, (natom, 3, 3))
|
||||
|
||||
from abipy.core.wyckoff import SiteSymmetries
|
||||
ss = SiteSymmetries(self.structure)
|
||||
return ss.check_site_symmetries(values, verbose=verbose)
|
||||
return self.structure.site_symmetries.check_site_symmetries(values_cart, verbose=verbose)
|
||||
|
||||
def _get_components(self, components):
|
||||
"""
|
||||
|
|
|
|||
|
|
@ -2833,7 +2833,7 @@ class PhdosFile(AbinitNcFile, Has_Structure, NotebookWriter):
|
|||
msqd_dos = None
|
||||
if msqd_dos is not None:
|
||||
yield msqd_dos.plot(units=units, show=False)
|
||||
yield msqd_dos.plot_bfactors(show=False)
|
||||
yield msqd_dos.plot_tensor(show=False)
|
||||
|
||||
def write_notebook(self, nbpath=None):
|
||||
"""
|
||||
|
|
|
|||
|
|
@ -13,7 +13,6 @@ from abipy.dfpt.phonons import PhononBands, get_dyn_mat_eigenvec, match_eigenvec
|
|||
from abipy.abio.inputs import AnaddbInput
|
||||
from abipy.tools.plotting import add_fig_kwargs, get_ax_fig_plt
|
||||
from pymatgen.core.units import bohr_to_angstrom, eV_to_Ha
|
||||
from pymatgen.symmetry.bandstructure import HighSymmKpath
|
||||
|
||||
|
||||
class SoundVelocity(Has_Structure, NotebookWriter):
|
||||
|
|
@ -108,8 +107,7 @@ class SoundVelocity(Has_Structure, NotebookWriter):
|
|||
)
|
||||
|
||||
if not directions:
|
||||
hs = HighSymmKpath(ddb.structure, symprec=1e-2)
|
||||
|
||||
hs = ddb.structure.hsym_kpath
|
||||
kpath = hs.kpath
|
||||
|
||||
directions = []
|
||||
|
|
|
|||
|
|
@ -227,9 +227,11 @@ codes), a looser tolerance of 0.1 (the value used in Materials Project) is often
|
|||
help="Extract/Compute Abinit space group from file with structure.")
|
||||
p_abispg.add_argument("-t", "--tolsym", type=float, default=None, help="""\
|
||||
Gives the tolerance on the atomic positions (reduced coordinates), primitive vectors, or magnetization,
|
||||
to be considered equivalent, thanks to symmetry operations. This is used in the recognition of the set
|
||||
to be considered equivalent, thanks to symmetry operations. This value is used by ABINIT in the recognition of the set
|
||||
of symmetries of the system, or the application of the symmetry operations to generate from a reduced set of atoms,
|
||||
the full set of atoms. Note that a value larger than 0.01 is considered to be unacceptable.""")
|
||||
The internal default is 1e-8. Setting tolsym to a value larger than 1e-8 will make Abinit detect the spacegroup within
|
||||
this tolerance and re-symmetrize the input structure. This option is useful if the structure has been taken from a CIF
|
||||
file that does not have enough significant digits.""")
|
||||
p_abispg.add_argument("-d", "--diff-mode", type=str, default="table", choices=["table", "diff"],
|
||||
help="Select diff output format.")
|
||||
|
||||
|
|
@ -283,12 +285,16 @@ Has to be all integers. Several options are possible:
|
|||
p_proto.add_argument("--angle-tol", default=5, type=float, help="angle tolerance.")
|
||||
|
||||
# Subparser for wyckoff.
|
||||
p_wyckoff = subparsers.add_parser('wyckoff', parents=[copts_parser, path_selector],
|
||||
p_wyckoff = subparsers.add_parser('wyckoff', parents=[copts_parser, spgopt_parser, path_selector],
|
||||
help="Print wyckoff positions. WARNING: still under development!")
|
||||
p_wyckoff.add_argument("--refine", default=False, action="store_true",
|
||||
help="Use spglib to refine structure before computation")
|
||||
|
||||
# Subparser for tensor_site.
|
||||
p_tensor_site = subparsers.add_parser('tensor_site', parents=[copts_parser, path_selector],
|
||||
p_tensor_site = subparsers.add_parser('tensor_site', parents=[copts_parser, spgopt_parser, path_selector],
|
||||
help="Print symmetry properties of tensors due to site-symmetries. WARNING: still under development!")
|
||||
p_tensor_site.add_argument("--refine", default=False, action="store_true",
|
||||
help="Use spglib to refine structure before computation")
|
||||
|
||||
# Subparser for neighbors.
|
||||
p_neighbors = subparsers.add_parser('neighbors', parents=[copts_parser, path_selector],
|
||||
|
|
@ -498,20 +504,25 @@ def main():
|
|||
elif options.command == "abispg":
|
||||
structure = abilab.Structure.from_file(options.filepath)
|
||||
check_ordered_structure(structure)
|
||||
spgrp = structure.abi_spacegroup
|
||||
abi_spg = structure.abi_spacegroup
|
||||
|
||||
if spgrp is not None:
|
||||
if abi_spg is not None and options.tolsym is None:
|
||||
print(structure.spget_summary(verbose=options.verbose))
|
||||
else:
|
||||
# Here we compare Abinit wrt spglib. If spgrp is None, we create a temporary
|
||||
# Here we compare Abinit wrt spglib. If abi_spg is None, we create a temporary
|
||||
# task to run the code in dry-run mode.
|
||||
print("FILE does not contain Abinit symmetry operations.")
|
||||
print("Calling Abinit in --dry-run mode with chkprim = 0 to get space group.")
|
||||
if abi_spg is None:
|
||||
print("FILE does not contain Abinit symmetry operations.")
|
||||
cprint("Calling Abinit in --dry-run mode with chkprim = 0 to get space group.")
|
||||
if options.tolsym is not None and options.tolsym > 1e-8:
|
||||
cprint("Crystal structure will be re-symmetrized by Abinit with tolsym: %s" % options.tolsym, "yellow")
|
||||
|
||||
from abipy.data.hgh_pseudos import HGH_TABLE
|
||||
gsinp = factories.gs_input(structure, HGH_TABLE, spin_mode="unpolarized")
|
||||
gsinp["chkprim"] = 0
|
||||
abistructure = gsinp.abiget_spacegroup(tolsym=options.tolsym)
|
||||
print(abistructure.spget_summary(verbose=options.verbose))
|
||||
print("")
|
||||
|
||||
diff_structures([structure, abistructure], mode=options.diff_mode,
|
||||
headers=["Input structure", "After Abinit symmetrization"], fmt="abivars")
|
||||
|
|
@ -659,17 +670,23 @@ def main():
|
|||
|
||||
elif options.command == "wyckoff":
|
||||
structure = abilab.Structure.from_file(options.filepath)
|
||||
from abipy.core.wyckoff import SiteSymmetries
|
||||
ss = SiteSymmetries(structure)
|
||||
df = ss.get_wyckoff_dataframe(view="all", select_symbols=None, verbose=options.verbose)
|
||||
abilab.print_dataframe(df)
|
||||
if options.refine:
|
||||
print("Refining structure with symprec: %s, angle_tolerance: %s" % (options.symprec, options.angle_tolerance))
|
||||
structure = structure.refine(symprec=options.symprec, angle_tolerance=options.angle_tolerance)
|
||||
print(structure.spget_summary(verbose=options.verbose))
|
||||
ss = structure.site_symmetries
|
||||
df = ss.get_wyckoff_dataframe(verbose=options.verbose)
|
||||
abilab.print_dataframe(df, title="\nWyckoff positions in reduced coordinates.")
|
||||
|
||||
elif options.command == "tensor_site":
|
||||
structure = abilab.Structure.from_file(options.filepath)
|
||||
from abipy.core.wyckoff import SiteSymmetries
|
||||
ss = SiteSymmetries(structure)
|
||||
df = ss.get_tensor_rank2_dataframe(view="all", select_symbols=None, verbose=options.verbose)
|
||||
abilab.print_dataframe(df)
|
||||
if options.refine:
|
||||
print("Refining structure with symprec: %s, angle_tolerance: %s" % (options.symprec, options.angle_tolerance))
|
||||
structure = structure.refine(symprec=options.symprec, angle_tolerance=options.angle_tolerance)
|
||||
print(structure.spget_summary(verbose=options.verbose))
|
||||
ss = structure.site_symmetries
|
||||
df = ss.get_tensor_rank2_dataframe(verbose=options.verbose)
|
||||
abilab.print_dataframe(df, title="\nTensor components in reduced coordinates (rank 2, symmetric)")
|
||||
|
||||
elif options.command == "neighbors":
|
||||
abilab.Structure.from_file(options.filepath).print_neighbors(radius=options.radius)
|
||||
|
|
@ -757,8 +774,11 @@ def main():
|
|||
|
||||
#elif options.command == "kmesh_jhu":
|
||||
# structure = abilab.Structure.from_file(options.filepath)
|
||||
# ksampling = structure.ksampling_from_jhudb(kppra=1000)
|
||||
# #print(ksampling)
|
||||
# from pymatgen.ext.jhu import get_kpoints
|
||||
# kpoints = get_kpoints(structure, min_distance=0, min_total_kpoints=1,
|
||||
# kppra=None, gap_distance=7, remove_symmetry=None,
|
||||
# include_gamma="auto", header="simple", incar=None)
|
||||
# #print(kpoints)
|
||||
|
||||
elif options.command == "lgk":
|
||||
structure = abilab.Structure.from_file(options.filepath)
|
||||
|
|
|
|||
|
|
@ -124,10 +124,10 @@ core Package
|
|||
:undoc-members:
|
||||
:show-inheritance:
|
||||
|
||||
:mod:`wyckoff` Module
|
||||
---------------------
|
||||
:mod:`site_symmetries` Module
|
||||
-----------------------------
|
||||
|
||||
.. automodule:: abipy.core.wyckoff
|
||||
.. automodule:: abipy.core.site_symmetries
|
||||
:members:
|
||||
:undoc-members:
|
||||
:show-inheritance:
|
||||
|
|
|
|||
|
|
@ -72,9 +72,9 @@ By default, the installer creates the ``anaconda`` directory in your home.
|
|||
Anaconda will add one line to your ``.bashrc`` to enable access to the anaconda executables.
|
||||
Once the installation is completed, execute::
|
||||
|
||||
source ~/anaconda/bin/activate root
|
||||
source ~/anaconda/bin/activate base
|
||||
|
||||
to activate the ``root`` environment.
|
||||
to activate the ``base`` environment.
|
||||
The output of ``which python`` should show that you are using the python interpreter provided by anaconda.
|
||||
|
||||
Use the conda_ command-line interface to install the packages not included in the official distribution.
|
||||
|
|
@ -86,17 +86,17 @@ Remember that if a package is not available in the official conda repository, yo
|
|||
download the package from one of the conda channels or use ``pip install`` if no conda package is available.
|
||||
|
||||
Fortunately there are conda channels providing all dependencies needed by AbiPy.
|
||||
To install the pymatgen_ package from the matsci_ channel, use::
|
||||
Now add ``conda-forge``, ``matsci`` and ``abinit`` to your conda channels with::
|
||||
|
||||
conda install pymatgen --channel matsci
|
||||
conda config --add channels conda-forge
|
||||
conda config --add channels matsci
|
||||
conda config --add channels abinit
|
||||
|
||||
then install Abipy from the abinit-channel_ with::
|
||||
These are the channels from which we will download pymatgen, abipy and abinit.
|
||||
Finally, install AbiPy from the abinit-channel_ with::
|
||||
|
||||
conda install abipy --channel abinit
|
||||
|
||||
Visit `materials.sh <http://materials.sh>`_ for instructions on how to use the
|
||||
matsci channel to install pymatgen and other packages.
|
||||
|
||||
Once you have completed the installation of AbiPy and pymatgen, open the ipython_ shell and type::
|
||||
|
||||
# make sure spglib library works
|
||||
|
|
|
|||
|
|
@ -64,7 +64,7 @@ Use::
|
|||
to read the structure from ``FILE`` and generate a CIF_ file (default behaviour).
|
||||
|
||||
The majority of the netcdf_ files produced by Abinit contain structural information
|
||||
so this command can be used with netcdf ouput files as well as Abinit input/output
|
||||
so this command can be used with netcdf output files as well as Abinit input/output
|
||||
files and all the other formats supported by pymatgen e.g. POSCAR files.
|
||||
Other formats can be specified with the ``-f`` option.
|
||||
For example::
|
||||
|
|
@ -101,7 +101,7 @@ to validate the input file with Abinit (requires ``manager.yml`` and, obviously,
|
|||
The script provides other options to invoke Abinit
|
||||
to get space group information, the list of k-points in the IBZ.
|
||||
the list of atomic perturbations for phonons or the list of autoparal configurations.
|
||||
See ``abinp.py --help`` for futher info.
|
||||
See ``abinp.py --help`` for further info.
|
||||
|
||||
Print the warnings in the log file
|
||||
----------------------------------
|
||||
|
|
@ -232,6 +232,24 @@ to generate a template with the input variables defining the k-path
|
|||
+0.62500 +0.25000 +0.62500 # U
|
||||
+0.50000 +0.00000 +0.50000 # X
|
||||
|
||||
|
||||
Re-symmetrize a structure when Abinit reports less symmetries than expected
|
||||
---------------------------------------------------------------------------
|
||||
|
||||
Crystalline structures saved in text format (e.g. CIF files downloaded from
|
||||
the Materials Project websites) may not have enough significant digits
|
||||
and Abinit may not find the same spacegroup as the one reported by the source
|
||||
as the default tolerance for symmetry detection in Abinit is tight (tolsym = 1e-8).
|
||||
|
||||
In this case, one can use the `abispg` option of abistruct.py to compute the spacegroup
|
||||
with Abinit and a tolerance larger that the default value::
|
||||
|
||||
abistruct.py abispg problematic.cif --tolsym=1e-6
|
||||
|
||||
Hopefully, the code will detect the correct spacegroup, will re-symmetrize
|
||||
the initial lattice vectors and atomic positions and print the new symmetrized structure to terminal.
|
||||
|
||||
|
||||
Get neighbors for each atom in the unit cell out to a distance radius
|
||||
---------------------------------------------------------------------
|
||||
|
||||
|
|
@ -402,7 +420,7 @@ Use::
|
|||
|
||||
abiview.py phbands out_PHBST.nc -web
|
||||
|
||||
to start a local webserver and open the HTML page inside the default browser
|
||||
to start a local web server and open the HTML page inside the default browser
|
||||
(the browser can be changed with the ``--browser`` option).
|
||||
|
||||
It is also possible to visualize the phonon modes starting directly from a DDB_ file with::
|
||||
|
|
@ -410,7 +428,7 @@ It is also possible to visualize the phonon modes starting directly from a DDB_
|
|||
abiview.py ddb -web
|
||||
|
||||
In this case, AbiPy will invoke anaddb to produce the ``PHBST.nc`` file on an automatically
|
||||
generated q-path and then start the webserver.
|
||||
generated q-path and then start the web server.
|
||||
|
||||
Visualize the results of a structural relaxation
|
||||
------------------------------------------------
|
||||
|
|
@ -461,7 +479,7 @@ the AbiPy scripts are quite handy for a quick analysis of the results.
|
|||
Compare multiple files
|
||||
----------------------
|
||||
|
||||
The :ref:`abicomp.py` script is explicitely designed for this kind of task.
|
||||
The :ref:`abicomp.py` script is explicitly designed for this kind of task.
|
||||
It operates on multiple files (usually files with the same extension) and
|
||||
either produces matplotlib_ plots or creates AbiPy robots providing methods
|
||||
to analyze the results, perform convergence studies and build pandas DataFrames_.
|
||||
|
|
|
|||
Loading…
Reference in New Issue