Skip to content

Elastic properties #16

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.

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
wants to merge 5 commits into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
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
22 changes: 22 additions & 0 deletions scripts/vasprun
Original file line number Diff line number Diff line change
Expand Up @@ -43,6 +43,22 @@ if __name__ == "__main__":
help="figure' dpi, default: 300", metavar="dpi")
parser.add_option("-S", "--save-bands", dest="saveBands", action='store_true',
help="saving of bands contributing to the plot, default: false", metavar="saveBands")
parser.add_option("-E", "--elastic", dest="runElasticCalculations", action='store_true',
help="run calculations for the elastic properties (requires IBRION >=5), default: false")
parser.add_option("-C", "--elasticsymmetry", dest="symmetry", default="Triclinic",
choices=['Triclinic', 'Monoclinic', 'Orthorhombic',
'Trigonal6', 'Trigonal7', 'Hexagonal',
'Tetragonal6', 'Tetragonal7', 'Cubic'],
help="symmetry group of the supercell, default: Triclinic;\
For calculation of elastic properties only!\
choices: \"Triclinic\", \"Monoclinic\", \"Orthorhombic\",\
\"Trigonal6\", \"Trigonal7\", \"Hexagonal\",\
\"Tetragonal6\", \"Tetragonal7\", and \"Cubic\"")
parser.add_option("-A", "--approximation", dest="approximation", default="voigt",
choices=['voigt','reuss','hill','all'],
help="approximation used for the calculation of the elastic moduli, default: voigt;\
For calculation of elastic properties only!\
choices: \"voigt\", \"reuss\", \"hill\", and \"all\"")

(options, args) = parser.parse_args()
if options.vasprun is None:
Expand Down Expand Up @@ -145,3 +161,9 @@ if __name__ == "__main__":
print("{:25s} {:12.3f} {:12.3f} {:12.3f}".format('DFPT', eps[0,0], eps[0,1], eps[0,2]))
print("{:25s} {:12.3f} {:12.3f} {:12.3f}".format('DFPT', eps[1,0], eps[1,1], eps[1,2]))
print("{:25s} {:12.3f} {:12.3f} {:12.3f}".format('DFPT', eps[2,0], eps[2,1], eps[2,2]))
elif options.runElasticCalculations:
from vasprun.elastic import elastic
solver = elastic.Calculator(options=options,verbosity=0)
solver() # running calculations
solver.print()

2 changes: 1 addition & 1 deletion setup.py
Original file line number Diff line number Diff line change
Expand Up @@ -17,7 +17,7 @@
long_description = long_description,
long_description_content_type = 'text/markdown',
url="https://github.com/qzhu2017/vasprun",
packages=['vasprun'],
packages=['vasprun','vasprun.elastic'],
scripts=['scripts/vasprun'],
classifiers=[
"Programming Language :: Python :: 3",
Expand Down
299 changes: 299 additions & 0 deletions vasprun/elastic/CMatrices.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,299 @@
"""
Functions allowing for calculation a 6x6 elastic-constant matrix
for 9 different point-symmetry space groups. Make sure that you use
the proper ordering of the Cij constants!
"""
import numpy as np

class ElasticTensor:
DegreesOfFreedom = {
'Triclinic' : 21,
'Monoclinic' : 13,
'Orthorhombic' : 9,
'Trigonal6' : 6,
'Trigonal7' : 7,
'Hexagonal' : 5,
'Tetragonal6' : 6,
'Tetragonal7' : 7,
'Cubic' : 3,
}
def __init__(self,symmetry=None):
self.CMatrix = {
'Triclinic': ElasticTensor.CMatrixTriclinic,
'Monoclinic': ElasticTensor.CMatrixMonoclinic,
'Orthorhombic': ElasticTensor.CMatrixOrthorhombic,
'Trigonal6': ElasticTensor.CMatrixTrigonal6,
'Trigonal7': ElasticTensor.CMatrixTrigonal7,
'Hexagonal': ElasticTensor.CMatrixHexagonal,
'Tetragonal6': ElasticTensor.CMatrixTetragonal6,
'Tetragonal7': ElasticTensor.CMatrixTetragonal7,
'Cubic': ElasticTensor.CMatrixCubic
}

self.symmetry = None
if isinstance(symmetry,str):
self.symmetry = symmetry

def __call__(self,*args,symmetry=None):
if symmetry is not None:
return self.CMatrix[symmetry](*args)
elif self.symmetry is not None:
return self.CMatrix[self.symmetry](*args)
else:
raise TypeError("You must define symmetry!")

@staticmethod
def CMatrixTriclinic(C11,C12,C13,C14,C15,C16,
C22,C23,C24,C25,C26,
C33,C34,C35,C36,
C44,C45,C46,
C55,C56,
C66):
"""
Triclinic cell - 21 constants
Order:
C11,C12,C13,C14,C15,C16,
C22,C23,C24,C25,C26,
C33,C34,C35,C36,
C44,C45,C46,
C55,C56,
C66
returns an array: [
[ C11, C12, C13, C14, C15, C16 ],
[ C12, C22, C23, C24, C25, C26 ],
[ C13, C23, C33, C34, C35, C36 ],
[ C14, C24, C34, C44, C45, C46 ],
[ C15, C25, C35, C45, C55, C56 ],
[ C16, C26, C36, C46, C56, C66 ],
]
"""
return np.array([
[ C11, C12, C13, C14, C15, C16 ],
[ C12, C22, C23, C24, C25, C26 ],
[ C13, C23, C33, C34, C35, C36 ],
[ C14, C24, C34, C44, C45, C46 ],
[ C15, C25, C35, C45, C55, C56 ],
[ C16, C26, C36, C46, C56, C66 ],
])

@staticmethod
def CMatrixMonoclinic(C11,C12,C13, C15,
C22,C23, C25,
C33, C35,
C44, C46,
C55,
C66):
"""
Monoclinic cell - 13 constants
Order:
C11,C12,C13,C15,
C22,C23,C25,
C33,C35,
C44,C46,
C55,
C66
returns an array: [
[ C11, C12, C13, 0.0, C15, 0.0 ],
[ C12, C22, C23, 0.0, C25, 0.0 ],
[ C13, C23, C33, 0.0, C35, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, C46 ],
[ C15, C25, C35, 0.0, C55, 0.0 ],
[ 0.0, 0.0, 0.0, C46, 0.0, C66 ],
]
"""
return np.array([
[ C11, C12, C13, 0.0, C15, 0.0 ],
[ C12, C22, C23, 0.0, C25, 0.0 ],
[ C13, C23, C33, 0.0, C35, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, C46 ],
[ C15, C25, C35, 0.0, C55, 0.0 ],
[ 0.0, 0.0, 0.0, C46, 0.0, C66 ],
])

@staticmethod
def CMatrixOrthorhombic(C11,C12,C13,
C22,C23,
C33,C44,C55,C66):
"""
Orthorhombic cell - 9 constants
Order:
C11,C12,C13,
C22,C23,
C33,C44,C55,C66
returns an array: [
[ C11, C12, C13, 0.0, 0.0, 0.0 ],
[ C12, C22, C23, 0.0, 0.0, 0.0 ],
[ C13, C23, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C55, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C66 ],
]
"""
return np.array([
[ C11, C12, C13, 0.0, 0.0, 0.0 ],
[ C12, C22, C23, 0.0, 0.0, 0.0 ],
[ C13, C23, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C55, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C66 ],
])

@staticmethod
def CMatrixTrigonal7(C11,C12,C13,C14,C15,
C33,C44):
"""
Trigonal cell - 7 constants
Order
C11,C12,C13,C14,C15,
C33,C44
returns an array: [
[ C11, C12, C13, C14, C15, 0.0 ],
[ C12, C11, C13,-C14,-C15, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ C14,-C14, 0.0, C44, 0.0,-C15 ],
[ C15,-C15, 0.0, 0.0, C44, C14 ],
[ 0.0, 0.0, 0.0,-C15, C14, C66 ],
]
where C66 = 0.5*(C11-C12)
"""
C66 = 0.5*(C11-C12)
return np.array([
[ C11, C12, C13, C14, C15, 0.0 ],
[ C12, C11, C13,-C14,-C15, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ C14,-C14, 0.0, C44, 0.0,-C15 ],
[ C15,-C15, 0.0, 0.0, C44, C14 ],
[ 0.0, 0.0, 0.0,-C15, C14, C66 ],
])

@staticmethod
def CMatrixTrigonal6(C11,C12,C13,C14,
C33,C44):
"""
Trigonal cell - 6 constants
Order
C11,C12,C13,C14,
C33,C44
returns an array: [
[ C11, C12, C13, C14, 0.0, 0.0 ],
[ C12, C11, C13,-C14, 0.0, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ C14,-C14, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, C14 ],
[ 0.0, 0.0, 0.0, 0.0, C14, C66 ],
]
where C66 = 0.5*(C11-C12)
"""
C66 = 0.5*(C11-C12)
return np.array([
[ C11, C12, C13, C14, 0.0, 0.0 ],
[ C12, C11, C13,-C14, 0.0, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ C14,-C14, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, C14 ],
[ 0.0, 0.0, 0.0, 0.0, C14, C66 ],
])

@staticmethod
def CMatrixHexagonal(C11,C12,C13,
C33,C44):
"""
Hexagonal cell - 5 constants
Order:
C11,C12,C13,
C33,C44
returns an array: [
[ C11, C12, C13, 0.0, 0.0, 0.0 ],
[ C12, C11, C13, 0.0, 0.0, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C66 ],
]
where C66 = 0.5*(C11-C12)
"""
C66 = 0.5*(C11-C12)
return np.array([
[ C11, C12, C13, 0.0, 0.0, 0.0 ],
[ C12, C11, C13, 0.0, 0.0, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C66 ],
])

@staticmethod
def CMatrixTetragonal7(C11,C12,C13, C16,
C33,C44,C66):
"""
Tetragonal cell - 7 constants
Order:
C11,C12,C13, C16,
C33,C44,C66
returns an array: [
[ C11, C12, C13, 0.0, 0.0, C16 ],
[ C12, C11, C13, 0.0, 0.0,-C16 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ C16,-C16, 0.0, 0.0, 0.0, C66 ],
]
"""
return np.array([
[ C11, C12, C13, 0.0, 0.0, C16 ],
[ C12, C11, C13, 0.0, 0.0,-C16 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ C16,-C16, 0.0, 0.0, 0.0, C66 ],
])

@staticmethod
def CMatrixTetragonal6(C11,C12,C13,
C33,C44,C66):
"""
Tetragonal cell - 6 constants
Order:
C11,C12,C13,
C33,C44,C66
returns an array: [
[ C11, C12, C13, 0.0, 0.0, 0.0 ],
[ C12, C11, C13, 0.0, 0.0, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C66 ],
]
"""
return np.array([
[ C11, C12, C13, 0.0, 0.0, 0.0 ],
[ C12, C11, C13, 0.0, 0.0, 0.0 ],
[ C13, C13, C33, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C66 ],
])

@staticmethod
def CMatrixCubic(C11,C12,C44):
"""
Cubic cell - 3 constants
Order:
C11,C12,C44
returns an array: [
[ C11, C12, C12, 0.0, 0.0, 0.0 ],
[ C12, C11, C12, 0.0, 0.0, 0.0 ],
[ C12, C12, C11, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C44 ],
]
"""
return np.array([
[ C11, C12, C12, 0.0, 0.0, 0.0 ],
[ C12, C11, C12, 0.0, 0.0, 0.0 ],
[ C12, C12, C11, 0.0, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, C44, 0.0, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, C44, 0.0 ],
[ 0.0, 0.0, 0.0, 0.0, 0.0, C44 ],
])
16 changes: 16 additions & 0 deletions vasprun/elastic/__init__.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,16 @@
"""
Submodule for calculation of the elastic properties
from the stress--strain relation based on VASP's ```vasprun.xml``` file.

Module CMatrix contains a class ```ElasticTensor```
with all 9 different symmetry-based constrains 6x6 elastic tensors.

Module ```elastic``` contains two classes:
(i) ```Calculation``` containing the scipy-based engine
fitting the elastic tensor C to the stress-strain relation.
(ii) ```Properites``` with the methods deriving the elastic properties
from the elastic tensor (within different approximations)

Module ```readdata``` contains a class ```XMLParser``` extracting
both stress and strain from the ```vasprun.xml```
"""
Loading