Skip to content

make compatible with Python 3 #13

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
jotelha wants to merge 1 commit into
base: master
Choose a base branch
from
Open

make compatible with Python 3 #13

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
152 changes: 84 additions & 68 deletions castep.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,10 @@
"""

import re
import scipy as S
import numpy as np
import logging

logger = logging.getLogger(__name__)

version = 0.1

Expand All @@ -32,35 +35,42 @@
dotcastep_poinggroup_RE = re.compile(r"^\s+Point group of crystal =\s+([\+\-]?\d+):")

def parse_dotcastep(seedname):
"""
Extract lattice and atom positions from a .castep
file. List of atoms may be empty (e.g. MgO)
"""
dotCastep = open(seedname+".castep","r")
# Find the lattice
latticeblock = dotcastep_latt_RE.findall(dotCastep.read())[-1] # Get the last block - handle concat restarts
lattice = []
lattice.append([float(latticeblock[0]), float(latticeblock[1]), float(latticeblock[2])])
lattice.append([float(latticeblock[6]), float(latticeblock[7]), float(latticeblock[8])])
lattice.append([float(latticeblock[12]), float(latticeblock[13]), float(latticeblock[14])])
# rewind and search for and final atomic positions (these will be absent if e.g. they are all on symmetry positions)
dotCastep.seek(0)
in_atoms = False
pointgroup = None
atoms = []
for line in dotCastep:
sym_line = dotcastep_poinggroup_RE.search(line)
atom_line = dotcastep_atomline_RE.search(line)
if (in_atoms and atom_line):
atoms.append([atom_line.group(1), float(atom_line.group(2)), \
float(atom_line.group(3)), float(atom_line.group(4))])
elif ((not in_atoms) and (dotcastep_infinal_RE.search(line))):
in_atoms = True
elif (sym_line):
pointgroup = int(sym_line.group(1))

dotCastep.close()
return (lattice, pointgroup, atoms)
"""
Extract lattice and atom positions from a .castep
file. List of atoms may be empty (e.g. MgO)
"""
with open(seedname+".castep","r") as dotCastep:
# Find the lattice
logger.debug("Find lattice block.")
latticeblock = dotcastep_latt_RE.findall(dotCastep.read())[-1] # Get the last block - handle concat restarts
logger.debug("Found lattice block '%s'.", latticeblock)
lattice = []
lattice.append([float(latticeblock[0]), float(latticeblock[1]), float(latticeblock[2])])
lattice.append([float(latticeblock[6]), float(latticeblock[7]), float(latticeblock[8])])
lattice.append([float(latticeblock[12]), float(latticeblock[13]), float(latticeblock[14])])
# rewind and search for and final atomic positions (these will be absent if e.g. they are all on symmetry positions)
dotCastep.seek(0)
in_atoms = False
pointgroup = None
atoms = []
logger.debug("Find atoms block.")
for line in dotCastep:
sym_line = dotcastep_poinggroup_RE.search(line)
atom_line = dotcastep_atomline_RE.search(line)
if (in_atoms and atom_line):
logger.debug("Found atom line '%s'.", line)
logger.debug("Append atom match '%s, %s ,%s, %s'.", atom_line.group(1), atom_line.group(2), atom_line.group(3), atom_line.group(4))
atoms.append([atom_line.group(1), float(atom_line.group(2)), \
float(atom_line.group(3)), float(atom_line.group(4))])
elif ((not in_atoms) and (dotcastep_infinal_RE.search(line))):
logger.debug("Found beginning of atoms block at '%s'.", line)
in_atoms = True
elif (sym_line):
logger.debug("Found symmetry line '%s'.", line.strip())
logger.debug("Append symmetry match '%s'.", sym_line.group(1))
pointgroup = int(sym_line.group(1))

return (lattice, pointgroup, atoms)


# Regular expressions to match a lattice block in a castep .cell file. Note that these
Expand All @@ -71,49 +81,55 @@ def parse_dotcastep(seedname):
dotcell_atoms_end_RE = re.compile(r"^\s*%ENDBLOCK\s+POSITIONS_(?:FRAC|ABS)", re.IGNORECASE)

def produce_dotcell(seedname, filename, defcell, atoms):
"""
produce_dotcell: reads <seedname>.cell (CASTEP cell file
and writes a new .cell file to <filename> replacing the
lattice block with a new crystalographic lattice <defcell>
(which should be supplied as a list of three lists, each with
three elements). Also adds command to fix cell during optimization.
"""
in_lattice = False
in_atoms = False
have_atoms = (atoms != []) # If we have an empty list, no atoms were optimized so just leave them in the .cell file.
inputfile = open(seedname+".cell", "r")
outputfile = open(filename, "w")
for line in inputfile:
if (dotcell_lattice_end_RE.search(line) and in_lattice):
in_lattice = False
elif (dotcell_lattice_start_RE.search(line) and not in_lattice):
outputfile.write("%block LATTICE_CART\n")
outputfile.write(str(defcell[0][0]) + " " + str(defcell[0][1]) + " " + str(defcell[0][2]) + "\n")
outputfile.write(str(defcell[1][0]) + " " + str(defcell[1][1]) + " " + str(defcell[1][2]) + "\n")
outputfile.write(str(defcell[2][0]) + " " + str(defcell[2][1]) + " " + str(defcell[2][2]) + "\n")
outputfile.write("%endblock LATTICE_CART\n")
outputfile.write("FIX_ALL_CELL true\n")
in_lattice = True
elif (dotcell_atoms_end_RE.search(line) and in_atoms and have_atoms):
in_atoms = False
elif ((dotcell_atoms_start_RE.search(line)) and (not in_atoms) and have_atoms):
outputfile.write("%block POSITIONS_FRAC\n")
for atom in atoms:
outputfile.write(" " + atom[0] + " " + str(atom[1]) + " " + str(atom[2]) + " " + str(atom[3]) + "\n")
outputfile.write("%endblock POSITIONS_FRAC\n")
in_atoms = True
elif(not (in_lattice or in_atoms)):
outputfile.write(line)
inputfile.close
outputfile.close
return()
"""
produce_dotcell: reads <seedname>.cell (CASTEP cell file
and writes a new .cell file to <filename> replacing the
lattice block with a new crystalographic lattice <defcell>
(which should be supplied as a list of three lists, each with
three elements). Also adds command to fix cell during optimization.
"""
in_lattice = False
in_atoms = False
have_atoms = (atoms != []) # If we have an empty list, no atoms were optimized so just leave them in the .cell file.
with open(seedname+".cell", "r") as inputfile, open(filename, "w") as outputfile:
for line in inputfile:
logger.debug("Process line '%s'", line.strip())
if (dotcell_lattice_end_RE.search(line) and in_lattice):
logger.debug("Reached end of LATTICE_CART block.")
in_lattice = False
elif (dotcell_lattice_start_RE.search(line) and not in_lattice):
logger.debug("Write LATTICE_CART block.")
outputfile.write("%block LATTICE_CART\n")
outputfile.write(str(defcell[0][0]) + " " + str(defcell[0][1]) + " " + str(defcell[0][2]) + "\n")
outputfile.write(str(defcell[1][0]) + " " + str(defcell[1][1]) + " " + str(defcell[1][2]) + "\n")
outputfile.write(str(defcell[2][0]) + " " + str(defcell[2][1]) + " " + str(defcell[2][2]) + "\n")
outputfile.write("%endblock LATTICE_CART\n")
outputfile.write("FIX_ALL_CELL true\n")
in_lattice = True
elif (dotcell_atoms_end_RE.search(line) and in_atoms and have_atoms):
logger.debug("Reached end of POSITIONS_FRAC block.")
in_atoms = False
elif ((dotcell_atoms_start_RE.search(line)) and (not in_atoms) and have_atoms):
logger.debug("Write POSITIONS_FRAC block.")
outputfile.write("%block POSITIONS_FRAC\n")
for atom in atoms:
logger.debug("Write atom '%s'.", atom)
outputfile.write(" " + atom[0] + " " + str(atom[1]) + " " + str(atom[2]) + " " + str(atom[3]) + "\n")
outputfile.write("%endblock POSITIONS_FRAC\n")
in_atoms = True
elif(not (in_lattice or in_atoms)):
logger.debug("Write line back to output file unchanged.")
outputfile.write(line)
else:
logger.debug("Do not write anything to output file.")
return()

# regular expression which matches the whole stress tensor block from a .castep file
stressRE = re.compile(r"\s\*+\s(?:Symmetrised\s)?Stress\sTensor\s\*+\n.+\n.+?\((\w+)\).+\n.+\n.+\n.+\n\s\*\s+x\s+([\+\-]?\d+.\d+)\s+([\+\-]?\d+.\d+)\s+([\+\-]?\d+.\d+)\s+\*\n\s\*\s+y\s+[\+\-]?\d+.\d+\s+([\+\-]?\d+.\d+)\s+([\+\-]?\d+.\d+)\s+\*\n\s\*\s+z\s+[\+\-]?\d+.\d+\s+[\+\-]?\d+.\d+\s+([\+\-]?\d+.\d+)\s+\*\n")

def get_stress_dotcastep(filename):
"""Extract the stress tensor from a .castep file

set_facecolor
Returns a tuple of (<units>, <stress>) where <units>
is a string representing the stress units and
<stress> is a numpy vector of the elements of the
Expand All @@ -124,7 +140,7 @@ def get_stress_dotcastep(filename):
stressData = stressRE.findall(dotCastep.read())[0]
dotCastep.close()
units = stressData[0]
stress = S.array([float(stressData[1]),float(stressData[4]),
stress = np.array([float(stressData[1]),float(stressData[4]),
float(stressData[6]),float(stressData[5]),
float(stressData[3]),float(stressData[2])])
return(units, stress)
Expand Down
Loading