From cc06fba34b8bcb969fc5ac4ac2156f4c67085433 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Thu, 19 Oct 2023 16:36:55 -0700 Subject: [PATCH 01/26] organize code a little bit. All material related files are nor in the materials folder. --- hexrd/__init__.py | 9 +- hexrd/materials/__init__.py | 0 hexrd/{ => materials}/crystallography.py | 0 hexrd/materials/jcpds.py | 296 +++++++++++++++++++++++ hexrd/{ => materials}/material.py | 9 +- hexrd/{ => materials}/mksupport.py | 4 +- hexrd/{ => materials}/spacegroup.py | 0 hexrd/{ => materials}/symbols.py | 0 hexrd/{ => materials}/symmetry.py | 0 hexrd/{ => materials}/unitcell.py | 7 +- hexrd/rotations.py | 6 +- hexrd/wppf/WPPF.py | 2 +- hexrd/wppf/phase.py | 6 +- hexrd/wppf/wppfsupport.py | 6 +- hexrd/xrdutil/utils.py | 2 +- 15 files changed, 326 insertions(+), 21 deletions(-) create mode 100644 hexrd/materials/__init__.py rename hexrd/{ => materials}/crystallography.py (100%) create mode 100644 hexrd/materials/jcpds.py rename hexrd/{ => materials}/material.py (99%) rename hexrd/{ => materials}/mksupport.py (99%) rename hexrd/{ => materials}/spacegroup.py (100%) rename hexrd/{ => materials}/symbols.py (100%) rename hexrd/{ => materials}/symmetry.py (100%) rename hexrd/{ => materials}/unitcell.py (99%) diff --git a/hexrd/__init__.py b/hexrd/__init__.py index 69e3be50d..e7346d721 100644 --- a/hexrd/__init__.py +++ b/hexrd/__init__.py @@ -1 +1,8 @@ -__path__ = __import__('pkgutil').extend_path(__path__, __name__) +from .materials import symbols +from .materials import symmetry +from .materials import spacegroup +from .materials import crystallography +from .materials import unitcell +from .materials import material +from .material import mksupport +from .materials import jcpds diff --git a/hexrd/materials/__init__.py b/hexrd/materials/__init__.py new file mode 100644 index 000000000..e69de29bb diff --git a/hexrd/crystallography.py b/hexrd/materials/crystallography.py similarity index 100% rename from hexrd/crystallography.py rename to hexrd/materials/crystallography.py diff --git a/hexrd/materials/jcpds.py b/hexrd/materials/jcpds.py new file mode 100644 index 000000000..91e2a2af9 --- /dev/null +++ b/hexrd/materials/jcpds.py @@ -0,0 +1,296 @@ +class JCPDS_extend(): + def __init__(self, filename=None): + if filename is None: + self.file = ' ' + self.name = ' ' + self.version = 0 + self.comments = '' + self.symmetry = '' + self.k0 = 0. + self.k0p = 0. # k0p at 298K + self.dk0dt = 0. + self.dk0pdt = 0. + self.thermal_expansion = 0. # alphat at 298K + self.thermal_expansion_dt = 0. + self.a0 = 0. + self.b0 = 0. + self.c0 = 0. + self.alpha0 = 0. + self.beta0 = 0. + self.gamma0 = 0. + self.v0 = 0. + else: + self.file = filename + self.dk0dt = 0. + self.dk0pdt = 0. + self.thermal_expansion_dt = 0. + self.b0 = 0. + self.c0 = 0. + self.alpha0 = 0. + self.beta0 = 0. + self.gamma0 = 0. + self.v0 = 0. + self.read_file(self.file) + self.a = 0. + self.b = 0. + self.c = 0. + self.alpha = 0. + self.beta = 0. + self.gamma = 0. + + self.version_status = '' + + def read_file(self, file): + """ + read a jcpds file + """ + self.file = file + # Construct base name = file without path and without extension + name = os.path.splitext(os.path.basename(self.file))[0] + self.name = name +# line = '', nd=0 + version = 0. + self.comments = [] + self.DiffLines = [] + + inp = open(file, 'r').readlines() +# my_list = [] # get all the text first and throw into my_list + + if inp[0][0] in ('2', '3', '4'): + version = int(inp[0]) # JCPDS version number + self.version = version + header = inp[1] # header + self.comments = header + + item = str.split(inp[2]) + crystal_system = int(item[0]) + if crystal_system == 1: + self.symmetry = 'cubic' + elif crystal_system == 2: + self.symmetry = 'hexagonal' + elif crystal_system == 3: + self.symmetry = 'tetragonal' + elif crystal_system == 4: + self.symmetry = 'orthorhombic' + elif crystal_system == 5: + self.symmetry = 'monoclinic' + elif crystal_system == 6: + self.symmetry = 'triclinic' + elif crystal_system == 7: + self.symmetry = 'manual' + # 1 cubic, 2 hexagonal, 3 tetragonal, 4 orthorhombic + # 5 monoclinic, 6 triclinic, 7 manual P, d-sp input + + k0 = float(item[1]) + k0p = float(item[2]) + self.k0 = k0 + self.k0p = k0p + + item = str.split(inp[3]) # line for unit-cell parameters + + if crystal_system == 1: # cubic + a = float(item[0]) + b = a + c = a + alpha = 90. + beta = 90. + gamma = 90. + elif crystal_system == 7: # P, d-sp input + a = float(item[0]) + b = a + c = a + alpha = 90. + beta = 90. + gamma = 90. + elif crystal_system == 2: # hexagonal + a = float(item[0]) + c = float(item[1]) + b = a + alpha = 90. + beta = 90. + gamma = 120. + elif crystal_system == 3: # tetragonal + a = float(item[0]) + c = float(item[1]) + b = a + alpha = 90. + beta = 90. + gamma = 90. + elif crystal_system == 4: # orthorhombic + a = float(item[0]) + b = float(item[1]) + c = float(item[2]) + alpha = 90. + beta = 90. + gamma = 90. + elif crystal_system == 5: # monoclinic + a = float(item[0]) + b = float(item[1]) + c = float(item[2]) + beta = float(item[3]) + alpha = 90. + gamma = 90. + elif crystal_system == 6: # triclinic + a = float(item[0]) + b = float(item[1]) + c = float(item[2]) + alpha = float(item[3]) + beta = float(item[4]) + gamma = float(item[5]) + + self.a0 = a + self.b0 = b + self.c0 = c + self.alpha0 = alpha + self.beta0 = beta + self.gamma0 = gamma + + item = str.split(inp[4]) + + if self.version == 3: + thermal_expansion = 0. + else: + thermal_expansion = float(item[0]) + self.thermal_expansion = thermal_expansion + + self.version_status = 'new' + + elif 'VERSION' in inp[0]: + jcpdsfile = open(file, 'r') + while True: + jcpdsline = jcpdsfile.readline() + if jcpdsline == '': + break + + jlinespl = jcpdsline.split() + + if jlinespl[0] == 'VERSION:': + version = int(jlinespl[1]) + self.version = version + + if jlinespl[0] == 'COMMENT:': + header = ' '.join(jlinespl[1:]) + self.comments = header + + if jlinespl[0] == 'K0:': + k0 = float(jlinespl[1]) + self.k0 = k0 + + if jlinespl[0] == 'K0P:': + k0p = float(jlinespl[1]) + self.k0p = k0p + + if jlinespl[0] == 'DK0DT:': + dk0dt = float(jlinespl[1]) + self.dk0dt = dk0dt + + if jlinespl[0] == 'DK0PDT:': + dk0pdt = float(jlinespl[1]) + self.dk0pdt = dk0pdt + + if jlinespl[0] == 'SYMMETRY:': + self.symmetry = jlinespl[1].lower() + + if jlinespl[0] == 'A:': + a = float(jlinespl[1]) + self.a0 = a + + if jlinespl[0] == 'B:': + b = float(jlinespl[1]) + self.b0 = b + + if jlinespl[0] == 'C:': + c = float(jlinespl[1]) + self.c0 = c + + if jlinespl[0] == 'ALPHA:': + alpha = float(jlinespl[1]) + self.alpha0 = alpha + + if jlinespl[0] == 'BETA:': + beta = float(jlinespl[1]) + self.beta0 = beta + + if jlinespl[0] == 'GAMMA:': + gamma = float(jlinespl[1]) + self.gamma0 = gamma + + if jlinespl[0] == 'VOLUME:': + v = float(jlinespl[1]) + self.v0 = v + + if jlinespl[0] == 'ALPHAT:': + alphat = float(jlinespl[1]) + self.thermal_expansion = alphat + + if jlinespl[0] == 'DALPHATDT:': + dalphatdt = float(jlinespl[1]) + self.thermal_expansion_dt = dalphatat + + if jlinespl[0] == 'DIHKL:': + pass + + if self.symmetry == 'cubic': + self.b0 = self.a0 + self.c0 = self.a0 + self.alpha0 = 90. + self.beta0 = 90. + self.gamma0 = 90. + elif self.symmetry == 'manual': + self.b0 = self.a0 + self.c0 = self.a0 + self.alpha0 = 90. + self.beta0 = 90. + self.gamma0 = 90. + elif self.symmetry == 'hexagonal' or self.symmetry == 'trigonal': + self.b0 = self.a0 + self.alpha0 = 90. + self.beta0 = 90. + self.gamma0 = 120. + elif self.symmetry == 'tetragonal': + self.b0 = self.a0 + self.alpha0 = 90. + self.beta0 = 90. + self.gamma0 = 90. + elif self.symmetry == 'orthorhombic': + self.alpha0 = 90. + self.beta0 = 90. + self.gamma0 = 90. + elif self.symmetry == 'monoclinic': + self.alpha0 = 90. + self.gamma0 = 90. + #elif self.symmetry == 'triclinic': + + jcpdsfile.close() + + self.version_status = 'new' + + else: + self.version_status = 'old' + + if self.version_status == 'new': + self.a = self.a0 + self.b = self.b0 + self.c = self.c0 + self.alpha = self.alpha0 + self.beta = self.beta0 + self.gamma = self.gamma0 + self.v = self.v0 + + def calc_pressure(self, volume=None, temperature=298): + '''calculate the pressure given the volume + and temperature using the third order + birch-murnaghan equation of state. + ''' + if volume is None: + return 0 + else: + v0 = self.v0 + alpha0 = self.thermal_expansion + alpha1 = self.thermal_expansion_dt + + if temperature == 298: + vt = v0 + else: + vt = v0*np.exp(alpha0*(temperature-298) + +0.5*alpha1*(temperature**2-298**2)) diff --git a/hexrd/material.py b/hexrd/materials/material.py similarity index 99% rename from hexrd/material.py rename to hexrd/materials/material.py index 029cbe52b..2ec2965e9 100644 --- a/hexrd/material.py +++ b/hexrd/materials/material.py @@ -34,7 +34,7 @@ from configparser import SafeConfigParser as Parser import numpy -from hexrd.crystallography import PlaneData as PData +from hexrd.materials.crystallography import PlaneData as PData from hexrd.valunits import valWUnit from hexrd import unitcell from hexrd.constants import ptable, ptableinverse, chargestate @@ -46,9 +46,10 @@ from CifFile import ReadCif import h5py from warnings import warn -from hexrd.mksupport import Write2H5File -from hexrd.symbols import xtal_sys_dict -from hexrd.symbols import Hall_to_sgnum, HM_to_sgnum +from hexrd.materials.mksupport import Write2H5File +from hexrd.materials.symbols import ( + xtal_sys_dict,Hall_to_sgnum, + HM_to_sgnum) from hexrd.utils.compatibility import h5py_read_string from hexrd.fitting.peakfunctions import _unit_gaussian diff --git a/hexrd/mksupport.py b/hexrd/materials/mksupport.py similarity index 99% rename from hexrd/mksupport.py rename to hexrd/materials/mksupport.py index 11d97199e..dc32b1c73 100644 --- a/hexrd/mksupport.py +++ b/hexrd/materials/mksupport.py @@ -1,4 +1,4 @@ -from hexrd.symbols import pstr_Elements, sitesym, \ +from hexrd.materials.symbols import pstr_Elements, sitesym, \ tworig, PrintPossibleSG, TRIG, pstr_spacegroup,\ pstr_mkxtal import h5py @@ -6,7 +6,7 @@ import numpy as np import datetime import getpass -from hexrd.unitcell import _StiffnessDict, _pgDict +from hexrd.materials.unitcell import _StiffnessDict, _pgDict def mk(filename, xtalname): diff --git a/hexrd/spacegroup.py b/hexrd/materials/spacegroup.py similarity index 100% rename from hexrd/spacegroup.py rename to hexrd/materials/spacegroup.py diff --git a/hexrd/symbols.py b/hexrd/materials/symbols.py similarity index 100% rename from hexrd/symbols.py rename to hexrd/materials/symbols.py diff --git a/hexrd/symmetry.py b/hexrd/materials/symmetry.py similarity index 100% rename from hexrd/symmetry.py rename to hexrd/materials/symmetry.py diff --git a/hexrd/unitcell.py b/hexrd/materials/unitcell.py similarity index 99% rename from hexrd/unitcell.py rename to hexrd/materials/unitcell.py index 094981dc2..3887b6439 100644 --- a/hexrd/unitcell.py +++ b/hexrd/materials/unitcell.py @@ -1,8 +1,9 @@ import importlib.resources import numpy as np from hexrd import constants -from hexrd import symmetry, symbols -from hexrd.spacegroup import Allowed_HKLs +from hexrd import symbols +from hexrd import symmetry +from hexrd import spacegroup from hexrd.ipfcolor import sphere_sector, colorspace from hexrd.valunits import valWUnit import hexrd.resources @@ -1033,7 +1034,7 @@ def getHKLs(self, dmin): for ik in np.arange(kmax, kmin, -1) for il in np.arange(lmax, lmin, -1)]) - hkl_allowed = Allowed_HKLs(self.sgnum, hkllist) + hkl_allowed = spacegroup.Allowed_HKLs(self.sgnum, hkllist) hkl = [] dsp = [] diff --git a/hexrd/rotations.py b/hexrd/rotations.py index c419cc29c..464e9a0af 100644 --- a/hexrd/rotations.py +++ b/hexrd/rotations.py @@ -49,7 +49,7 @@ columnNorm, unitVector, \ skewMatrixOfVector, findDuplicateVectors, \ multMatArray, nullSpace -from hexrd import symmetry + from hexrd.utils.decorators import numba_njit_if_available # ============================================================================= @@ -445,7 +445,7 @@ def quatAverageCluster(q_in, qsym): q_in) # second, re-cast to FR - qrot = symmetry.toFundamentalRegion(qrot.squeeze(), crysSym=qsym) + qrot = toFundamentalRegion(qrot.squeeze(), crysSym=qsym) # compute arithmetic average q_bar = unitVector(average(qrot, axis=1).reshape(4, 1)) @@ -456,7 +456,7 @@ def quatAverageCluster(q_in, qsym): q_bar) # re-map - q_bar = symmetry.toFundamentalRegion(q_bar, crysSym=qsym) + q_bar = toFundamentalRegion(q_bar, crysSym=qsym) return q_bar diff --git a/hexrd/wppf/WPPF.py b/hexrd/wppf/WPPF.py index ceeb5fe97..171a09463 100644 --- a/hexrd/wppf/WPPF.py +++ b/hexrd/wppf/WPPF.py @@ -19,7 +19,7 @@ # ------------- from hexrd import constants from hexrd.imageutil import snip1d_quad -from hexrd.material import Material +from hexrd.materials.material import Material from hexrd.utils.multiprocess_generic import GenericMultiprocessing from hexrd.valunits import valWUnit from hexrd.wppf.peakfunctions import ( diff --git a/hexrd/wppf/phase.py b/hexrd/wppf/phase.py index 8ef597208..7412ad61f 100644 --- a/hexrd/wppf/phase.py +++ b/hexrd/wppf/phase.py @@ -1,9 +1,9 @@ import numpy as np from hexrd.valunits import valWUnit -from hexrd.spacegroup import Allowed_HKLs, SpaceGroup +from hexrd.materials.spacegroup import Allowed_HKLs, SpaceGroup from hexrd import symmetry, symbols, constants -from hexrd.material import Material -from hexrd.unitcell import _rqpDict +from hexrd.materials.material import Material +from hexrd.materials.unitcell import _rqpDict from hexrd.wppf import wppfsupport from hexrd.wppf.xtal import _calc_dspacing, _get_tth, _calcxrsf,\ _calc_extinction_factor, _calc_absorption_factor diff --git a/hexrd/wppf/wppfsupport.py b/hexrd/wppf/wppfsupport.py index 07e320173..4a52b6aea 100644 --- a/hexrd/wppf/wppfsupport.py +++ b/hexrd/wppf/wppfsupport.py @@ -32,12 +32,12 @@ classes are put here to minimize code duplication. Some examples include initialize background, generate_default_parameter list etc. """ -from hexrd.symbols import pstr_spacegroup +from hexrd.materials.symbols import pstr_spacegroup from hexrd.wppf.parameters import Parameters from lmfit import Parameters as Parameters_lmfit from hexrd.wppf.phase import Phases_LeBail, Phases_Rietveld -from hexrd.material import Material -from hexrd.unitcell import _rqpDict +from hexrd.materials.material import Material +from hexrd.materials.unitcell import _rqpDict import hexrd import numpy as np from hexrd import constants diff --git a/hexrd/xrdutil/utils.py b/hexrd/xrdutil/utils.py index 376fbf40d..139c22d31 100644 --- a/hexrd/xrdutil/utils.py +++ b/hexrd/xrdutil/utils.py @@ -32,7 +32,7 @@ from hexrd import matrixutil as mutil from hexrd import gridutil as gutil -from hexrd.crystallography import processWavelength, PlaneData +from hexrd.materials.crystallography import processWavelength, PlaneData from hexrd.transforms import xf from hexrd.transforms import xfcapi From 3d174d66a759d1b7d2d336153c014666c21c38e2 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Thu, 19 Oct 2023 16:57:12 -0700 Subject: [PATCH 02/26] compute pressure. --- hexrd/materials/jcpds.py | 25 +++++++++++++++++++++---- 1 file changed, 21 insertions(+), 4 deletions(-) diff --git a/hexrd/materials/jcpds.py b/hexrd/materials/jcpds.py index 91e2a2af9..40def56cd 100644 --- a/hexrd/materials/jcpds.py +++ b/hexrd/materials/jcpds.py @@ -1,3 +1,6 @@ +import os +import numpy as np + class JCPDS_extend(): def __init__(self, filename=None): if filename is None: @@ -277,7 +280,7 @@ def read_file(self, file): self.gamma = self.gamma0 self.v = self.v0 - def calc_pressure(self, volume=None, temperature=298): + def calc_pressure(self, volume=None, temperature=None): '''calculate the pressure given the volume and temperature using the third order birch-murnaghan equation of state. @@ -289,8 +292,22 @@ def calc_pressure(self, volume=None, temperature=298): alpha0 = self.thermal_expansion alpha1 = self.thermal_expansion_dt - if temperature == 298: + beta0 = self.dk0dt + beta1 = self.dk0pdt + + k0 = self.k0 + k0p = self.k0p + if temperature is None: vt = v0 + kt = k0 + ktp = k0p else: - vt = v0*np.exp(alpha0*(temperature-298) - +0.5*alpha1*(temperature**2-298**2)) + delT = (temperature-298) + delT2 = (temperature**2-298**2) + vt = v0*np.exp(alpha0*delT + +0.5*alpha1*delT2) + kt = k0 + beta0*delT + ktp = k0p + beta1*delT + f = 0.5*((vt/volume)**(2./3.) - 1) + + return 3.0*kt*f*(1 - 1.5*(4 - ktp)*f)*(1 + 2*f)**2.5 From 720bcd1e3e0746c415eef08e7c2dac07a26d5182 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Fri, 20 Oct 2023 11:18:51 -0700 Subject: [PATCH 03/26] more import path corrections. --- hexrd/instrument/detector.py | 2 +- hexrd/instrument/hedm_instrument.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/hexrd/instrument/detector.py b/hexrd/instrument/detector.py index fbd3763e1..58b5aa703 100644 --- a/hexrd/instrument/detector.py +++ b/hexrd/instrument/detector.py @@ -9,7 +9,7 @@ from hexrd import matrixutil as mutil from hexrd import xrdutil -from hexrd.crystallography import PlaneData +from hexrd.materials.crystallography import PlaneData from hexrd.transforms.xfcapi import ( detectorXYToGvec, diff --git a/hexrd/instrument/hedm_instrument.py b/hexrd/instrument/hedm_instrument.py index 1ecbdb5ae..fe4f76386 100644 --- a/hexrd/instrument/hedm_instrument.py +++ b/hexrd/instrument/hedm_instrument.py @@ -68,7 +68,7 @@ unitRowVector, ) from hexrd import xrdutil -from hexrd.crystallography import PlaneData +from hexrd.materials.crystallography import PlaneData from hexrd import constants as ct from hexrd.rotations import ( angleAxisOfRotMat, From e9eb841fdeeca5c55b25d76e0113a8b5618788de Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Fri, 20 Oct 2023 16:45:55 -0700 Subject: [PATCH 04/26] finish volume computation. --- hexrd/materials/jcpds.py | 103 ++++++++++++++++++++++++++++++--------- 1 file changed, 81 insertions(+), 22 deletions(-) diff --git a/hexrd/materials/jcpds.py b/hexrd/materials/jcpds.py index 40def56cd..1c3e3dbca 100644 --- a/hexrd/materials/jcpds.py +++ b/hexrd/materials/jcpds.py @@ -32,8 +32,8 @@ def __init__(self, filename=None): self.alpha0 = 0. self.beta0 = 0. self.gamma0 = 0. - self.v0 = 0. self.read_file(self.file) + self.v0 = self.calc_volume_unitcell() self.a = 0. self.b = 0. self.c = 0. @@ -278,7 +278,21 @@ def read_file(self, file): self.alpha = self.alpha0 self.beta = self.beta0 self.gamma = self.gamma0 - self.v = self.v0 + self.v = self.calc_volume_unitcell() + + def calc_volume_unitcell(self): + '''calculate volume of the unitcell + Ref: Introduction to conventional TEM + Marc De Graef, Appendix 1 pg. 662 + ''' + ca = np.cos(np.radians(self.alpha0)) + cb = np.cos(np.radians(self.beta0)) + cg = np.cos(np.radians(self.gamma0)) + + v0 = self.a0*self.b0*self.c0 + f = np.sqrt(1 - ca**2 - cb**2 - cg**2 + +2*ca*cb*cg) + return v0*f def calc_pressure(self, volume=None, temperature=None): '''calculate the pressure given the volume @@ -288,26 +302,71 @@ def calc_pressure(self, volume=None, temperature=None): if volume is None: return 0 else: - v0 = self.v0 - alpha0 = self.thermal_expansion - alpha1 = self.thermal_expansion_dt - - beta0 = self.dk0dt - beta1 = self.dk0pdt - - k0 = self.k0 - k0p = self.k0p - if temperature is None: - vt = v0 - kt = k0 - ktp = k0p - else: - delT = (temperature-298) - delT2 = (temperature**2-298**2) - vt = v0*np.exp(alpha0*delT - +0.5*alpha1*delT2) - kt = k0 + beta0*delT - ktp = k0p + beta1*delT + vt = self.vt(temperature=temperature) + kt = self.kt(temperature=temperature) + ktp = self.ktp(temperature=temperature) f = 0.5*((vt/volume)**(2./3.) - 1) return 3.0*kt*f*(1 - 1.5*(4 - ktp)*f)*(1 + 2*f)**2.5 + + def calc_volume(self, pressure=None, temperature=None): + '''solve for volume in the birch-murnaghan EoS to + compute the volume. this number will be propagated + to the Material object as updated lattice constants. + ''' + vt = self.vt(temperature=temperature) + kt = self.kt(temperature=temperature) + ktp = self.ktp(temperature=temperature) + + if pressure is None: + return vt + else: + alpha = 0.75*(ktp - 4) + p = np.zeros([10,]) + p[0] = 1.0 + p[2] = (1 - 2*alpha)/alpha + p[4] = (alpha-1)/alpha + p[9] = -2*pressure/3/kt/alpha + res = np.roots(p) + res = res[np.isreal(res)] + res = 1/np.real(res)**3 + + mask = np.logical_and(res >= 0., res <= 1.0) + return res[mask] * vt + + def vt(self, temperature=None): + '''calculate volume at high + temperature + ''' + alpha0 = self.thermal_expansion + alpha1 = self.thermal_expansion_dt + if temperature is None: + vt = self.v0 + else: + delT = (temperature-298) + delT2 = (temperature**2-298**2) + vt = self.v0*np.exp(alpha0*delT + +0.5*alpha1*delT2) + return vt + + def kt(self, temperature=None): + '''calculate bulk modulus for + high temperature + ''' + k0 = self.k0 + if temperature is None: + return k0 + else: + delT = (temperature-298) + return k0 + self.dk0dt*delT + + def ktp(self, temperature=None): + '''calculate bulk modulus derivative + for high temperature + ''' + k0p = self.k0p + if temperature is None: + return k0p + else: + delT = (temperature-298) + return k0p + self.dk0pdt*delT \ No newline at end of file From 290935d064a252fc0e03dd0664801c5783f5594a Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Fri, 20 Oct 2023 17:14:51 -0700 Subject: [PATCH 05/26] compute the scaling factor for lattice parameter lengths. --- hexrd/materials/jcpds.py | 89 +++++++++++++++++++++++----------------- 1 file changed, 52 insertions(+), 37 deletions(-) diff --git a/hexrd/materials/jcpds.py b/hexrd/materials/jcpds.py index 1c3e3dbca..259495812 100644 --- a/hexrd/materials/jcpds.py +++ b/hexrd/materials/jcpds.py @@ -294,6 +294,43 @@ def calc_volume_unitcell(self): +2*ca*cb*cg) return v0*f + def vt(self, temperature=None): + '''calculate volume at high + temperature + ''' + alpha0 = self.thermal_expansion + alpha1 = self.thermal_expansion_dt + if temperature is None: + vt = self.v0 + else: + delT = (temperature-298) + delT2 = (temperature**2-298**2) + vt = self.v0*np.exp(alpha0*delT + +0.5*alpha1*delT2) + return vt + + def kt(self, temperature=None): + '''calculate bulk modulus for + high temperature + ''' + k0 = self.k0 + if temperature is None: + return k0 + else: + delT = (temperature-298) + return k0 + self.dk0dt*delT + + def ktp(self, temperature=None): + '''calculate bulk modulus derivative + for high temperature + ''' + k0p = self.k0p + if temperature is None: + return k0p + else: + delT = (temperature-298) + return k0p + self.dk0pdt*delT + def calc_pressure(self, volume=None, temperature=None): '''calculate the pressure given the volume and temperature using the third order @@ -332,41 +369,19 @@ def calc_volume(self, pressure=None, temperature=None): res = 1/np.real(res)**3 mask = np.logical_and(res >= 0., res <= 1.0) - return res[mask] * vt - - def vt(self, temperature=None): - '''calculate volume at high - temperature + res = res[mask] + if len(res) != 1: + msg = (f'more than one physically ' + f'reasonable solution found!') + raise ValueError(msg) + return res[0] * vt + + def calc_lp_factor(self, pressure=None, temperature=None): + '''calculate the factor to multiply the lattice + constants by. only the lengths will be modified, the + angles will be kept constant. ''' - alpha0 = self.thermal_expansion - alpha1 = self.thermal_expansion_dt - if temperature is None: - vt = self.v0 - else: - delT = (temperature-298) - delT2 = (temperature**2-298**2) - vt = self.v0*np.exp(alpha0*delT - +0.5*alpha1*delT2) - return vt - - def kt(self, temperature=None): - '''calculate bulk modulus for - high temperature - ''' - k0 = self.k0 - if temperature is None: - return k0 - else: - delT = (temperature-298) - return k0 + self.dk0dt*delT - - def ktp(self, temperature=None): - '''calculate bulk modulus derivative - for high temperature - ''' - k0p = self.k0p - if temperature is None: - return k0p - else: - delT = (temperature-298) - return k0p + self.dk0pdt*delT \ No newline at end of file + vt = self.vt(temperature=temperature) + vpt = self.calc_volume(pressure=pressure, + temperature=temperature) + return (vpt/vt)**(1./3.) \ No newline at end of file From 994f5a47eb19e2de614de7cdbcfbf8c0b78a8444 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Fri, 20 Oct 2023 17:22:13 -0700 Subject: [PATCH 06/26] calculate lattice parameter as a function of pressure and temperature using the Birch-Murnaghan EOS. --- hexrd/materials/jcpds.py | 13 ++++++++++++- 1 file changed, 12 insertions(+), 1 deletion(-) diff --git a/hexrd/materials/jcpds.py b/hexrd/materials/jcpds.py index 259495812..7b97859d6 100644 --- a/hexrd/materials/jcpds.py +++ b/hexrd/materials/jcpds.py @@ -384,4 +384,15 @@ def calc_lp_factor(self, pressure=None, temperature=None): vt = self.vt(temperature=temperature) vpt = self.calc_volume(pressure=pressure, temperature=temperature) - return (vpt/vt)**(1./3.) \ No newline at end of file + return (vpt/vt)**(1./3.) + + def calc_lp_at_PT(self, pressure=None, temperature=None): + '''calculate the lattice parameters for a given + pressure and temperature using the BM EoS. This + is the main function which will be called from + the GUI. + ''' + f = self.calc_lp_factor(pressure=pressure, + temperature=temperature) + return np.array([f*self.a0, f*self.b0, f*self.c0, + self.alpha0, self.beta0, self.gamma0,]) \ No newline at end of file From bcc809f6a50981bae28ff1209cccef244994a911 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 2 Nov 2023 16:45:27 -0500 Subject: [PATCH 07/26] Rename "materials" package back to "material" There are still some compatibility issues, but not as many. Signed-off-by: Patrick Avery --- hexrd/__init__.py | 13 ++++++------- hexrd/config/material.py | 2 +- hexrd/instrument/detector.py | 2 +- hexrd/instrument/hedm_instrument.py | 2 +- hexrd/material/__init__.py | 7 +++++++ hexrd/{materials => material}/crystallography.py | 0 hexrd/{materials => material}/jcpds.py | 0 hexrd/{materials => material}/material.py | 11 +++++------ hexrd/{materials => material}/mksupport.py | 4 ++-- hexrd/{materials => material}/spacegroup.py | 3 ++- hexrd/{materials => material}/symbols.py | 0 hexrd/{materials => material}/symmetry.py | 0 hexrd/{materials => material}/unitcell.py | 4 +--- hexrd/materials/__init__.py | 0 hexrd/wppf/WPPF.py | 2 +- hexrd/wppf/phase.py | 10 +++++----- hexrd/wppf/wppfsupport.py | 6 +++--- hexrd/xrdutil/utils.py | 2 +- 18 files changed, 36 insertions(+), 32 deletions(-) create mode 100644 hexrd/material/__init__.py rename hexrd/{materials => material}/crystallography.py (100%) rename hexrd/{materials => material}/jcpds.py (100%) rename hexrd/{materials => material}/material.py (99%) rename hexrd/{materials => material}/mksupport.py (99%) rename hexrd/{materials => material}/spacegroup.py (99%) rename hexrd/{materials => material}/symbols.py (100%) rename hexrd/{materials => material}/symmetry.py (100%) rename hexrd/{materials => material}/unitcell.py (99%) delete mode 100644 hexrd/materials/__init__.py diff --git a/hexrd/__init__.py b/hexrd/__init__.py index e7346d721..f9da2f348 100644 --- a/hexrd/__init__.py +++ b/hexrd/__init__.py @@ -1,8 +1,7 @@ -from .materials import symbols -from .materials import symmetry -from .materials import spacegroup -from .materials import crystallography -from .materials import unitcell -from .materials import material +from .material import crystallography +from .material import jcpds from .material import mksupport -from .materials import jcpds +from .material import spacegroup +from .material import symbols +from .material import symmetry +from .material import unitcell diff --git a/hexrd/config/material.py b/hexrd/config/material.py index fc73dfa3a..24bbce64a 100644 --- a/hexrd/config/material.py +++ b/hexrd/config/material.py @@ -109,5 +109,5 @@ def beam_energy(self): def beam_energy(self, x): if not isinstance(x, valWUnit): x = valWUnit("beam energy", "energy", x, "keV") - for matl in self.materials.values(): + for matl in self.material.values(): matl.beamEnergy = x diff --git a/hexrd/instrument/detector.py b/hexrd/instrument/detector.py index 58b5aa703..6cfd8d834 100644 --- a/hexrd/instrument/detector.py +++ b/hexrd/instrument/detector.py @@ -9,7 +9,7 @@ from hexrd import matrixutil as mutil from hexrd import xrdutil -from hexrd.materials.crystallography import PlaneData +from hexrd.material.crystallography import PlaneData from hexrd.transforms.xfcapi import ( detectorXYToGvec, diff --git a/hexrd/instrument/hedm_instrument.py b/hexrd/instrument/hedm_instrument.py index fe4f76386..8f5b7e1ca 100644 --- a/hexrd/instrument/hedm_instrument.py +++ b/hexrd/instrument/hedm_instrument.py @@ -68,7 +68,7 @@ unitRowVector, ) from hexrd import xrdutil -from hexrd.materials.crystallography import PlaneData +from hexrd.material.crystallography import PlaneData from hexrd import constants as ct from hexrd.rotations import ( angleAxisOfRotMat, diff --git a/hexrd/material/__init__.py b/hexrd/material/__init__.py new file mode 100644 index 000000000..8ce0c9625 --- /dev/null +++ b/hexrd/material/__init__.py @@ -0,0 +1,7 @@ +from .material import ( + _angstroms, + _kev, + load_materials_hdf5, + Material, + save_materials_hdf5, +) diff --git a/hexrd/materials/crystallography.py b/hexrd/material/crystallography.py similarity index 100% rename from hexrd/materials/crystallography.py rename to hexrd/material/crystallography.py diff --git a/hexrd/materials/jcpds.py b/hexrd/material/jcpds.py similarity index 100% rename from hexrd/materials/jcpds.py rename to hexrd/material/jcpds.py diff --git a/hexrd/materials/material.py b/hexrd/material/material.py similarity index 99% rename from hexrd/materials/material.py rename to hexrd/material/material.py index 2ec2965e9..2313609d7 100644 --- a/hexrd/materials/material.py +++ b/hexrd/material/material.py @@ -34,11 +34,10 @@ from configparser import SafeConfigParser as Parser import numpy -from hexrd.materials.crystallography import PlaneData as PData +from hexrd.material.crystallography import PlaneData as PData +from hexrd.material import symmetry, unitcell from hexrd.valunits import valWUnit -from hexrd import unitcell from hexrd.constants import ptable, ptableinverse, chargestate -from hexrd import symmetry import copy from os import path @@ -46,9 +45,9 @@ from CifFile import ReadCif import h5py from warnings import warn -from hexrd.materials.mksupport import Write2H5File -from hexrd.materials.symbols import ( - xtal_sys_dict,Hall_to_sgnum, +from hexrd.material.mksupport import Write2H5File +from hexrd.material.symbols import ( + xtal_sys_dict,Hall_to_sgnum, HM_to_sgnum) from hexrd.utils.compatibility import h5py_read_string from hexrd.fitting.peakfunctions import _unit_gaussian diff --git a/hexrd/materials/mksupport.py b/hexrd/material/mksupport.py similarity index 99% rename from hexrd/materials/mksupport.py rename to hexrd/material/mksupport.py index dc32b1c73..b575cf66b 100644 --- a/hexrd/materials/mksupport.py +++ b/hexrd/material/mksupport.py @@ -1,4 +1,4 @@ -from hexrd.materials.symbols import pstr_Elements, sitesym, \ +from hexrd.material.symbols import pstr_Elements, sitesym, \ tworig, PrintPossibleSG, TRIG, pstr_spacegroup,\ pstr_mkxtal import h5py @@ -6,7 +6,7 @@ import numpy as np import datetime import getpass -from hexrd.materials.unitcell import _StiffnessDict, _pgDict +from hexrd.material.unitcell import _StiffnessDict, _pgDict def mk(filename, xtalname): diff --git a/hexrd/materials/spacegroup.py b/hexrd/material/spacegroup.py similarity index 99% rename from hexrd/materials/spacegroup.py rename to hexrd/material/spacegroup.py index c1a559392..3efb7212c 100644 --- a/hexrd/materials/spacegroup.py +++ b/hexrd/material/spacegroup.py @@ -72,7 +72,8 @@ from collections import OrderedDict from math import sqrt, floor -from hexrd import symbols, constants, symmetry +from hexrd import constants +from hexrd.material import symbols, symmetry import numpy as np # __all__ = ['SpaceGroup'] diff --git a/hexrd/materials/symbols.py b/hexrd/material/symbols.py similarity index 100% rename from hexrd/materials/symbols.py rename to hexrd/material/symbols.py diff --git a/hexrd/materials/symmetry.py b/hexrd/material/symmetry.py similarity index 100% rename from hexrd/materials/symmetry.py rename to hexrd/material/symmetry.py diff --git a/hexrd/materials/unitcell.py b/hexrd/material/unitcell.py similarity index 99% rename from hexrd/materials/unitcell.py rename to hexrd/material/unitcell.py index 3887b6439..2acbb22b3 100644 --- a/hexrd/materials/unitcell.py +++ b/hexrd/material/unitcell.py @@ -1,9 +1,7 @@ import importlib.resources import numpy as np from hexrd import constants -from hexrd import symbols -from hexrd import symmetry -from hexrd import spacegroup +from hexrd.material import spacegroup, symbols, symmetry from hexrd.ipfcolor import sphere_sector, colorspace from hexrd.valunits import valWUnit import hexrd.resources diff --git a/hexrd/materials/__init__.py b/hexrd/materials/__init__.py deleted file mode 100644 index e69de29bb..000000000 diff --git a/hexrd/wppf/WPPF.py b/hexrd/wppf/WPPF.py index 171a09463..ceeb5fe97 100644 --- a/hexrd/wppf/WPPF.py +++ b/hexrd/wppf/WPPF.py @@ -19,7 +19,7 @@ # ------------- from hexrd import constants from hexrd.imageutil import snip1d_quad -from hexrd.materials.material import Material +from hexrd.material import Material from hexrd.utils.multiprocess_generic import GenericMultiprocessing from hexrd.valunits import valWUnit from hexrd.wppf.peakfunctions import ( diff --git a/hexrd/wppf/phase.py b/hexrd/wppf/phase.py index 7412ad61f..ce6a878bd 100644 --- a/hexrd/wppf/phase.py +++ b/hexrd/wppf/phase.py @@ -1,9 +1,9 @@ import numpy as np from hexrd.valunits import valWUnit -from hexrd.materials.spacegroup import Allowed_HKLs, SpaceGroup +from hexrd.material.spacegroup import Allowed_HKLs, SpaceGroup from hexrd import symmetry, symbols, constants -from hexrd.materials.material import Material -from hexrd.materials.unitcell import _rqpDict +from hexrd.material import Material +from hexrd.material.unitcell import _rqpDict from hexrd.wppf import wppfsupport from hexrd.wppf.xtal import _calc_dspacing, _get_tth, _calcxrsf,\ _calc_extinction_factor, _calc_absorption_factor @@ -21,7 +21,7 @@ class Material_LeBail: 09/14/2020 SS 1.1 class can now be initialized using a material.Material class instance >> @DETAILS: Material_LeBail class is a stripped down version of the - materials.Material class.this is done to keep the class lightweight + material.Material class.this is done to keep the class lightweight and make sure only the information necessary for the lebail fit is kept ========================================================================================= @@ -605,7 +605,7 @@ class Material_Rietveld: >> @DATE: 05/18/2020 SS 1.0 original 02/01/2021 SS 1.1 class can now be initialized using a material.Material class instance - >> @DETAILS: Material_LeBail class is a stripped down version of the materials.Material + >> @DETAILS: Material_LeBail class is a stripped down version of the material.Material class.this is done to keep the class lightweight and make sure only the information necessary for the Rietveld fit is kept =========================================================================================== diff --git a/hexrd/wppf/wppfsupport.py b/hexrd/wppf/wppfsupport.py index 4a52b6aea..0720e56f7 100644 --- a/hexrd/wppf/wppfsupport.py +++ b/hexrd/wppf/wppfsupport.py @@ -32,12 +32,12 @@ classes are put here to minimize code duplication. Some examples include initialize background, generate_default_parameter list etc. """ -from hexrd.materials.symbols import pstr_spacegroup +from hexrd.material.symbols import pstr_spacegroup from hexrd.wppf.parameters import Parameters from lmfit import Parameters as Parameters_lmfit from hexrd.wppf.phase import Phases_LeBail, Phases_Rietveld -from hexrd.materials.material import Material -from hexrd.materials.unitcell import _rqpDict +from hexrd.material import Material +from hexrd.material.unitcell import _rqpDict import hexrd import numpy as np from hexrd import constants diff --git a/hexrd/xrdutil/utils.py b/hexrd/xrdutil/utils.py index 139c22d31..768e05ddc 100644 --- a/hexrd/xrdutil/utils.py +++ b/hexrd/xrdutil/utils.py @@ -32,7 +32,7 @@ from hexrd import matrixutil as mutil from hexrd import gridutil as gutil -from hexrd.materials.crystallography import processWavelength, PlaneData +from hexrd.material.crystallography import processWavelength, PlaneData from hexrd.transforms import xf from hexrd.transforms import xfcapi From bf95d2149bfdda75adc218874ba24d6ebafead91 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Fri, 3 Nov 2023 14:08:46 -0500 Subject: [PATCH 08/26] Fix mistaken renaming Signed-off-by: Patrick Avery --- hexrd/config/material.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/hexrd/config/material.py b/hexrd/config/material.py index 24bbce64a..fc73dfa3a 100644 --- a/hexrd/config/material.py +++ b/hexrd/config/material.py @@ -109,5 +109,5 @@ def beam_energy(self): def beam_energy(self, x): if not isinstance(x, valWUnit): x = valWUnit("beam energy", "energy", x, "keV") - for matl in self.material.values(): + for matl in self.materials.values(): matl.beamEnergy = x From 993f1d1024820306b24de2dd957c0b7b9dfee476 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Tue, 7 Nov 2023 10:19:06 -0600 Subject: [PATCH 09/26] Make changes needed to work with HEXRDGUI Signed-off-by: Patrick Avery --- hexrd/material/jcpds.py | 148 +++++++++++++++++++++++++--------------- 1 file changed, 92 insertions(+), 56 deletions(-) diff --git a/hexrd/material/jcpds.py b/hexrd/material/jcpds.py index 7b97859d6..0b515dccb 100644 --- a/hexrd/material/jcpds.py +++ b/hexrd/material/jcpds.py @@ -1,47 +1,36 @@ import os import numpy as np + class JCPDS_extend(): def __init__(self, filename=None): - if filename is None: - self.file = ' ' - self.name = ' ' - self.version = 0 - self.comments = '' - self.symmetry = '' - self.k0 = 0. - self.k0p = 0. # k0p at 298K - self.dk0dt = 0. - self.dk0pdt = 0. - self.thermal_expansion = 0. # alphat at 298K - self.thermal_expansion_dt = 0. - self.a0 = 0. - self.b0 = 0. - self.c0 = 0. - self.alpha0 = 0. - self.beta0 = 0. - self.gamma0 = 0. - self.v0 = 0. - else: + self.a0 = 0 + self.b0 = 0 + self.c0 = 0 + self.alpha0 = 0 + self.beta0 = 0 + self.gamma0 = 0 + self.v0 = 0 + + self.k0 = 100 + self.k0p = 5 # k0p at 298K + + self.dk0dt = 0 + self.dk0pdt = 0 + + self.symmetry = '' + self.thermal_expansion = 0 # alphat at 298K + self.thermal_expansion_dt = 0 + + self.file = ' ' + self.name = ' ' + self.version = 0 + self.comments = '' + + if filename: self.file = filename - self.dk0dt = 0. - self.dk0pdt = 0. - self.thermal_expansion_dt = 0. - self.b0 = 0. - self.c0 = 0. - self.alpha0 = 0. - self.beta0 = 0. - self.gamma0 = 0. self.read_file(self.file) - self.v0 = self.calc_volume_unitcell() - self.a = 0. - self.b = 0. - self.c = 0. - self.alpha = 0. - self.beta = 0. - self.gamma = 0. - - self.version_status = '' + self.update_v0() def read_file(self, file): """ @@ -56,6 +45,8 @@ def read_file(self, file): self.comments = [] self.DiffLines = [] + version_status = '' + inp = open(file, 'r').readlines() # my_list = [] # get all the text first and throw into my_list @@ -156,7 +147,7 @@ def read_file(self, file): thermal_expansion = float(item[0]) self.thermal_expansion = thermal_expansion - self.version_status = 'new' + version_status = 'new' elif 'VERSION' in inp[0]: jcpdsfile = open(file, 'r') @@ -228,7 +219,7 @@ def read_file(self, file): if jlinespl[0] == 'DALPHATDT:': dalphatdt = float(jlinespl[1]) - self.thermal_expansion_dt = dalphatat + self.thermal_expansion_dt = dalphatdt if jlinespl[0] == 'DIHKL:': pass @@ -266,20 +257,61 @@ def read_file(self, file): jcpdsfile.close() - self.version_status = 'new' + version_status = 'new' else: - self.version_status = 'old' - - if self.version_status == 'new': - self.a = self.a0 - self.b = self.b0 - self.c = self.c0 - self.alpha = self.alpha0 - self.beta = self.beta0 - self.gamma = self.gamma0 + version_status = 'old' + + if version_status == 'new': self.v = self.calc_volume_unitcell() + def verify_symmetry_match(self, mat): + if not self.symmetry_matches(mat): + msg = ( + f'The JCPDS symmetry "{self.symmetry}" does not match the ' + f'symmetry of the material "{mat.latticeType}"!' + ) + raise SymmetryMismatch(msg) + + def symmetry_matches(self, mat): + return mat.latticeType == self.symmetry + + @property + def _lat_param_names(self): + return ['a0', 'b0', 'c0', 'alpha0', 'beta0', 'gamma0'] + + @property + def lattice_params(self): + return [getattr(self, x) for x in self._lat_param_names] + + @lattice_params.setter + def lattice_params(self, v): + for name, val in zip(self._lat_param_names, v): + setattr(self, name, val) + + self.update_v0() + + def matches_material(self, mat): + self.verify_symmetry_match(mat) + mat_lp = [x.value for x in mat.latticeParameters] + for v1, v2 in zip(mat_lp, self.lattice_params): + if not np.isclose(v1, v2): + return False + + def load_from_material(self, mat): + if self.symmetry: + self.verify_symmetry_match(mat) + + lp = [x.value for x in mat.latticeParameters] + self.lattice_params = lp + + def write_to_material(self, mat): + self.verify_symmetry_match(mat) + mat.latticeParameters = self.lattice_params + + def update_v0(self): + self.v0 = self.calc_volume_unitcell() + def calc_volume_unitcell(self): '''calculate volume of the unitcell Ref: Introduction to conventional TEM @@ -333,7 +365,7 @@ def ktp(self, temperature=None): def calc_pressure(self, volume=None, temperature=None): '''calculate the pressure given the volume - and temperature using the third order + and temperature using the third order birch-murnaghan equation of state. ''' if volume is None: @@ -348,7 +380,7 @@ def calc_pressure(self, volume=None, temperature=None): def calc_volume(self, pressure=None, temperature=None): '''solve for volume in the birch-murnaghan EoS to - compute the volume. this number will be propagated + compute the volume. this number will be propagated to the Material object as updated lattice constants. ''' vt = self.vt(temperature=temperature) @@ -371,14 +403,14 @@ def calc_volume(self, pressure=None, temperature=None): mask = np.logical_and(res >= 0., res <= 1.0) res = res[mask] if len(res) != 1: - msg = (f'more than one physically ' - f'reasonable solution found!') + msg = ('more than one physically ' + 'reasonable solution found!') raise ValueError(msg) return res[0] * vt def calc_lp_factor(self, pressure=None, temperature=None): '''calculate the factor to multiply the lattice - constants by. only the lengths will be modified, the + constants by. only the lengths will be modified, the angles will be kept constant. ''' vt = self.vt(temperature=temperature) @@ -388,11 +420,15 @@ def calc_lp_factor(self, pressure=None, temperature=None): def calc_lp_at_PT(self, pressure=None, temperature=None): '''calculate the lattice parameters for a given - pressure and temperature using the BM EoS. This + pressure and temperature using the BM EoS. This is the main function which will be called from the GUI. ''' f = self.calc_lp_factor(pressure=pressure, temperature=temperature) return np.array([f*self.a0, f*self.b0, f*self.c0, - self.alpha0, self.beta0, self.gamma0,]) \ No newline at end of file + self.alpha0, self.beta0, self.gamma0,]) + + +class SymmetryMismatch(Exception): + pass From 43fd6b96d1c9879f3e7ec7e3324281de04cde2b9 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Tue, 7 Nov 2023 16:28:26 -0800 Subject: [PATCH 10/26] import path fix --- hexrd/wppf/phase.py | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/hexrd/wppf/phase.py b/hexrd/wppf/phase.py index ce6a878bd..b6bd5f58b 100644 --- a/hexrd/wppf/phase.py +++ b/hexrd/wppf/phase.py @@ -1,7 +1,8 @@ import numpy as np from hexrd.valunits import valWUnit from hexrd.material.spacegroup import Allowed_HKLs, SpaceGroup -from hexrd import symmetry, symbols, constants +from hexrd import constants +from hexrd.material import symmetry, symbols from hexrd.material import Material from hexrd.material.unitcell import _rqpDict from hexrd.wppf import wppfsupport From 106286d362e8e4d889ec5d61eca9869467775ce8 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Wed, 8 Nov 2023 11:39:59 -0800 Subject: [PATCH 11/26] Hande I/O of 'k0', 'k0p' and 'v0'. --- hexrd/material/material.py | 37 +++++++++++++++++++++++++++++++++++++ hexrd/material/mksupport.py | 8 ++++++++ 2 files changed, 45 insertions(+) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index 2313609d7..a15169e54 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -259,6 +259,9 @@ def _newUnitcell(self): if pdata.laueGroup != laue or pdata.lparms != reduced_lparms: pdata.set_laue_and_lparms(laue, reduced_lparms) + if not hasattr(self, 'v0'): + self.v0 = self.vol + def _hkls_changed(self): # Call this when something happens that changes the hkls... self._newPdata() @@ -704,6 +707,28 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): else: self.stiffness = numpy.zeros([6, 6]) + self.k0 = 100.0 + if('k0' in gid): + # this is the isotropic bulk modulus + k0 = numpy.array(gid.get('k0'), + dtype=numpy.float64).item() + self.k0 = k0 + + self.k0p = 0.0 + if('k0p' in gid): + # this is the pressure derivation of + # the isotropic bulk modulus + k0p = numpy.array(gid.get('k0p'), + dtype=numpy.float64).item() + self.k0p = k0p + + if('v0' in gid): + # this is the isotropic ambient unitcell + # volume + v0 = numpy.array(gid.get('v0'), + dtype=numpy.float64).item() + self.v0 = v0 + # if dmin is present in file: if('dmin' in gid): # if dmin is present in the HDF5 file, then use that @@ -761,6 +786,18 @@ def dump_material(self, file, path=None): AtomInfo['tThWidth'] = numpy.degrees(Material.DFLT_TTH) else: AtomInfo['tThWidth'] = numpy.degrees(self.planeData.tThWidth) + + AtomInfo['k0'] = 100.0 + if hasattr(self, 'k0'): + AtomInfo['k0'] = self.k0 + + AtomInfo['k0p'] = 0.0 + if hasattr(self, 'k0p'): + AtomInfo['k0p'] = self.k0p + + AtomInfo['v0'] = self.vol + if hasattr(self, 'v0'): + AtomInfo['v0'] = self.v0 ''' lattice parameters ''' diff --git a/hexrd/material/mksupport.py b/hexrd/material/mksupport.py index b575cf66b..f2f0a10cc 100644 --- a/hexrd/material/mksupport.py +++ b/hexrd/material/mksupport.py @@ -515,6 +515,14 @@ def WriteH5Data(fid, AtomInfo, lat_param, path=None): did = gid.create_dataset("stiffness", (6, 6), dtype=np.float64) did.write_direct(np.array(AtomInfo['stiffness'], dtype=np.float64)) + did = gid.create_dataset("k0", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['k0']], dtype=np.float64)) + + did = gid.create_dataset("k0p", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['k0p']], dtype=np.float64)) + + did = gid.create_dataset("v0", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['v0']], dtype=np.float64)) if "tThWidth" in AtomInfo: did = gid.create_dataset("tThWidth", (1,), dtype=np.float64) From 3ab929061de52008ed6ab3c8438dfc02d88fd690 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Wed, 8 Nov 2023 13:13:19 -0800 Subject: [PATCH 12/26] add temperature variables to the material class --- hexrd/material/material.py | 24 ++++++++++++++++++++++++ hexrd/material/mksupport.py | 6 ++++++ 2 files changed, 30 insertions(+) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index a15169e54..f15078ea9 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -722,6 +722,22 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): dtype=numpy.float64).item() self.k0p = k0p + self.dk0dt = 0.0 + if('dk0dt' in gid): + # this is the temperature derivation of + # the isotropic bulk modulus + dk0dt = numpy.array(gid.get('dk0dt'), + dtype=numpy.float64).item() + self.dk0dt = dk0dt + + self.dk0pdt = 0.0 + if('dk0pdt' in gid): + # this is the temperature derivation of + # the pressure derivative of isotropic bulk modulus + dk0pdt = numpy.array(gid.get('dk0pdt'), + dtype=numpy.float64).item() + self.dk0pdt = dk0pdt + if('v0' in gid): # this is the isotropic ambient unitcell # volume @@ -798,6 +814,14 @@ def dump_material(self, file, path=None): AtomInfo['v0'] = self.vol if hasattr(self, 'v0'): AtomInfo['v0'] = self.v0 + + AtomInfo['dk0dt'] = 0.0 + if hasattr(self, 'dk0dt'): + AtomInfo['dk0dt'] = self.dk0dt + + AtomInfo['dk0pdt'] = 0.0 + if hasattr(self, 'dk0pdt'): + AtomInfo['v0'] = self.dk0pdt ''' lattice parameters ''' diff --git a/hexrd/material/mksupport.py b/hexrd/material/mksupport.py index f2f0a10cc..3f6a30416 100644 --- a/hexrd/material/mksupport.py +++ b/hexrd/material/mksupport.py @@ -521,6 +521,12 @@ def WriteH5Data(fid, AtomInfo, lat_param, path=None): did = gid.create_dataset("k0p", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['k0p']], dtype=np.float64)) + did = gid.create_dataset("dk0dt", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['dk0dt']], dtype=np.float64)) + + did = gid.create_dataset("dk0pdt", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['dk0pdt']], dtype=np.float64)) + did = gid.create_dataset("v0", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['v0']], dtype=np.float64)) From 093b8e32675031e55ff464fea8404889949f7dc6 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Wed, 8 Nov 2023 14:26:46 -0800 Subject: [PATCH 13/26] I/O for material class with JCPDS parameters. Moved BM EOS functions to material class --- hexrd/material/jcpds.py | 103 ---------------------- hexrd/material/material.py | 171 +++++++++++++++++++++++++++++++++++- hexrd/material/mksupport.py | 6 ++ 3 files changed, 176 insertions(+), 104 deletions(-) diff --git a/hexrd/material/jcpds.py b/hexrd/material/jcpds.py index 0b515dccb..39b21a868 100644 --- a/hexrd/material/jcpds.py +++ b/hexrd/material/jcpds.py @@ -326,109 +326,6 @@ def calc_volume_unitcell(self): +2*ca*cb*cg) return v0*f - def vt(self, temperature=None): - '''calculate volume at high - temperature - ''' - alpha0 = self.thermal_expansion - alpha1 = self.thermal_expansion_dt - if temperature is None: - vt = self.v0 - else: - delT = (temperature-298) - delT2 = (temperature**2-298**2) - vt = self.v0*np.exp(alpha0*delT - +0.5*alpha1*delT2) - return vt - - def kt(self, temperature=None): - '''calculate bulk modulus for - high temperature - ''' - k0 = self.k0 - if temperature is None: - return k0 - else: - delT = (temperature-298) - return k0 + self.dk0dt*delT - - def ktp(self, temperature=None): - '''calculate bulk modulus derivative - for high temperature - ''' - k0p = self.k0p - if temperature is None: - return k0p - else: - delT = (temperature-298) - return k0p + self.dk0pdt*delT - - def calc_pressure(self, volume=None, temperature=None): - '''calculate the pressure given the volume - and temperature using the third order - birch-murnaghan equation of state. - ''' - if volume is None: - return 0 - else: - vt = self.vt(temperature=temperature) - kt = self.kt(temperature=temperature) - ktp = self.ktp(temperature=temperature) - f = 0.5*((vt/volume)**(2./3.) - 1) - - return 3.0*kt*f*(1 - 1.5*(4 - ktp)*f)*(1 + 2*f)**2.5 - - def calc_volume(self, pressure=None, temperature=None): - '''solve for volume in the birch-murnaghan EoS to - compute the volume. this number will be propagated - to the Material object as updated lattice constants. - ''' - vt = self.vt(temperature=temperature) - kt = self.kt(temperature=temperature) - ktp = self.ktp(temperature=temperature) - - if pressure is None: - return vt - else: - alpha = 0.75*(ktp - 4) - p = np.zeros([10,]) - p[0] = 1.0 - p[2] = (1 - 2*alpha)/alpha - p[4] = (alpha-1)/alpha - p[9] = -2*pressure/3/kt/alpha - res = np.roots(p) - res = res[np.isreal(res)] - res = 1/np.real(res)**3 - - mask = np.logical_and(res >= 0., res <= 1.0) - res = res[mask] - if len(res) != 1: - msg = ('more than one physically ' - 'reasonable solution found!') - raise ValueError(msg) - return res[0] * vt - - def calc_lp_factor(self, pressure=None, temperature=None): - '''calculate the factor to multiply the lattice - constants by. only the lengths will be modified, the - angles will be kept constant. - ''' - vt = self.vt(temperature=temperature) - vpt = self.calc_volume(pressure=pressure, - temperature=temperature) - return (vpt/vt)**(1./3.) - - def calc_lp_at_PT(self, pressure=None, temperature=None): - '''calculate the lattice parameters for a given - pressure and temperature using the BM EoS. This - is the main function which will be called from - the GUI. - ''' - f = self.calc_lp_factor(pressure=pressure, - temperature=temperature) - return np.array([f*self.a0, f*self.b0, f*self.c0, - self.alpha0, self.beta0, self.gamma0,]) - class SymmetryMismatch(Exception): pass diff --git a/hexrd/material/material.py b/hexrd/material/material.py index f15078ea9..8e62bcba2 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -193,6 +193,13 @@ def __init__(self, name=None, self._dmin = Material.DFLT_DMIN self._beamEnergy = Material.DFLT_KEV + self.k0 = 100.0 + self.k0p = 0.0 + self.dk0dt = 0.0 + self.dk0pdt = 0.0 + self.alphaT = 0.0 + self.dalphaTdT = 0.0 + # If these were specified, they override any other method of # obtaining them (including loading them from files). if dmin is not None: @@ -405,6 +412,110 @@ def remove_duplicate_atoms(self): self.charge = self.unitcell.chargestates self._hkls_changed() + + def vt(self, temperature=None): + '''calculate volume at high + temperature + ''' + alpha0 = self.thermal_expansion + alpha1 = self.thermal_expansion_dt + if temperature is None: + vt = self.v0 + else: + delT = (temperature-298) + delT2 = (temperature**2-298**2) + vt = self.v0*np.exp(alpha0*delT + +0.5*alpha1*delT2) + return vt + + def kt(self, temperature=None): + '''calculate bulk modulus for + high temperature + ''' + k0 = self.k0 + if temperature is None: + return k0 + else: + delT = (temperature-298) + return k0 + self.dk0dt*delT + + def ktp(self, temperature=None): + '''calculate bulk modulus derivative + for high temperature + ''' + k0p = self.k0p + if temperature is None: + return k0p + else: + delT = (temperature-298) + return k0p + self.dk0pdt*delT + + def calc_pressure(self, volume=None, temperature=None): + '''calculate the pressure given the volume + and temperature using the third order + birch-murnaghan equation of state. + ''' + if volume is None: + return 0 + else: + vt = self.vt(temperature=temperature) + kt = self.kt(temperature=temperature) + ktp = self.ktp(temperature=temperature) + f = 0.5*((vt/volume)**(2./3.) - 1) + + return 3.0*kt*f*(1 - 1.5*(4 - ktp)*f)*(1 + 2*f)**2.5 + + def calc_volume(self, pressure=None, temperature=None): + '''solve for volume in the birch-murnaghan EoS to + compute the volume. this number will be propagated + to the Material object as updated lattice constants. + ''' + vt = self.vt(temperature=temperature) + kt = self.kt(temperature=temperature) + ktp = self.ktp(temperature=temperature) + + if pressure is None: + return vt + else: + alpha = 0.75*(ktp - 4) + p = np.zeros([10,]) + p[0] = 1.0 + p[2] = (1 - 2*alpha)/alpha + p[4] = (alpha-1)/alpha + p[9] = -2*pressure/3/kt/alpha + res = np.roots(p) + res = res[np.isreal(res)] + res = 1/np.real(res)**3 + + mask = np.logical_and(res >= 0., res <= 1.0) + res = res[mask] + if len(res) != 1: + msg = ('more than one physically ' + 'reasonable solution found!') + raise ValueError(msg) + return res[0] * vt + + def calc_lp_factor(self, pressure=None, temperature=None): + '''calculate the factor to multiply the lattice + constants by. only the lengths will be modified, the + angles will be kept constant. + ''' + vt = self.vt(temperature=temperature) + vpt = self.calc_volume(pressure=pressure, + temperature=temperature) + return (vpt/vt)**(1./3.) + + def calc_lp_at_PT(self, pressure=None, temperature=None): + '''calculate the lattice parameters for a given + pressure and temperature using the BM EoS. This + is the main function which will be called from + the GUI. + ''' + f = self.calc_lp_factor(pressure=pressure, + temperature=temperature) + return np.array([f*self.a0, f*self.b0, f*self.c0, + self.alpha0, self.beta0, self.gamma0,]) + def _readCif(self, fcif=DFLT_NAME+'.cif'): """ @@ -627,6 +738,18 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): self._beamEnergy = Material.DFLT_KEV self._tThWidth = Material.DFLT_TTH + '''set the Birch-Murnaghan equation of state + parameters to default values. These values can + be updated by user or by reading a JCPDS file + ''' + self.k0 = 100.0 + self.k0p = 0.0 + self.dk0dt = 0.0 + self.dk0pdt = 0.0 + self.alphaT = 0.0 + self.dalphaTdT = 0.0 + + def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): """ >> @AUTHOR: Saransh Singh, Lawrence Livermore National Lab, @@ -707,6 +830,9 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): else: self.stiffness = numpy.zeros([6, 6]) + '''start reading the Birch-Murnaghan equation of state + parameters + ''' self.k0 = 100.0 if('k0' in gid): # this is the isotropic bulk modulus @@ -738,6 +864,25 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): dtype=numpy.float64).item() self.dk0pdt = dk0pdt + self.alphaT = 0.0 + if('alphaT' in gid): + # this is the temperature derivation of + # the pressure derivative of isotropic bulk modulus + alphaT = numpy.array(gid.get('alphaT'), + dtype=numpy.float64).item() + self.alphaT = alphaT + + self.dalphaTdT = 0.0 + if('dalphaTdT' in gid): + # this is the temperature derivation of + # the pressure derivative of isotropic bulk modulus + dalphaTdT = numpy.array(gid.get('dalphaTdT'), + dtype=numpy.float64).item() + self.dalphaTdT = dalphaTdT + + '''Finished with the BM EOS + ''' + if('v0' in gid): # this is the isotropic ambient unitcell # volume @@ -821,7 +966,15 @@ def dump_material(self, file, path=None): AtomInfo['dk0pdt'] = 0.0 if hasattr(self, 'dk0pdt'): - AtomInfo['v0'] = self.dk0pdt + AtomInfo['dk0pdt'] = self.dk0pdt + + AtomInfo['alphaT'] = 0.0 + if hasattr(self, 'alphaT'): + AtomInfo['alphaT'] = self.alphaT + + AtomInfo['dalphaTdT'] = 0.0 + if hasattr(self, 'dalphaTdT'): + AtomInfo['dalphaTdT'] = self.dalphaTdT ''' lattice parameters ''' @@ -1081,6 +1234,22 @@ def _set_atomtype(self, v): _get_atomtype, _set_atomtype, None, "Information about atomic types") + @property + def thermal_expansion(self): + return self.alphaT + + @thermal_expansion.setter + def thermal_expansion(self, val): + self.alphaT = val + + @property + def thermal_expansion_dt(self): + return self.dalphaTdT + + @thermal_expansion_dt.setter + def thermal_expansion_dt(self, val): + self.dalphaTdT = val + def _set_atomdata(self, atomtype, atominfo, U, charge): """ sometimes the number of atom types and their diff --git a/hexrd/material/mksupport.py b/hexrd/material/mksupport.py index 3f6a30416..c3b44c2de 100644 --- a/hexrd/material/mksupport.py +++ b/hexrd/material/mksupport.py @@ -527,6 +527,12 @@ def WriteH5Data(fid, AtomInfo, lat_param, path=None): did = gid.create_dataset("dk0pdt", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['dk0pdt']], dtype=np.float64)) + did = gid.create_dataset("alphaT", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['alphaT']], dtype=np.float64)) + + did = gid.create_dataset("dalphaTdT", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['dalphaTdT']], dtype=np.float64)) + did = gid.create_dataset("v0", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['v0']], dtype=np.float64)) From c812f78694f7a2566f1b994aa5e01c1bf59f5007 Mon Sep 17 00:00:00 2001 From: Saransh Singh Date: Wed, 8 Nov 2023 14:37:00 -0800 Subject: [PATCH 14/26] change np --> numpy --- hexrd/material/material.py | 14 +++++++------- 1 file changed, 7 insertions(+), 7 deletions(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index 8e62bcba2..426ab81a3 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -424,7 +424,7 @@ def vt(self, temperature=None): else: delT = (temperature-298) delT2 = (temperature**2-298**2) - vt = self.v0*np.exp(alpha0*delT + vt = self.v0*numpy.exp(alpha0*delT +0.5*alpha1*delT2) return vt @@ -478,16 +478,16 @@ def calc_volume(self, pressure=None, temperature=None): return vt else: alpha = 0.75*(ktp - 4) - p = np.zeros([10,]) + p = numpy.zeros([10,]) p[0] = 1.0 p[2] = (1 - 2*alpha)/alpha p[4] = (alpha-1)/alpha p[9] = -2*pressure/3/kt/alpha - res = np.roots(p) - res = res[np.isreal(res)] - res = 1/np.real(res)**3 + res = numpy.roots(p) + res = res[numpy.isreal(res)] + res = 1/numpy.real(res)**3 - mask = np.logical_and(res >= 0., res <= 1.0) + mask = numpy.logical_and(res >= 0., res <= 1.0) res = res[mask] if len(res) != 1: msg = ('more than one physically ' @@ -513,7 +513,7 @@ def calc_lp_at_PT(self, pressure=None, temperature=None): ''' f = self.calc_lp_factor(pressure=pressure, temperature=temperature) - return np.array([f*self.a0, f*self.b0, f*self.c0, + return numpy.array([f*self.a0, f*self.b0, f*self.c0, self.alpha0, self.beta0, self.gamma0,]) def _readCif(self, fcif=DFLT_NAME+'.cif'): From e31bb162aa92d3abe92bf77760e5be328ab71f1c Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Fri, 10 Nov 2023 12:02:06 -0600 Subject: [PATCH 15/26] Memoize the result of GenerateSGSym This is a time consuming function, its input is simple, and its output doesn't take up too much memory. It's the perfect candidate for memoization. Signed-off-by: Patrick Avery --- hexrd/material/symmetry.py | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/hexrd/material/symmetry.py b/hexrd/material/symmetry.py index 27e420083..a6a566eb4 100644 --- a/hexrd/material/symmetry.py +++ b/hexrd/material/symmetry.py @@ -37,7 +37,7 @@ # from hexrd.rotations import quatOfAngleAxis, quatProductMatrix, fixQuat from hexrd import rotations as rot from hexrd import constants -from hexrd.utils.decorators import numba_njit_if_available +from hexrd.utils.decorators import memoize, numba_njit_if_available # ============================================================================= @@ -417,6 +417,8 @@ def SYM_fillgen(t): mat = np.broadcast_to(mat, [1,4,4]) return mat + +@memoize(maxsize=20) def GenerateSGSym(sgnum, setting=0): ''' From 9d299141b1cc287916a5335015abd8da112fbb88 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 30 Nov 2023 09:19:03 -0600 Subject: [PATCH 16/26] Change all "numpy" to "np" in material.py "np" is the standard, so let's use that. Signed-off-by: Patrick Avery --- hexrd/material/material.py | 180 ++++++++++++++++++------------------- 1 file changed, 90 insertions(+), 90 deletions(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index 426ab81a3..f7b247798 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -32,7 +32,7 @@ materials are defined by name in materialDict. """ from configparser import SafeConfigParser as Parser -import numpy +import numpy as np from hexrd.material.crystallography import PlaneData as PData from hexrd.material import symmetry, unitcell @@ -98,7 +98,7 @@ class Material(object): DFLT_KEV = valWUnit('wavelength', 'energy', 80.725e0, 'keV') DFLT_STR = 0.0025 - DFLT_TTH = numpy.radians(0.25) + DFLT_TTH = np.radians(0.25) DFLT_TTHMAX = None """ ATOMINFO Fractional Atom Position of an atom in the @@ -108,10 +108,10 @@ class Material(object): ATOMTYPE atomic number of all the different species in the unitcell """ - DFLT_ATOMINFO = numpy.array([[0., 0., 0., 1.]]) - DFLT_U = numpy.array([6.33E-3]) - DFLT_ATOMTYPE = numpy.array([28]) - DFLT_CHARGE = numpy.array(["0"]) + DFLT_ATOMINFO = np.array([[0., 0., 0., 1.]]) + DFLT_U = np.array([6.33E-3]) + DFLT_ATOMTYPE = np.array([28]) + DFLT_CHARGE = np.array(["0"]) ''' the dmin parameter is used to figure out the maximum sampling for g-vectors @@ -122,7 +122,7 @@ class Material(object): ''' default stiffness tensor in voight notation ''' - DFLT_STIFFNESS = numpy.eye(6) + DFLT_STIFFNESS = np.eye(6) ''' some materials have more than one space group setting. for ex @@ -298,10 +298,10 @@ def _newPdata(self): old_indices, new_indices = map_hkls_to_new(old_pdata, self._pData) # Map the new exclusions to the old ones. New ones default to True. - exclusions = numpy.ones(hkls.shape[1], dtype=bool) + exclusions = np.ones(hkls.shape[1], dtype=bool) exclusions[new_indices] = old_pdata.exclusions[old_indices] - if numpy.all(exclusions): + if np.all(exclusions): # If they are all excluded, just set the default exclusions self.set_default_exclusions() else: @@ -328,7 +328,7 @@ def set_default_exclusions(self): def enable_hkls_from_file(self): # Enable hkls from the file # 'hkl_from_file' must be an attribute on `self` - exclusions = numpy.ones_like(self._pData.exclusions, dtype=bool) + exclusions = np.ones_like(self._pData.exclusions, dtype=bool) for i, g in enumerate(self._pData.hklDataList): if g['hkl'].tolist() in self.hkl_from_file.tolist(): exclusions[i] = False @@ -337,7 +337,7 @@ def enable_hkls_from_file(self): def enable_hkls_below_index(self, index=5): # Enable hkls with indices less than @index - exclusions = numpy.ones_like(self._pData.exclusions, dtype=bool) + exclusions = np.ones_like(self._pData.exclusions, dtype=bool) exclusions[:index] = False self._pData.exclusions = exclusions @@ -345,12 +345,12 @@ def enable_hkls_below_tth(self, tth_threshold=90.0): ''' enable reflections with two-theta less than @tth_threshold degrees ''' - tth_threshold = numpy.radians(tth_threshold) + tth_threshold = np.radians(tth_threshold) - tth = numpy.array([hkldata['tTheta'] + tth = np.array([hkldata['tTheta'] for hkldata in self._pData.hklDataList]) - dflt_excl = numpy.ones(tth.shape, dtype=numpy.bool) - dflt_excl2 = numpy.ones(tth.shape, dtype=numpy.bool) + dflt_excl = np.ones(tth.shape, dtype=np.bool) + dflt_excl2 = np.ones(tth.shape, dtype=np.bool) if(hasattr(self, 'hkl_from_file')): """ @@ -361,16 +361,16 @@ def enable_hkls_below_tth(self, tth_threshold=90.0): if(g['hkl'].tolist() in self.hkl_from_file.tolist()): dflt_excl[i] = False - dflt_excl2[~numpy.isnan(tth)] = \ - ~((tth[~numpy.isnan(tth)] >= 0.0) & - (tth[~numpy.isnan(tth)] <= tth_threshold)) + dflt_excl2[~np.isnan(tth)] = \ + ~((tth[~np.isnan(tth)] >= 0.0) & + (tth[~np.isnan(tth)] <= tth_threshold)) - dflt_excl = numpy.logical_or(dflt_excl, dflt_excl2) + dflt_excl = np.logical_or(dflt_excl, dflt_excl2) else: - dflt_excl[~numpy.isnan(tth)] = \ - ~((tth[~numpy.isnan(tth)] >= 0.0) & - (tth[~numpy.isnan(tth)] <= tth_threshold)) + dflt_excl[~np.isnan(tth)] = \ + ~((tth[~np.isnan(tth)] >= 0.0) & + (tth[~np.isnan(tth)] <= tth_threshold)) dflt_excl[0] = False self._pData.exclusions = dflt_excl @@ -381,7 +381,7 @@ def update_structure_factor(self): self.planeData.set_structFact(sf) def compute_powder_overlay(self, - ttharray=numpy.linspace(0, 80, 2000), + ttharray=np.linspace(0, 80, 2000), fwhm=0.25, scale=1.0): """ @@ -393,9 +393,9 @@ def compute_powder_overlay(self, requested feature from Amy Jenei """ - tth = numpy.degrees(self.planeData.getTTh()) # convert to degrees + tth = np.degrees(self.planeData.getTTh()) # convert to degrees Ip = self.planeData.powder_intensity - self.powder_overlay = numpy.zeros_like(ttharray) + self.powder_overlay = np.zeros_like(ttharray) for t, I in zip(tth, Ip): p = [t, fwhm] self.powder_overlay += scale*I*_unit_gaussian(p, ttharray) @@ -424,7 +424,7 @@ def vt(self, temperature=None): else: delT = (temperature-298) delT2 = (temperature**2-298**2) - vt = self.v0*numpy.exp(alpha0*delT + vt = self.v0*np.exp(alpha0*delT +0.5*alpha1*delT2) return vt @@ -478,16 +478,16 @@ def calc_volume(self, pressure=None, temperature=None): return vt else: alpha = 0.75*(ktp - 4) - p = numpy.zeros([10,]) + p = np.zeros([10,]) p[0] = 1.0 p[2] = (1 - 2*alpha)/alpha p[4] = (alpha-1)/alpha p[9] = -2*pressure/3/kt/alpha - res = numpy.roots(p) - res = res[numpy.isreal(res)] - res = 1/numpy.real(res)**3 + res = np.roots(p) + res = res[np.isreal(res)] + res = 1/np.real(res)**3 - mask = numpy.logical_and(res >= 0., res <= 1.0) + mask = np.logical_and(res >= 0., res <= 1.0) res = res[mask] if len(res) != 1: msg = ('more than one physically ' @@ -513,7 +513,7 @@ def calc_lp_at_PT(self, pressure=None, temperature=None): ''' f = self.calc_lp_factor(pressure=pressure, temperature=temperature) - return numpy.array([f*self.a0, f*self.b0, f*self.c0, + return np.array([f*self.a0, f*self.b0, f*self.c0, self.alpha0, self.beta0, self.gamma0,]) def _readCif(self, fcif=DFLT_NAME+'.cif'): @@ -629,8 +629,8 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): sometimes cif files have negative values so need to bring them back to fractional coordinates between 0-1 ''' - pos = numpy.asarray(pos).astype(numpy.float64) - pos, _ = numpy.modf(pos+100.0) + pos = np.asarray(pos).astype(np.float64) + pos, _ = np.modf(pos+100.0) atompos.append(pos) """note that the vibration amplitude, U is just the amplitude (in A) @@ -645,7 +645,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): if(not pocc): warn('occupation fraction not present. setting it to 1') - occ = numpy.ones(atompos[0].shape) + occ = np.ones(atompos[0].shape) atompos.append(occ) else: slist = cifdata[occ_U[0]] @@ -658,17 +658,17 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): else: occ.append(p) - chkstr = numpy.asarray([isinstance(x, str) for x in occ]) - occstr = numpy.array(occ) + chkstr = np.asarray([isinstance(x, str) for x in occ]) + occstr = np.array(occ) occstr[chkstr] = 1.0 - atompos.append(numpy.asarray(occstr).astype(numpy.float64)) + atompos.append(np.asarray(occstr).astype(np.float64)) if(not pU): msg = (f"'Debye-Waller factors not present. " f"setting to same values for all atoms.'") warn(msg) - U = self.DFLT_U[0] * numpy.ones(atompos[0].shape) + U = self.DFLT_U[0] * np.ones(atompos[0].shape) self._U = U else: if(occ_U[1] in cifdata.keys()): @@ -686,7 +686,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): else: U.append(p) - chkstr = numpy.asarray([isinstance(x, str) for x in U]) + chkstr = np.asarray([isinstance(x, str) for x in U]) for ii, x in enumerate(chkstr): if x: @@ -695,11 +695,11 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): except ValueError: U[ii] = self.DFLT_U[0] - self._U = numpy.asarray(U) + self._U = np.asarray(U) ''' format everything in the right shape etc. ''' - self._atominfo = numpy.asarray(atompos).T + self._atominfo = np.asarray(atompos).T ''' get atome types here i.e. the atomic number of atoms at each site @@ -730,7 +730,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): atomtype.append(ptable[ss]) charge.append(c) - self._atomtype = numpy.asarray(atomtype).astype(numpy.int32) + self._atomtype = np.asarray(atomtype).astype(np.int32) self._charge = charge self._sgsetting = 0 @@ -781,8 +781,8 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): gid = root_gid.get(xtal) - sgnum = numpy.array(gid.get('SpaceGroupNumber'), - dtype=numpy.int32).item() + sgnum = np.array(gid.get('SpaceGroupNumber'), + dtype=np.int32).item() """ IMPORTANT NOTE: note that the latice parameters is nm by default @@ -804,31 +804,31 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): self.sgnum = sgnum # the U factors are related to B by the relation B = 8pi^2 U - self._atominfo = numpy.transpose(numpy.array( - gid.get('AtomData'), dtype=numpy.float64)) - self._U = numpy.transpose(numpy.array( - gid.get('U'), dtype=numpy.float64)) + self._atominfo = np.transpose(np.array( + gid.get('AtomData'), dtype=np.float64)) + self._U = np.transpose(np.array( + gid.get('U'), dtype=np.float64)) # read atom types (by atomic number, Z) - self._atomtype = numpy.array(gid.get('Atomtypes'), dtype=numpy.int32) + self._atomtype = np.array(gid.get('Atomtypes'), dtype=np.int32) if 'ChargeStates' in gid: self._charge = h5py_read_string(gid['ChargeStates']) else: self._charge = ['0']*self._atomtype.shape[0] self._atom_ntype = self._atomtype.shape[0] - self._sgsetting = numpy.array(gid.get('SpaceGroupSetting'), - dtype=numpy.int32).item() + self._sgsetting = np.array(gid.get('SpaceGroupSetting'), + dtype=np.int32).item() if('stiffness' in gid): # we're assuming the stiffness is in units of GPa - self.stiffness = numpy.array(gid.get('stiffness')) + self.stiffness = np.array(gid.get('stiffness')) elif('compliance' in gid): # we're assuming the compliance is in units of TPa^-1 - self.stiffness = numpy.linalg.inv( - numpy.array(gid.get('compliance'))) * 1.e3 + self.stiffness = np.linalg.inv( + np.array(gid.get('compliance'))) * 1.e3 else: - self.stiffness = numpy.zeros([6, 6]) + self.stiffness = np.zeros([6, 6]) '''start reading the Birch-Murnaghan equation of state parameters @@ -836,48 +836,48 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): self.k0 = 100.0 if('k0' in gid): # this is the isotropic bulk modulus - k0 = numpy.array(gid.get('k0'), - dtype=numpy.float64).item() + k0 = np.array(gid.get('k0'), + dtype=np.float64).item() self.k0 = k0 self.k0p = 0.0 if('k0p' in gid): # this is the pressure derivation of # the isotropic bulk modulus - k0p = numpy.array(gid.get('k0p'), - dtype=numpy.float64).item() + k0p = np.array(gid.get('k0p'), + dtype=np.float64).item() self.k0p = k0p self.dk0dt = 0.0 if('dk0dt' in gid): # this is the temperature derivation of # the isotropic bulk modulus - dk0dt = numpy.array(gid.get('dk0dt'), - dtype=numpy.float64).item() + dk0dt = np.array(gid.get('dk0dt'), + dtype=np.float64).item() self.dk0dt = dk0dt self.dk0pdt = 0.0 if('dk0pdt' in gid): # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - dk0pdt = numpy.array(gid.get('dk0pdt'), - dtype=numpy.float64).item() + dk0pdt = np.array(gid.get('dk0pdt'), + dtype=np.float64).item() self.dk0pdt = dk0pdt self.alphaT = 0.0 if('alphaT' in gid): # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - alphaT = numpy.array(gid.get('alphaT'), - dtype=numpy.float64).item() + alphaT = np.array(gid.get('alphaT'), + dtype=np.float64).item() self.alphaT = alphaT self.dalphaTdT = 0.0 if('dalphaTdT' in gid): # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - dalphaTdT = numpy.array(gid.get('dalphaTdT'), - dtype=numpy.float64).item() + dalphaTdT = np.array(gid.get('dalphaTdT'), + dtype=np.float64).item() self.dalphaTdT = dalphaTdT '''Finished with the BM EOS @@ -886,39 +886,39 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): if('v0' in gid): # this is the isotropic ambient unitcell # volume - v0 = numpy.array(gid.get('v0'), - dtype=numpy.float64).item() + v0 = np.array(gid.get('v0'), + dtype=np.float64).item() self.v0 = v0 # if dmin is present in file: if('dmin' in gid): # if dmin is present in the HDF5 file, then use that - dmin = numpy.array(gid.get('dmin'), - dtype=numpy.float64).item() + dmin = np.array(gid.get('dmin'), + dtype=np.float64).item() self._dmin = _angstroms(dmin*10.) else: self._dmin = Material.DFLT_DMIN # if kev is present in file if('kev' in gid): - kev = numpy.array(gid.get('kev'), - dtype=numpy.float64).item() + kev = np.array(gid.get('kev'), + dtype=np.float64).item() self._beamEnergy = _kev(kev) else: self._beamEnergy = Material.DFLT_KEV if 'tThWidth' in gid: - tThWidth = numpy.array(gid.get('tThWidth'), - dtype=numpy.float64).item() - tThWidth = numpy.radians(tThWidth) + tThWidth = np.array(gid.get('tThWidth'), + dtype=np.float64).item() + tThWidth = np.radians(tThWidth) else: tThWidth = Material.DFLT_TTH self._tThWidth = tThWidth if('hkls' in gid): - self.hkl_from_file = numpy.array(gid.get('hkls'), - dtype=numpy.int32) + self.hkl_from_file = np.array(gid.get('hkls'), + dtype=np.int32) if isinstance(fhdf, (Path, str)): # Close the file... @@ -944,9 +944,9 @@ def dump_material(self, file, path=None): AtomInfo['dmin'] = self.unitcell.dmin AtomInfo['kev'] = self.beamEnergy.getVal("keV") if self.planeData.tThWidth is None: - AtomInfo['tThWidth'] = numpy.degrees(Material.DFLT_TTH) + AtomInfo['tThWidth'] = np.degrees(Material.DFLT_TTH) else: - AtomInfo['tThWidth'] = numpy.degrees(self.planeData.tThWidth) + AtomInfo['tThWidth'] = np.degrees(self.planeData.tThWidth) AtomInfo['k0'] = 100.0 if hasattr(self, 'k0'): @@ -1014,7 +1014,7 @@ def vol(self): @property def lparms(self): - return numpy.array([x.getVal("nm") if x.isLength() else + return np.array([x.getVal("nm") if x.isLength() else x.getVal("degrees") for x in self._lparms]) @property @@ -1043,7 +1043,7 @@ def U(self): @U.setter def U(self, Uarr): - Uarr = numpy.array(Uarr) + Uarr = np.array(Uarr) self._U = Uarr # property: sgnum @@ -1245,7 +1245,7 @@ def thermal_expansion(self, val): @property def thermal_expansion_dt(self): return self.dalphaTdT - + @thermal_expansion_dt.setter def thermal_expansion_dt(self, val): self.dalphaTdT = val @@ -1269,9 +1269,9 @@ def _set_atomdata(self, atomtype, atominfo, U, charge): """ # check for consistency of sizes here - atomtype = numpy.array(atomtype) - atominfo = numpy.array(atominfo) - U = numpy.array(U) + atomtype = np.array(atomtype) + atominfo = np.array(atominfo) + U = np.array(U) if atomtype.shape[0] != atominfo.shape[0]: msg = (f"inconsistent shapes: number of atoms " @@ -1377,8 +1377,8 @@ def save_materials_hdf5(f, materials, path=None): def hkls_match(a, b): # Check if hkls match. Expects inputs to have shape (x, 3). def sorted_hkls(x): - return x[numpy.lexsort((x[:, 2], x[:, 1], x[:, 0]))] - return numpy.array_equal(sorted_hkls(a), sorted_hkls(b)) + return x[np.lexsort((x[:, 2], x[:, 1], x[:, 0]))] + return np.array_equal(sorted_hkls(a), sorted_hkls(b)) def map_hkls_to_new(old_pdata, new_pdata): @@ -1392,7 +1392,7 @@ def get_hkl_strings(pdata): 'ar2': get_hkl_strings(new_pdata), 'return_indices': True, } - return numpy.intersect1d(**kwargs)[1:] + return np.intersect1d(**kwargs)[1:] # # ============================== Executable section for testing From 716ffaee99e798167ff7306707d593e80e8df269 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 30 Nov 2023 09:27:15 -0600 Subject: [PATCH 17/26] Add config for Black This is so that we can run black on files without specifying the options manually. Signed-off-by: Patrick Avery --- pyproject.toml | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/pyproject.toml b/pyproject.toml index f08833fec..eb28bae5d 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,2 +1,6 @@ [build-system] requires = ["setuptools", "wheel", "numpy<1.27", "setuptools_scm[toml]", "pybind11>=2.11.0"] + +[tool.black] +line-length = 79 +skip-string-normalization = true From 0571c59a0bb34e83a3dc69671e6fe771bf2c6c88 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 30 Nov 2023 09:27:37 -0600 Subject: [PATCH 18/26] Run "black" on material.py Signed-off-by: Patrick Avery --- hexrd/material/material.py | 475 +++++++++++++++++++++---------------- 1 file changed, 265 insertions(+), 210 deletions(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index f7b247798..9abfbe520 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -38,7 +38,6 @@ from hexrd.material import symmetry, unitcell from hexrd.valunits import valWUnit from hexrd.constants import ptable, ptableinverse, chargestate -import copy from os import path from pathlib import Path @@ -47,8 +46,10 @@ from warnings import warn from hexrd.material.mksupport import Write2H5File from hexrd.material.symbols import ( - xtal_sys_dict,Hall_to_sgnum, - HM_to_sgnum) + xtal_sys_dict, + Hall_to_sgnum, + HM_to_sgnum, +) from hexrd.utils.compatibility import h5py_read_string from hexrd.fitting.peakfunctions import _unit_gaussian @@ -61,15 +62,15 @@ def _angstroms(x): - return valWUnit('lp', 'length', x, 'angstrom') + return valWUnit('lp', 'length', x, 'angstrom') def _degrees(x): - return valWUnit('lp', 'angle', x, 'degrees') + return valWUnit('lp', 'angle', x, 'degrees') def _kev(x): - return valWUnit('xrayenergy', 'energy', x, 'keV') + return valWUnit('xrayenergy', 'energy', x, 'keV') def _key(x): @@ -88,12 +89,19 @@ class Material(object): space group data. default data is for nickel, but name is material """ + DFLT_NAME = 'material.xtal' DFLT_XTAL = 'Ni' DFLT_SGNUM = 225 - DFLT_LPARMS = [_angstroms(3.61), _angstroms(3.61), _angstroms(3.61), - _degrees(90.0), _degrees(90.0), _degrees(90.0)] + DFLT_LPARMS = [ + _angstroms(3.61), + _angstroms(3.61), + _angstroms(3.61), + _degrees(90.0), + _degrees(90.0), + _degrees(90.0), + ] DFLT_SSMAX = 100 DFLT_KEV = valWUnit('wavelength', 'energy', 80.725e0, 'keV') @@ -108,8 +116,8 @@ class Material(object): ATOMTYPE atomic number of all the different species in the unitcell """ - DFLT_ATOMINFO = np.array([[0., 0., 0., 1.]]) - DFLT_U = np.array([6.33E-3]) + DFLT_ATOMINFO = np.array([[0.0, 0.0, 0.0, 1.0]]) + DFLT_U = np.array([6.33e-3]) DFLT_ATOMTYPE = np.array([28]) DFLT_CHARGE = np.array(["0"]) @@ -134,11 +142,14 @@ class Material(object): ''' DFLT_SGSETTING = 0 - def __init__(self, name=None, - material_file=None, - dmin=None, - kev=None, - sgsetting=DFLT_SGSETTING): + def __init__( + self, + name=None, + material_file=None, + dmin=None, + kev=None, + sgsetting=DFLT_SGSETTING, + ): """Constructor for Material name -- (str) name of crystal @@ -157,7 +168,6 @@ def __init__(self, name=None, self.sgsetting = sgsetting if material_file: - # >> @ date 08/20/2020 SS removing dependence on hklmax if isinstance(material_file, (Path, str)): form = Path(material_file).suffix[1:] @@ -234,7 +244,7 @@ def _reset_lparms(self): lparms = unitcell._rqpDict[ltype][1](lparms) lparms_vu = [] for i in range(6): - if(i < 3): + if i < 3: lparms_vu.append(_angstroms(lparms[i])) else: lparms_vu.append(_degrees(lparms[i])) @@ -250,10 +260,16 @@ def _newUnitcell(self): """ self._reset_lparms() self._unitcell = unitcell.unitcell( - self._lparms, self.sgnum, self._atomtype, self._charge, - self._atominfo, self._U, - self._dmin.getVal('nm'), self._beamEnergy.value, - self._sgsetting) + self._lparms, + self.sgnum, + self._atomtype, + self._charge, + self._atominfo, + self._U, + self._dmin.getVal('nm'), + self._beamEnergy.value, + self._sgsetting, + ) if hasattr(self, 'stiffness'): self._unitcell.stiffness = self.stiffness @@ -310,10 +326,15 @@ def _newPdata(self): # Make the PlaneData object from scratch... lprm = self.reduced_lattice_parameters laue = self.unitcell._laueGroup - self._pData = PData(hkls, lprm, laue, - self._beamEnergy, Material.DFLT_STR, - tThWidth=self._tThWidth, - tThMax=Material.DFLT_TTHMAX) + self._pData = PData( + hkls, + lprm, + laue, + self._beamEnergy, + Material.DFLT_STR, + tThWidth=self._tThWidth, + tThMax=Material.DFLT_TTHMAX, + ) self.set_default_exclusions() @@ -343,34 +364,37 @@ def enable_hkls_below_index(self, index=5): def enable_hkls_below_tth(self, tth_threshold=90.0): ''' - enable reflections with two-theta less than @tth_threshold degrees + enable reflections with two-theta less than @tth_threshold degrees ''' tth_threshold = np.radians(tth_threshold) - tth = np.array([hkldata['tTheta'] - for hkldata in self._pData.hklDataList]) + tth = np.array( + [hkldata['tTheta'] for hkldata in self._pData.hklDataList] + ) dflt_excl = np.ones(tth.shape, dtype=np.bool) dflt_excl2 = np.ones(tth.shape, dtype=np.bool) - if(hasattr(self, 'hkl_from_file')): + if hasattr(self, 'hkl_from_file'): """ hkls were read from the file so the exclusions will be set based on what hkls are present """ for i, g in enumerate(self._pData.hklDataList): - if(g['hkl'].tolist() in self.hkl_from_file.tolist()): + if g['hkl'].tolist() in self.hkl_from_file.tolist(): dflt_excl[i] = False - dflt_excl2[~np.isnan(tth)] = \ - ~((tth[~np.isnan(tth)] >= 0.0) & - (tth[~np.isnan(tth)] <= tth_threshold)) + dflt_excl2[~np.isnan(tth)] = ~( + (tth[~np.isnan(tth)] >= 0.0) + & (tth[~np.isnan(tth)] <= tth_threshold) + ) dflt_excl = np.logical_or(dflt_excl, dflt_excl2) else: - dflt_excl[~np.isnan(tth)] = \ - ~((tth[~np.isnan(tth)] >= 0.0) & - (tth[~np.isnan(tth)] <= tth_threshold)) + dflt_excl[~np.isnan(tth)] = ~( + (tth[~np.isnan(tth)] >= 0.0) + & (tth[~np.isnan(tth)] <= tth_threshold) + ) dflt_excl[0] = False self._pData.exclusions = dflt_excl @@ -380,10 +404,9 @@ def update_structure_factor(self): sf = self.unitcell.CalcXRSF(hkls) self.planeData.set_structFact(sf) - def compute_powder_overlay(self, - ttharray=np.linspace(0, 80, 2000), - fwhm=0.25, - scale=1.0): + def compute_powder_overlay( + self, ttharray=np.linspace(0, 80, 2000), fwhm=0.25, scale=1.0 + ): """ this function computes a simulated spectra for using in place of lines for the powder @@ -398,7 +421,7 @@ def compute_powder_overlay(self, self.powder_overlay = np.zeros_like(ttharray) for t, I in zip(tth, Ip): p = [t, fwhm] - self.powder_overlay += scale*I*_unit_gaussian(p, ttharray) + self.powder_overlay += scale * I * _unit_gaussian(p, ttharray) def remove_duplicate_atoms(self): """ @@ -412,7 +435,6 @@ def remove_duplicate_atoms(self): self.charge = self.unitcell.chargestates self._hkls_changed() - def vt(self, temperature=None): '''calculate volume at high temperature @@ -422,10 +444,9 @@ def vt(self, temperature=None): if temperature is None: vt = self.v0 else: - delT = (temperature-298) - delT2 = (temperature**2-298**2) - vt = self.v0*np.exp(alpha0*delT - +0.5*alpha1*delT2) + delT = temperature - 298 + delT2 = temperature**2 - 298**2 + vt = self.v0 * np.exp(alpha0 * delT + 0.5 * alpha1 * delT2) return vt def kt(self, temperature=None): @@ -436,8 +457,8 @@ def kt(self, temperature=None): if temperature is None: return k0 else: - delT = (temperature-298) - return k0 + self.dk0dt*delT + delT = temperature - 298 + return k0 + self.dk0dt * delT def ktp(self, temperature=None): '''calculate bulk modulus derivative @@ -447,23 +468,25 @@ def ktp(self, temperature=None): if temperature is None: return k0p else: - delT = (temperature-298) - return k0p + self.dk0pdt*delT + delT = temperature - 298 + return k0p + self.dk0pdt * delT def calc_pressure(self, volume=None, temperature=None): '''calculate the pressure given the volume - and temperature using the third order - birch-murnaghan equation of state. + and temperature using the third order + birch-murnaghan equation of state. ''' if volume is None: return 0 else: - vt = self.vt(temperature=temperature) - kt = self.kt(temperature=temperature) + vt = self.vt(temperature=temperature) + kt = self.kt(temperature=temperature) ktp = self.ktp(temperature=temperature) - f = 0.5*((vt/volume)**(2./3.) - 1) + f = 0.5 * ((vt / volume) ** (2.0 / 3.0) - 1) - return 3.0*kt*f*(1 - 1.5*(4 - ktp)*f)*(1 + 2*f)**2.5 + return ( + 3.0 * kt * f * (1 - 1.5 * (4 - ktp) * f) * (1 + 2 * f) ** 2.5 + ) def calc_volume(self, pressure=None, temperature=None): '''solve for volume in the birch-murnaghan EoS to @@ -471,27 +494,30 @@ def calc_volume(self, pressure=None, temperature=None): to the Material object as updated lattice constants. ''' vt = self.vt(temperature=temperature) - kt = self.kt(temperature=temperature) + kt = self.kt(temperature=temperature) ktp = self.ktp(temperature=temperature) if pressure is None: return vt else: - alpha = 0.75*(ktp - 4) - p = np.zeros([10,]) + alpha = 0.75 * (ktp - 4) + p = np.zeros( + [ + 10, + ] + ) p[0] = 1.0 - p[2] = (1 - 2*alpha)/alpha - p[4] = (alpha-1)/alpha - p[9] = -2*pressure/3/kt/alpha + p[2] = (1 - 2 * alpha) / alpha + p[4] = (alpha - 1) / alpha + p[9] = -2 * pressure / 3 / kt / alpha res = np.roots(p) res = res[np.isreal(res)] - res = 1/np.real(res)**3 + res = 1 / np.real(res) ** 3 - mask = np.logical_and(res >= 0., res <= 1.0) + mask = np.logical_and(res >= 0.0, res <= 1.0) res = res[mask] if len(res) != 1: - msg = ('more than one physically ' - 'reasonable solution found!') + msg = 'more than one physically ' 'reasonable solution found!' raise ValueError(msg) return res[0] * vt @@ -500,10 +526,9 @@ def calc_lp_factor(self, pressure=None, temperature=None): constants by. only the lengths will be modified, the angles will be kept constant. ''' - vt = self.vt(temperature=temperature) - vpt = self.calc_volume(pressure=pressure, - temperature=temperature) - return (vpt/vt)**(1./3.) + vt = self.vt(temperature=temperature) + vpt = self.calc_volume(pressure=pressure, temperature=temperature) + return (vpt / vt) ** (1.0 / 3.0) def calc_lp_at_PT(self, pressure=None, temperature=None): '''calculate the lattice parameters for a given @@ -511,13 +536,19 @@ def calc_lp_at_PT(self, pressure=None, temperature=None): is the main function which will be called from the GUI. ''' - f = self.calc_lp_factor(pressure=pressure, - temperature=temperature) - return np.array([f*self.a0, f*self.b0, f*self.c0, - self.alpha0, self.beta0, self.gamma0,]) - - def _readCif(self, fcif=DFLT_NAME+'.cif'): - + f = self.calc_lp_factor(pressure=pressure, temperature=temperature) + return np.array( + [ + f * self.a0, + f * self.b0, + f * self.c0, + self.alpha0, + self.beta0, + self.gamma0, + ] + ) + + def _readCif(self, fcif=DFLT_NAME + '.cif'): """ >> @AUTHOR: Saransh Singh, Lawrence Livermore National Lab, saransh1@llnl.gov @@ -536,7 +567,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): # read the file for k in cif.keys(): - if('_cell_length_a' in cif[k]): + if '_cell_length_a' in cif[k]: m = k break cifdata = cif[m] @@ -544,51 +575,58 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): # make sure the space group is present in the cif file, either as # international table number, hermann-maguain or hall symbol - sgkey = ['_space_group_IT_number', - '_symmetry_space_group_name_h-m', - '_symmetry_space_group_name_hall', - '_symmetry_Int_Tables_number'] + sgkey = [ + '_space_group_IT_number', + '_symmetry_space_group_name_h-m', + '_symmetry_space_group_name_hall', + '_symmetry_Int_Tables_number', + ] sgdata = False for key in sgkey: sgdata = sgdata or (key in cifdata) - if(sgdata): + if sgdata: skey = key break - if(not(sgdata)): + if not (sgdata): raise RuntimeError(' No space group information in CIF file! ') sgnum = 0 if skey is sgkey[0]: sgnum = int(cifdata[sgkey[0]]) - elif (skey is sgkey[1]): + elif skey is sgkey[1]: HM = cifdata[sgkey[1]] HM = HM.replace(" ", "") HM = HM.replace("_", "") sgnum = HM_to_sgnum[HM] - elif (skey is sgkey[2]): + elif skey is sgkey[2]: hall = cifdata[sgkey[2]] hall = hall.replace(" ", "") sgnum = Hall_to_sgnum[HM] - elif(skey is sgkey[3]): + elif skey is sgkey[3]: sgnum = int(cifdata[sgkey[3]]) # lattice parameters lparms = [] - lpkey = ['_cell_length_a', '_cell_length_b', - '_cell_length_c', '_cell_angle_alpha', - '_cell_angle_beta', '_cell_angle_gamma'] + lpkey = [ + '_cell_length_a', + '_cell_length_b', + '_cell_length_c', + '_cell_angle_alpha', + '_cell_angle_beta', + '_cell_angle_gamma', + ] for key in lpkey: n = cifdata[key].find('(') - if(n != -1): + if n != -1: lparms.append(float(cifdata[key][:n])) else: lparms.append(float(cifdata[key])) for i in range(6): - if(i < 3): + if i < 3: lparms[i] = _angstroms(lparms[i]) else: lparms[i] = _degrees(lparms[i]) @@ -597,20 +635,27 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): self.sgnum = sgnum # fractional atomic site, occ and vibration amplitude - fracsitekey = ['_atom_site_fract_x', '_atom_site_fract_y', - '_atom_site_fract_z', ] - - occ_U = ['_atom_site_occupancy', - '_atom_site_u_iso_or_equiv', '_atom_site_U_iso_or_equiv'] + fracsitekey = [ + '_atom_site_fract_x', + '_atom_site_fract_y', + '_atom_site_fract_z', + ] + + occ_U = [ + '_atom_site_occupancy', + '_atom_site_u_iso_or_equiv', + '_atom_site_U_iso_or_equiv', + ] sitedata = True for key in fracsitekey: sitedata = sitedata and (key in cifdata) - if(not(sitedata)): + if not (sitedata): raise RuntimeError( ' fractional site position is not present \ - or incomplete in the CIF file! ') + or incomplete in the CIF file! ' + ) atompos = [] for key in fracsitekey: @@ -620,7 +665,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): for p in slist: n = p.find('(') - if(n != -1): + if n != -1: pos.append(p[:n]) else: pos.append(p) @@ -630,7 +675,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): bring them back to fractional coordinates between 0-1 ''' pos = np.asarray(pos).astype(np.float64) - pos, _ = np.modf(pos+100.0) + pos, _ = np.modf(pos + 100.0) atompos.append(pos) """note that the vibration amplitude, U is just the amplitude (in A) @@ -640,10 +685,10 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): this will be done here so we dont have to worry about it later """ - pocc = (occ_U[0] in cifdata.keys()) + pocc = occ_U[0] in cifdata.keys() pU = (occ_U[1] in cifdata.keys()) or (occ_U[2] in cifdata.keys()) - if(not pocc): + if not pocc: warn('occupation fraction not present. setting it to 1') occ = np.ones(atompos[0].shape) atompos.append(occ) @@ -653,7 +698,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): for p in slist: n = p.find('(') - if(n != -1): + if n != -1: occ.append(p[:n]) else: occ.append(p) @@ -664,14 +709,16 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): atompos.append(np.asarray(occstr).astype(np.float64)) - if(not pU): - msg = (f"'Debye-Waller factors not present. " - f"setting to same values for all atoms.'") + if not pU: + msg = ( + "'Debye-Waller factors not present. " + "setting to same values for all atoms.'" + ) warn(msg) U = self.DFLT_U[0] * np.ones(atompos[0].shape) self._U = U else: - if(occ_U[1] in cifdata.keys()): + if occ_U[1] in cifdata.keys(): k = occ_U[1] else: k = occ_U[2] @@ -681,7 +728,7 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): for p in slist: n = p.find('(') - if(n != -1): + if n != -1: U.append(p[:n]) else: U.append(p) @@ -705,8 +752,8 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): get atome types here i.e. the atomic number of atoms at each site ''' atype = '_atom_site_type_symbol' - patype = (atype in cifdata) - if(not patype): + patype = atype in cifdata + if not patype: raise RuntimeError('atom types not defined in cif file.') satype = cifdata[atype] @@ -749,7 +796,6 @@ def _readCif(self, fcif=DFLT_NAME+'.cif'): self.alphaT = 0.0 self.dalphaTdT = 0.0 - def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): """ >> @AUTHOR: Saransh Singh, Lawrence Livermore National Lab, @@ -777,12 +823,11 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): raise Exception(f'Unknown type for fhdf: {fhdf}') if xtal not in root_gid: - raise IOError('crystal doesn''t exist in material file.') + raise IOError('crystal doesn' 't exist in material file.') gid = root_gid.get(xtal) - sgnum = np.array(gid.get('SpaceGroupNumber'), - dtype=np.int32).item() + sgnum = np.array(gid.get('SpaceGroupNumber'), dtype=np.int32).item() """ IMPORTANT NOTE: note that the latice parameters is nm by default @@ -793,8 +838,8 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): lparms = list(gid.get('LatticeParameters')) for i in range(6): - if(i < 3): - lparms[i] = _angstroms(lparms[i]*10.0) + if i < 3: + lparms[i] = _angstroms(lparms[i] * 10.0) else: lparms[i] = _degrees(lparms[i]) @@ -804,29 +849,31 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): self.sgnum = sgnum # the U factors are related to B by the relation B = 8pi^2 U - self._atominfo = np.transpose(np.array( - gid.get('AtomData'), dtype=np.float64)) - self._U = np.transpose(np.array( - gid.get('U'), dtype=np.float64)) + self._atominfo = np.transpose( + np.array(gid.get('AtomData'), dtype=np.float64) + ) + self._U = np.transpose(np.array(gid.get('U'), dtype=np.float64)) # read atom types (by atomic number, Z) self._atomtype = np.array(gid.get('Atomtypes'), dtype=np.int32) if 'ChargeStates' in gid: self._charge = h5py_read_string(gid['ChargeStates']) else: - self._charge = ['0']*self._atomtype.shape[0] + self._charge = ['0'] * self._atomtype.shape[0] self._atom_ntype = self._atomtype.shape[0] - self._sgsetting = np.array(gid.get('SpaceGroupSetting'), - dtype=np.int32).item() + self._sgsetting = np.array( + gid.get('SpaceGroupSetting'), dtype=np.int32 + ).item() - if('stiffness' in gid): + if 'stiffness' in gid: # we're assuming the stiffness is in units of GPa self.stiffness = np.array(gid.get('stiffness')) - elif('compliance' in gid): + elif 'compliance' in gid: # we're assuming the compliance is in units of TPa^-1 - self.stiffness = np.linalg.inv( - np.array(gid.get('compliance'))) * 1.e3 + self.stiffness = ( + np.linalg.inv(np.array(gid.get('compliance'))) * 1.0e3 + ) else: self.stiffness = np.zeros([6, 6]) @@ -834,91 +881,80 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): parameters ''' self.k0 = 100.0 - if('k0' in gid): + if 'k0' in gid: # this is the isotropic bulk modulus - k0 = np.array(gid.get('k0'), - dtype=np.float64).item() + k0 = np.array(gid.get('k0'), dtype=np.float64).item() self.k0 = k0 self.k0p = 0.0 - if('k0p' in gid): + if 'k0p' in gid: # this is the pressure derivation of # the isotropic bulk modulus - k0p = np.array(gid.get('k0p'), - dtype=np.float64).item() + k0p = np.array(gid.get('k0p'), dtype=np.float64).item() self.k0p = k0p self.dk0dt = 0.0 - if('dk0dt' in gid): + if 'dk0dt' in gid: # this is the temperature derivation of # the isotropic bulk modulus - dk0dt = np.array(gid.get('dk0dt'), - dtype=np.float64).item() + dk0dt = np.array(gid.get('dk0dt'), dtype=np.float64).item() self.dk0dt = dk0dt self.dk0pdt = 0.0 - if('dk0pdt' in gid): + if 'dk0pdt' in gid: # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - dk0pdt = np.array(gid.get('dk0pdt'), - dtype=np.float64).item() + dk0pdt = np.array(gid.get('dk0pdt'), dtype=np.float64).item() self.dk0pdt = dk0pdt self.alphaT = 0.0 - if('alphaT' in gid): + if 'alphaT' in gid: # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - alphaT = np.array(gid.get('alphaT'), - dtype=np.float64).item() + alphaT = np.array(gid.get('alphaT'), dtype=np.float64).item() self.alphaT = alphaT self.dalphaTdT = 0.0 - if('dalphaTdT' in gid): + if 'dalphaTdT' in gid: # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - dalphaTdT = np.array(gid.get('dalphaTdT'), - dtype=np.float64).item() + dalphaTdT = np.array(gid.get('dalphaTdT'), dtype=np.float64).item() self.dalphaTdT = dalphaTdT '''Finished with the BM EOS ''' - if('v0' in gid): + if 'v0' in gid: # this is the isotropic ambient unitcell # volume - v0 = np.array(gid.get('v0'), - dtype=np.float64).item() + v0 = np.array(gid.get('v0'), dtype=np.float64).item() self.v0 = v0 # if dmin is present in file: - if('dmin' in gid): + if 'dmin' in gid: # if dmin is present in the HDF5 file, then use that - dmin = np.array(gid.get('dmin'), - dtype=np.float64).item() - self._dmin = _angstroms(dmin*10.) + dmin = np.array(gid.get('dmin'), dtype=np.float64).item() + self._dmin = _angstroms(dmin * 10.0) else: self._dmin = Material.DFLT_DMIN # if kev is present in file - if('kev' in gid): - kev = np.array(gid.get('kev'), - dtype=np.float64).item() + if 'kev' in gid: + kev = np.array(gid.get('kev'), dtype=np.float64).item() self._beamEnergy = _kev(kev) else: self._beamEnergy = Material.DFLT_KEV if 'tThWidth' in gid: - tThWidth = np.array(gid.get('tThWidth'), - dtype=np.float64).item() + tThWidth = np.array(gid.get('tThWidth'), dtype=np.float64).item() tThWidth = np.radians(tThWidth) else: tThWidth = Material.DFLT_TTH self._tThWidth = tThWidth - if('hkls' in gid): - self.hkl_from_file = np.array(gid.get('hkls'), - dtype=np.int32) + if 'hkls' in gid: + self.hkl_from_file = np.array(gid.get('hkls'), dtype=np.int32) if isinstance(fhdf, (Path, str)): # Close the file... @@ -978,12 +1014,14 @@ def dump_material(self, file, path=None): ''' lattice parameters ''' - lat_param = {'a': self.unitcell.a, - 'b': self.unitcell.b, - 'c': self.unitcell.c, - 'alpha': self.unitcell.alpha, - 'beta': self.unitcell.beta, - 'gamma': self.unitcell.gamma} + lat_param = { + 'a': self.unitcell.a, + 'b': self.unitcell.b, + 'c': self.unitcell.c, + 'alpha': self.unitcell.alpha, + 'beta': self.unitcell.beta, + 'gamma': self.unitcell.gamma, + } Write2H5File(AtomInfo, lat_param, path) @@ -1010,12 +1048,16 @@ def spaceGroup(self): @property def vol(self): - return self.unitcell.vol*1e3 + return self.unitcell.vol * 1e3 @property def lparms(self): - return np.array([x.getVal("nm") if x.isLength() else - x.getVal("degrees") for x in self._lparms]) + return np.array( + [ + x.getVal("nm") if x.isLength() else x.getVal("degrees") + for x in self._lparms + ] + ) @property def latticeType(self): @@ -1064,15 +1106,12 @@ def _set_sgnum(self, v): self._newUnitcell() self._hkls_changed() - sgnum = property(_get_sgnum, _set_sgnum, None, - "Space group number") + sgnum = property(_get_sgnum, _set_sgnum, None, "Space group number") @property def pgnum(self): return self.unitcell.pgnum - # property: beamEnergy - def _get_beamEnergy(self): """Get method for beamEnergy""" return self._beamEnergy @@ -1085,7 +1124,7 @@ def _set_beamEnergy(self, keV): float arguments. Also can take a valWUnit instance """ - if(isinstance(keV, valWUnit)): + if isinstance(keV, valWUnit): self._beamEnergy = keV else: self._beamEnergy = valWUnit('kev', 'energy', keV, 'keV') @@ -1095,16 +1134,18 @@ def _set_beamEnergy(self, keV): @TODO make voltage valWUnit instance so this that we dont have to track this manually inn the future ''' - self.unitcell.voltage = self.beamEnergy.value*1e3 + self.unitcell.voltage = self.beamEnergy.value * 1e3 self.planeData.wavelength = keV self._hkls_changed() - beamEnergy = property(_get_beamEnergy, _set_beamEnergy, None, - "Beam energy in keV") + beamEnergy = property( + _get_beamEnergy, _set_beamEnergy, None, "Beam energy in keV" + ) """ 03/11/2021 SS 1.0 original """ + @property def unitcell(self): return self._unitcell @@ -1122,7 +1163,7 @@ def latticeParameters(self): @latticeParameters.setter def latticeParameters(self, v): """Set method for latticeParameters""" - if(len(v) != 6): + if len(v) != 6: v = unitcell._rqpDict[self.unitcell.latticeType][1](v) lp = [_angstroms(v[i]) for i in range(3)] for i in range(3, 6): @@ -1146,8 +1187,7 @@ def _set_name(self, v): return - name = property(_get_name, _set_name, None, - "Name of material") + name = property(_get_name, _set_name, None, "Name of material") @property def dmin(self): @@ -1175,8 +1215,10 @@ def charge(self, vals): first make sure the lengths are correct """ if len(vals) != len(self.atomtype): - msg = (f"incorrect size of charge. " - f"must be same size as number of aom types.") + msg = ( + "incorrect size of charge. " + "must be same size as number of aom types." + ) raise ValueError(msg) """ now check if charge states are actually allowed @@ -1185,9 +1227,11 @@ def charge(self, vals): elem = ptableinverse[z] cs = chargestate[elem] if not vals[ii] in cs: - msg = (f"element {elem} does not allow " - f"charge state of {vals[ii]}. " - f"allowed charge states are {cs}.") + msg = ( + f"element {elem} does not allow " + f"charge state of {vals[ii]}. " + f"allowed charge states are {cs}." + ) raise ValueError(msg) self._charge = vals @@ -1218,8 +1262,11 @@ def _set_atominfo(self, v): self._atominfo = v atominfo = property( - _get_atominfo, _set_atominfo, None, - "Information about atomic positions and electron number") + _get_atominfo, + _set_atominfo, + None, + "Information about atomic positions and electron number", + ) # property: "atomtype" def _get_atomtype(self): @@ -1231,8 +1278,8 @@ def _set_atomtype(self, v): self._atomtype = v atomtype = property( - _get_atomtype, _set_atomtype, None, - "Information about atomic types") + _get_atomtype, _set_atomtype, None, "Information about atomic types" + ) @property def thermal_expansion(self): @@ -1274,27 +1321,35 @@ def _set_atomdata(self, atomtype, atominfo, U, charge): U = np.array(U) if atomtype.shape[0] != atominfo.shape[0]: - msg = (f"inconsistent shapes: number of atoms " - f"types passed = {atomtype.shape[0]} \n" - f" number of atom positions passed = {atominfo.shape[0]}") + msg = ( + f"inconsistent shapes: number of atoms " + f"types passed = {atomtype.shape[0]} \n" + f" number of atom positions passed = {atominfo.shape[0]}" + ) raise ValueError(msg) if atomtype.shape[0] != U.shape[0]: - msg = (f"inconsistent shapes: number of atoms " - f"types passed = {atomtype.shape[0]} \n" - f" U passed for {U.shape[0]} atoms.") + msg = ( + f"inconsistent shapes: number of atoms " + f"types passed = {atomtype.shape[0]} \n" + f" U passed for {U.shape[0]} atoms." + ) raise ValueError(msg) if atominfo.shape[0] != U.shape[0]: - msg = (f"inconsistent shapes: number of atom " - f"positions passed = {atominfo.shape[0]} \n" - f"U passed for {U.shape[0]} atoms.") + msg = ( + f"inconsistent shapes: number of atom " + f"positions passed = {atominfo.shape[0]} \n" + f"U passed for {U.shape[0]} atoms." + ) raise ValueError(msg) if len(charge) != atomtype.shape[0]: - msg = (f"inconsistent shapes: number of atoms " - f"types passed = {atomtype.shape[0]} \n" - f"charge value passed for {len(charge)} atoms.") + msg = ( + f"inconsistent shapes: number of atoms " + f"types passed = {atomtype.shape[0]} \n" + f"charge value passed for {len(charge)} atoms." + ) raise ValueError(msg) self.atomtype = atomtype @@ -1304,6 +1359,7 @@ def _set_atomdata(self, atomtype, atominfo, U, charge): self._newUnitcell() self.update_structure_factor() + # # ========== Methods # @@ -1321,8 +1377,7 @@ def _set_atomdata(self, atomtype, atominfo, U, charge): def loadMaterialList(cfgFile): """Load a list of materials from a file - The file uses the config file format. See ConfigParser module. -""" + The file uses the config file format. See ConfigParser module.""" p = Parser() p.read(cfgFile) # @@ -1336,8 +1391,9 @@ def loadMaterialList(cfgFile): return matList -def load_materials_hdf5(f, dmin=None, kev=None, - sgsetting=Material.DFLT_SGSETTING): +def load_materials_hdf5( + f, dmin=None, kev=None, sgsetting=Material.DFLT_SGSETTING +): """Load materials from an HDF5 file The file uses the HDF5 file format. @@ -1362,10 +1418,7 @@ def load_materials_hdf5(f, dmin=None, kev=None, 'kev': kev, 'sgsetting': sgsetting, } - return { - name: Material(name, **kwargs) - for name in names - } + return {name: Material(name, **kwargs) for name in names} def save_materials_hdf5(f, materials, path=None): @@ -1378,6 +1431,7 @@ def hkls_match(a, b): # Check if hkls match. Expects inputs to have shape (x, 3). def sorted_hkls(x): return x[np.lexsort((x[:, 2], x[:, 1], x[:, 0]))] + return np.array_equal(sorted_hkls(a), sorted_hkls(b)) @@ -1394,6 +1448,7 @@ def get_hkl_strings(pdata): } return np.intersect1d(**kwargs)[1:] + # # ============================== Executable section for testing # From eea778b700145f088938358a56eb4bb95e346491 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 30 Nov 2023 09:32:27 -0600 Subject: [PATCH 19/26] Initialize v0 in the __init__ method if needed It is better to initialize it here instead of in the newUnitCell method, because it should only need to happen once: when the material is initialized. Signed-off-by: Patrick Avery --- hexrd/material/material.py | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index 9abfbe520..68088b318 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -219,6 +219,10 @@ def __init__( self._beamEnergy = kev self._newUnitcell() + + if not hasattr(self, 'v0'): + self.v0 = self.vol + self._newPdata() self.update_structure_factor() @@ -282,9 +286,6 @@ def _newUnitcell(self): if pdata.laueGroup != laue or pdata.lparms != reduced_lparms: pdata.set_laue_and_lparms(laue, reduced_lparms) - if not hasattr(self, 'v0'): - self.v0 = self.vol - def _hkls_changed(self): # Call this when something happens that changes the hkls... self._newPdata() From bdd749a0050060da1db2741fe097182fb354311a Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 30 Nov 2023 10:22:22 -0600 Subject: [PATCH 20/26] Several fixes and improvements to use with the GUI Signed-off-by: Patrick Avery --- hexrd/material/jcpds.py | 43 ++++++++---------- hexrd/material/material.py | 87 ++++++++++++++++++++++++------------- hexrd/material/mksupport.py | 19 +++++--- 3 files changed, 90 insertions(+), 59 deletions(-) diff --git a/hexrd/material/jcpds.py b/hexrd/material/jcpds.py index 39b21a868..f7e672890 100644 --- a/hexrd/material/jcpds.py +++ b/hexrd/material/jcpds.py @@ -19,8 +19,8 @@ def __init__(self, filename=None): self.dk0pdt = 0 self.symmetry = '' - self.thermal_expansion = 0 # alphat at 298K - self.thermal_expansion_dt = 0 + self.alpha_t = 0 # alphat at 298K + self.dalpha_t_dt = 0 self.file = ' ' self.name = ' ' @@ -142,10 +142,10 @@ def read_file(self, file): item = str.split(inp[4]) if self.version == 3: - thermal_expansion = 0. + alpha_t = 0. else: - thermal_expansion = float(item[0]) - self.thermal_expansion = thermal_expansion + alpha_t = float(item[0]) + self.alpha_t = alpha_t version_status = 'new' @@ -215,11 +215,11 @@ def read_file(self, file): if jlinespl[0] == 'ALPHAT:': alphat = float(jlinespl[1]) - self.thermal_expansion = alphat + self.alpha_t = alphat if jlinespl[0] == 'DALPHATDT:': dalphatdt = float(jlinespl[1]) - self.thermal_expansion_dt = dalphatdt + self.dalpha_t_dt = dalphatdt if jlinespl[0] == 'DIHKL:': pass @@ -253,7 +253,7 @@ def read_file(self, file): elif self.symmetry == 'monoclinic': self.alpha0 = 90. self.gamma0 = 90. - #elif self.symmetry == 'triclinic': + # elif self.symmetry == 'triclinic': jcpdsfile.close() @@ -284,13 +284,6 @@ def _lat_param_names(self): def lattice_params(self): return [getattr(self, x) for x in self._lat_param_names] - @lattice_params.setter - def lattice_params(self, v): - for name, val in zip(self._lat_param_names, v): - setattr(self, name, val) - - self.update_v0() - def matches_material(self, mat): self.verify_symmetry_match(mat) mat_lp = [x.value for x in mat.latticeParameters] @@ -298,16 +291,18 @@ def matches_material(self, mat): if not np.isclose(v1, v2): return False - def load_from_material(self, mat): - if self.symmetry: - self.verify_symmetry_match(mat) - - lp = [x.value for x in mat.latticeParameters] - self.lattice_params = lp + def write_lattice_params_to_material(self, mat): + self.verify_symmetry_match(mat) + mat.latticeParameters0 = self.lattice_params - def write_to_material(self, mat): + def write_pt_params_to_material(self, mat): self.verify_symmetry_match(mat) - mat.latticeParameters = self.lattice_params + mat.k0 = self.k0 + mat.k0p = self.k0p + mat.dk0dt = self.dk0dt + mat.dk0pdt = self.dk0pdt + mat.alpha_t = self.alpha_t + mat.dalpha_t_dt = self.dalpha_t_dt def update_v0(self): self.v0 = self.calc_volume_unitcell() @@ -323,7 +318,7 @@ def calc_volume_unitcell(self): v0 = self.a0*self.b0*self.c0 f = np.sqrt(1 - ca**2 - cb**2 - cg**2 - +2*ca*cb*cg) + + 2 * ca * cb * cg) return v0*f diff --git a/hexrd/material/material.py b/hexrd/material/material.py index 68088b318..d6ce17b53 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -203,12 +203,15 @@ def __init__( self._dmin = Material.DFLT_DMIN self._beamEnergy = Material.DFLT_KEV + self.pressure = 0 + self.temperature = 298 + self.pt_lp_factor = 0 self.k0 = 100.0 self.k0p = 0.0 self.dk0dt = 0.0 self.dk0pdt = 0.0 - self.alphaT = 0.0 - self.dalphaTdT = 0.0 + self.alpha_t = 0.0 + self.dalpha_t_dt = 0.0 # If these were specified, they override any other method of # obtaining them (including loading them from files). @@ -221,7 +224,7 @@ def __init__( self._newUnitcell() if not hasattr(self, 'v0'): - self.v0 = self.vol + self.update_v0() self._newPdata() self.update_structure_factor() @@ -339,6 +342,9 @@ def _newPdata(self): self.set_default_exclusions() + def update_v0(self): + self.v0 = self.vol + def set_default_exclusions(self): if hasattr(self, 'hkl_from_file'): # If we loaded hkls from the file, use those @@ -531,7 +537,7 @@ def calc_lp_factor(self, pressure=None, temperature=None): vpt = self.calc_volume(pressure=pressure, temperature=temperature) return (vpt / vt) ** (1.0 / 3.0) - def calc_lp_at_PT(self, pressure=None, temperature=None): + def calc_lp_at_PT(self, lparms0, pressure=None, temperature=None): '''calculate the lattice parameters for a given pressure and temperature using the BM EoS. This is the main function which will be called from @@ -540,12 +546,8 @@ def calc_lp_at_PT(self, pressure=None, temperature=None): f = self.calc_lp_factor(pressure=pressure, temperature=temperature) return np.array( [ - f * self.a0, - f * self.b0, - f * self.c0, - self.alpha0, - self.beta0, - self.gamma0, + *(f * lparms0[:3]), + *lparms0[3:], ] ) @@ -790,12 +792,15 @@ def _readCif(self, fcif=DFLT_NAME + '.cif'): parameters to default values. These values can be updated by user or by reading a JCPDS file ''' + self.pressure = 0 + self.temperature = 298 + self.pt_lp_factor = 1 self.k0 = 100.0 self.k0p = 0.0 self.dk0dt = 0.0 self.dk0pdt = 0.0 - self.alphaT = 0.0 - self.dalphaTdT = 0.0 + self.alpha_t = 0.0 + self.dalpha_t_dt = 0.0 def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): """ @@ -881,6 +886,21 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): '''start reading the Birch-Murnaghan equation of state parameters ''' + self.pressure = 0 + if 'pressure' in gid: + self.pressure = np.array(gid.get('pressure'), + dtype=np.float64).item() + + self.temperature = 298 + if 'temperature' in gid: + self.temperature = np.array(gid.get('temperature'), + dtype=np.float64).item() + + self.pt_lp_factor = 1 + if 'pt_lp_factor' in gid: + self.pt_lp_factor = np.array(gid.get('pt_lp_factor'), + dtype=np.float64).item() + self.k0 = 100.0 if 'k0' in gid: # this is the isotropic bulk modulus @@ -908,19 +928,19 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): dk0pdt = np.array(gid.get('dk0pdt'), dtype=np.float64).item() self.dk0pdt = dk0pdt - self.alphaT = 0.0 - if 'alphaT' in gid: + self.alpha_t = 0.0 + if 'alpha_t' in gid: # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - alphaT = np.array(gid.get('alphaT'), dtype=np.float64).item() - self.alphaT = alphaT + alpha_t = np.array(gid.get('alpha_t'), dtype=np.float64).item() + self.alpha_t = alpha_t - self.dalphaTdT = 0.0 - if 'dalphaTdT' in gid: + self.dalpha_t_dt = 0.0 + if 'dalpha_t_dt' in gid: # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - dalphaTdT = np.array(gid.get('dalphaTdT'), dtype=np.float64).item() - self.dalphaTdT = dalphaTdT + dalpha_t_dt = np.array(gid.get('dalpha_t_dt'), dtype=np.float64).item() + self.dalpha_t_dt = dalpha_t_dt '''Finished with the BM EOS ''' @@ -985,6 +1005,10 @@ def dump_material(self, file, path=None): else: AtomInfo['tThWidth'] = np.degrees(self.planeData.tThWidth) + AtomInfo['pressure'] = self.pressure + AtomInfo['temperature'] = self.temperature + AtomInfo['pt_lp_factor'] = self.pt_lp_factor + AtomInfo['k0'] = 100.0 if hasattr(self, 'k0'): AtomInfo['k0'] = self.k0 @@ -1005,13 +1029,13 @@ def dump_material(self, file, path=None): if hasattr(self, 'dk0pdt'): AtomInfo['dk0pdt'] = self.dk0pdt - AtomInfo['alphaT'] = 0.0 - if hasattr(self, 'alphaT'): - AtomInfo['alphaT'] = self.alphaT + AtomInfo['alpha_t'] = 0.0 + if hasattr(self, 'alpha_t'): + AtomInfo['alpha_t'] = self.alpha_t - AtomInfo['dalphaTdT'] = 0.0 - if hasattr(self, 'dalphaTdT'): - AtomInfo['dalphaTdT'] = self.dalphaTdT + AtomInfo['dalpha_t_dt'] = 0.0 + if hasattr(self, 'dalpha_t_dt'): + AtomInfo['dalpha_t_dt'] = self.dalpha_t_dt ''' lattice parameters ''' @@ -1284,19 +1308,19 @@ def _set_atomtype(self, v): @property def thermal_expansion(self): - return self.alphaT + return self.alpha_t @thermal_expansion.setter def thermal_expansion(self, val): - self.alphaT = val + self.alpha_t = val @property def thermal_expansion_dt(self): - return self.dalphaTdT + return self.dalpha_t_dt @thermal_expansion_dt.setter def thermal_expansion_dt(self, val): - self.dalphaTdT = val + self.dalpha_t_dt = val def _set_atomdata(self, atomtype, atominfo, U, charge): """ @@ -1433,6 +1457,9 @@ def hkls_match(a, b): def sorted_hkls(x): return x[np.lexsort((x[:, 2], x[:, 1], x[:, 0]))] + if a.shape != b.shape: + return False + return np.array_equal(sorted_hkls(a), sorted_hkls(b)) diff --git a/hexrd/material/mksupport.py b/hexrd/material/mksupport.py index c3b44c2de..edf0223cc 100644 --- a/hexrd/material/mksupport.py +++ b/hexrd/material/mksupport.py @@ -515,6 +515,15 @@ def WriteH5Data(fid, AtomInfo, lat_param, path=None): did = gid.create_dataset("stiffness", (6, 6), dtype=np.float64) did.write_direct(np.array(AtomInfo['stiffness'], dtype=np.float64)) + did = gid.create_dataset("pressure", (1,), dtype=np.float64) + did.write_direct(np.array(AtomInfo['pressure'], dtype=np.float64)) + + did = gid.create_dataset("temperature", (1,), dtype=np.float64) + did.write_direct(np.array(AtomInfo['temperature'], dtype=np.float64)) + + did = gid.create_dataset("pt_lp_factor", (1,), dtype=np.float64) + did.write_direct(np.array(AtomInfo['pt_lp_factor'], dtype=np.float64)) + did = gid.create_dataset("k0", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['k0']], dtype=np.float64)) @@ -527,16 +536,16 @@ def WriteH5Data(fid, AtomInfo, lat_param, path=None): did = gid.create_dataset("dk0pdt", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['dk0pdt']], dtype=np.float64)) - did = gid.create_dataset("alphaT", (1,), dtype=np.float64) - did.write_direct(np.array([AtomInfo['alphaT']], dtype=np.float64)) + did = gid.create_dataset("alpha_t", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['alpha_t']], dtype=np.float64)) - did = gid.create_dataset("dalphaTdT", (1,), dtype=np.float64) - did.write_direct(np.array([AtomInfo['dalphaTdT']], dtype=np.float64)) + did = gid.create_dataset("dalpha_t_dt", (1,), dtype=np.float64) + did.write_direct(np.array([AtomInfo['dalpha_t_dt']], dtype=np.float64)) did = gid.create_dataset("v0", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['v0']], dtype=np.float64)) - if "tThWidth" in AtomInfo: + if "tThWidth" in AtomInfo: did = gid.create_dataset("tThWidth", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['tThWidth']], dtype=np.float64)) From 6e2e9bd9042c4a2e81292c9cbe547739c96b68fd Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Tue, 19 Dec 2023 13:27:19 -0600 Subject: [PATCH 21/26] Set default LP factor to 1 Signed-off-by: Patrick Avery --- hexrd/material/material.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index d6ce17b53..127da6938 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -205,7 +205,7 @@ def __init__( self.pressure = 0 self.temperature = 298 - self.pt_lp_factor = 0 + self.pt_lp_factor = 1 self.k0 = 100.0 self.k0p = 0.0 self.dk0dt = 0.0 From 8baa384fc0c6868b0237be5d9bb80b41585c65c7 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Tue, 19 Dec 2023 16:02:14 -0600 Subject: [PATCH 22/26] Compute lparms0 and pt_lp_factor These can be computed and don't need to be stored. Signed-off-by: Patrick Avery --- hexrd/material/material.py | 24 +++++++++++++++--------- hexrd/material/mksupport.py | 3 --- 2 files changed, 15 insertions(+), 12 deletions(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index 127da6938..cd3fa1e26 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -205,7 +205,6 @@ def __init__( self.pressure = 0 self.temperature = 298 - self.pt_lp_factor = 1 self.k0 = 100.0 self.k0p = 0.0 self.dk0dt = 0.0 @@ -478,6 +477,19 @@ def ktp(self, temperature=None): delT = temperature - 298 return k0p + self.dk0pdt * delT + @property + def pt_lp_factor(self): + return (self.unitcell.vol * 1e3 / self.v0) ** (1 / 3) + + @property + def lparms0(self): + # Get the lattice parameters for 0 pressure and temperature (at v0) + lparms = self.lparms + return np.array([ + *(lparms[:3] / self.pt_lp_factor), + *lparms[3:], + ]) + def calc_pressure(self, volume=None, temperature=None): '''calculate the pressure given the volume and temperature using the third order @@ -537,13 +549,14 @@ def calc_lp_factor(self, pressure=None, temperature=None): vpt = self.calc_volume(pressure=pressure, temperature=temperature) return (vpt / vt) ** (1.0 / 3.0) - def calc_lp_at_PT(self, lparms0, pressure=None, temperature=None): + def calc_lp_at_PT(self, pressure=None, temperature=None): '''calculate the lattice parameters for a given pressure and temperature using the BM EoS. This is the main function which will be called from the GUI. ''' f = self.calc_lp_factor(pressure=pressure, temperature=temperature) + lparms0 = self.lparms0 return np.array( [ *(f * lparms0[:3]), @@ -794,7 +807,6 @@ def _readCif(self, fcif=DFLT_NAME + '.cif'): ''' self.pressure = 0 self.temperature = 298 - self.pt_lp_factor = 1 self.k0 = 100.0 self.k0p = 0.0 self.dk0dt = 0.0 @@ -896,11 +908,6 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): self.temperature = np.array(gid.get('temperature'), dtype=np.float64).item() - self.pt_lp_factor = 1 - if 'pt_lp_factor' in gid: - self.pt_lp_factor = np.array(gid.get('pt_lp_factor'), - dtype=np.float64).item() - self.k0 = 100.0 if 'k0' in gid: # this is the isotropic bulk modulus @@ -1007,7 +1014,6 @@ def dump_material(self, file, path=None): AtomInfo['pressure'] = self.pressure AtomInfo['temperature'] = self.temperature - AtomInfo['pt_lp_factor'] = self.pt_lp_factor AtomInfo['k0'] = 100.0 if hasattr(self, 'k0'): diff --git a/hexrd/material/mksupport.py b/hexrd/material/mksupport.py index edf0223cc..083272624 100644 --- a/hexrd/material/mksupport.py +++ b/hexrd/material/mksupport.py @@ -521,9 +521,6 @@ def WriteH5Data(fid, AtomInfo, lat_param, path=None): did = gid.create_dataset("temperature", (1,), dtype=np.float64) did.write_direct(np.array(AtomInfo['temperature'], dtype=np.float64)) - did = gid.create_dataset("pt_lp_factor", (1,), dtype=np.float64) - did.write_direct(np.array(AtomInfo['pt_lp_factor'], dtype=np.float64)) - did = gid.create_dataset("k0", (1,), dtype=np.float64) did.write_direct(np.array([AtomInfo['k0']], dtype=np.float64)) From 96ff679444bd120af7bc3143d69980aa7d9996f3 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Tue, 19 Dec 2023 16:10:09 -0600 Subject: [PATCH 23/26] Fix import paths for tests Signed-off-by: Patrick Avery --- tests/find_orientations_testing.py | 2 +- tests/test_find_orientations.py | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/find_orientations_testing.py b/tests/find_orientations_testing.py index b9ff117fa..d8a0907f3 100755 --- a/tests/find_orientations_testing.py +++ b/tests/find_orientations_testing.py @@ -11,7 +11,7 @@ import numpy as np -from hexrd.crystallography import PlaneData +from hexrd.material.crystallography import PlaneData from hexrd.rotations import misorientation diff --git a/tests/test_find_orientations.py b/tests/test_find_orientations.py index bb4a472fe..c6208f012 100644 --- a/tests/test_find_orientations.py +++ b/tests/test_find_orientations.py @@ -10,7 +10,7 @@ from hexrd.findorientations import find_orientations, generate_eta_ome_maps from hexrd import config -from hexrd.crystallography import PlaneData +from hexrd.material.crystallography import PlaneData import find_orientations_testing as test_utils From 9f036b816f141049d2d0284c221296238dfc1486 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Thu, 21 Dec 2023 20:05:24 -0500 Subject: [PATCH 24/26] Fix pep8 issues Signed-off-by: Patrick Avery --- hexrd/material/material.py | 3 ++- hexrd/wppf/phase.py | 10 ++++++---- 2 files changed, 8 insertions(+), 5 deletions(-) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index cd3fa1e26..dba3b652b 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -946,7 +946,8 @@ def _readHDFxtal(self, fhdf=DFLT_NAME, xtal=DFLT_NAME): if 'dalpha_t_dt' in gid: # this is the temperature derivation of # the pressure derivative of isotropic bulk modulus - dalpha_t_dt = np.array(gid.get('dalpha_t_dt'), dtype=np.float64).item() + dalpha_t_dt = np.array(gid.get('dalpha_t_dt'), + dtype=np.float64).item() self.dalpha_t_dt = dalpha_t_dt '''Finished with the BM EOS diff --git a/hexrd/wppf/phase.py b/hexrd/wppf/phase.py index b6bd5f58b..b90ca018c 100644 --- a/hexrd/wppf/phase.py +++ b/hexrd/wppf/phase.py @@ -22,7 +22,8 @@ class Material_LeBail: 09/14/2020 SS 1.1 class can now be initialized using a material.Material class instance >> @DETAILS: Material_LeBail class is a stripped down version of the - material.Material class.this is done to keep the class lightweight + material.Material class.this is done to keep the class + lightweight and make sure only the information necessary for the lebail fit is kept ========================================================================================= @@ -606,9 +607,10 @@ class Material_Rietveld: >> @DATE: 05/18/2020 SS 1.0 original 02/01/2021 SS 1.1 class can now be initialized using a material.Material class instance - >> @DETAILS: Material_LeBail class is a stripped down version of the material.Material - class.this is done to keep the class lightweight and make sure only the - information necessary for the Rietveld fit is kept + >> @DETAILS: Material_LeBail class is a stripped down version of the + material.Material class. This is done to keep the class + lightweight and make sure only the information necessary + for the Rietveld fit is kept =========================================================================================== ========================================================================================== """ From 7c28c6b74587c23247cfec9bc7c3adf9de553e2f Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Fri, 29 Dec 2023 15:38:45 -0600 Subject: [PATCH 25/26] Remove setting latticeParameters0 This property no longer exists, and we should just set the lattice parameters instead. Signed-off-by: Patrick Avery --- hexrd/material/jcpds.py | 3 ++- hexrd/material/material.py | 4 ++-- 2 files changed, 4 insertions(+), 3 deletions(-) diff --git a/hexrd/material/jcpds.py b/hexrd/material/jcpds.py index f7e672890..0affaa8c2 100644 --- a/hexrd/material/jcpds.py +++ b/hexrd/material/jcpds.py @@ -293,7 +293,8 @@ def matches_material(self, mat): def write_lattice_params_to_material(self, mat): self.verify_symmetry_match(mat) - mat.latticeParameters0 = self.lattice_params + mat.latticeParameters = self.lattice_params + mat.reset_v0() def write_pt_params_to_material(self, mat): self.verify_symmetry_match(mat) diff --git a/hexrd/material/material.py b/hexrd/material/material.py index dba3b652b..abb139f08 100644 --- a/hexrd/material/material.py +++ b/hexrd/material/material.py @@ -223,7 +223,7 @@ def __init__( self._newUnitcell() if not hasattr(self, 'v0'): - self.update_v0() + self.reset_v0() self._newPdata() self.update_structure_factor() @@ -341,7 +341,7 @@ def _newPdata(self): self.set_default_exclusions() - def update_v0(self): + def reset_v0(self): self.v0 = self.vol def set_default_exclusions(self): From 23d0db7094545cc1f9bcbddcf6b171681a2c3783 Mon Sep 17 00:00:00 2001 From: Patrick Avery Date: Wed, 3 Jan 2024 05:36:59 -0600 Subject: [PATCH 26/26] Create aliases to import paths to fix old scripts These import path aliases should keep old scripts from breaking. For example, you can now still do all of the following: ```python from hexrd.crystallography import PlaneData from hexrd import crystallography crystallography.PlaneData ``` And using the new import paths, you can do the following: ```python from hexrd.material.crystallography import PlaneData from hexrd.material import crystallography crystallography.PlaneData ``` We also added a check to verify that the import path alias is not a real file path. If it is, then an exception is raised. This will keep us from having real file paths overlapping with the alias import paths. Signed-off-by: Patrick Avery --- hexrd/__init__.py | 25 +++++++++++++++++++++++++ 1 file changed, 25 insertions(+) diff --git a/hexrd/__init__.py b/hexrd/__init__.py index f9da2f348..77ca9b19a 100644 --- a/hexrd/__init__.py +++ b/hexrd/__init__.py @@ -1,3 +1,6 @@ +import importlib +import sys + from .material import crystallography from .material import jcpds from .material import mksupport @@ -5,3 +8,25 @@ from .material import symbols from .material import symmetry from .material import unitcell + +# These are aliases for import paths, so we don't break old HEXRD scripts. +# We will verify the alias files *do not* exist, to avoid confusion. +module_aliases = { + 'hexrd.crystallography': crystallography, + 'hexrd.mksupport': mksupport, + 'hexrd.spacegroup': spacegroup, + 'hexrd.symbols': symbols, + 'hexrd.symmetry': symmetry, + 'hexrd.unitcell': unitcell, +} + +for alias, module in module_aliases.items(): + try: + file_exists = importlib.import_module(alias).__name__ == alias + except ImportError: + file_exists = False + + if file_exists: + raise Exception(f'"{alias}" is an alias path and should not exist') + + sys.modules[alias] = module