diff --git a/ImageD11/eps_sig_solver.py b/ImageD11/depreciated/eps_sig_solver.py similarity index 100% rename from ImageD11/eps_sig_solver.py rename to ImageD11/depreciated/eps_sig_solver.py diff --git a/ImageD11/eps_sig_solver_NEW.py b/ImageD11/stress.py similarity index 94% rename from ImageD11/eps_sig_solver_NEW.py rename to ImageD11/stress.py index 2837d9fc..2a3af829 100644 --- a/ImageD11/eps_sig_solver_NEW.py +++ b/ImageD11/stress.py @@ -40,14 +40,13 @@ def __init__(self, ----------- cij (float) : elastic constants of the crystal Cij_symmetry (str) : symmetry considered for the Stiffness and Compliance matrices. Should be one of the following: - 'cubic', 'trigonal_high', 'trigonal_low', 'tetragonal', 'hexagonal','orthorombic', 'monoclinic', 'triclinic' + 'cubic', 'trigonal_high', 'trigonal_low', 'tetragonal_high', 'tetragonal_low', 'hexagonal', 'orthorhombic', 'monoclinic', 'triclinic' - dzero_unitcell (array_like) : Unstrained unit cell parameters [a, b, c, alpha,beta, gamma] + unitcell (array_like) : Unstrained unit cell parameters [a, b, c, alpha,beta, gamma] UBI_list (list of 3x3 arrays) : List of real-space unit cell vectors (ubi in ImageD11). """ - assert symmetry in ['cubic', 'trigonal_high','trigonal_low', 'tetragonal', - 'hexagonal','orthorombic', 'monoclinic', 'triclinic'], 'symmetry not recognized!' + assert symmetry in Cij_symmetry.keys(), 'symmetry not recognized!' self.phase_name = name self.unitcell = unitcell @@ -80,16 +79,17 @@ def __init__(self, self.stress_unit = 'GPa' self.Cij_symmetry = Cij_symmetry[self.symmetry] self.Cij = self.form_stiffness_tensor() - self.Sij = np.linalg.inv(self.Cij) + if not np.alltrue(self.Cij == 0): + self.Sij = np.linalg.inv(self.Cij) self.UBIs = UBI_list self.F_list = None - def __str__(self): - return f"EpsSigSolver:\n phase name: {self.phase_name}\n reference unitcell: {self.unitcell}\n symmetry:" +\ - f"{self.symmetry}\n unit:{self.stress_unit}\n Stiffness:\n {self.Cij}\n n ubis: {len(self.UBIs)}" + # def __str__(self): + # return f"EpsSigSolver:\n phase name: {self.phase_name}\n reference unitcell: {self.unitcell}\n symmetry:" +\ + # f"{self.symmetry}\n unit:{self.stress_unit}\n Stiffness:\n {self.Cij}\n n ubis: {len(self.UBIs)}" # Load / save / update functions for parameters (from former eps_sig_solver.py) ######################################## @@ -138,7 +138,8 @@ def form_stiffness_tensor(self): Cij : 6x6 matrix containing the elastic components """ Cij = np.zeros((6,6)) - pattern = self.Cij_symmetry # pattern for the stiffness matrix + pattern = self.Cij_symmetry # pattern for the stiffness matrix. From Mouhat & Coudert (2014). Necessary and sufficient elastic stability conditions in various crystal systems. Physical Review + parlist = 'c11,c22,c33,c44,c55,c66,c12,c13,c14,c15,c16,c23,c24,c25,c26,c34,c35,c36,c45,c46,c56'.split(',') @@ -146,6 +147,9 @@ def form_stiffness_tensor(self): v = self.parameterobj.parameters[parname] if v is None: continue + if self.symmetry in ['hexagonal','trigonal_high','trigonal_low'] and parname == 'c66': + v = 0.5 * (self.parameterobj.parameters['c11'] - self.parameterobj.parameters['c12']) + c_ij = np.where( np.abs(pattern) == i+1, np.sign(pattern) * v, 0 ) Cij += c_ij return Cij @@ -274,7 +278,7 @@ def convert_tensors_to_vecs(self, output_format = 'voigt', debug=0): # select all strain and stress tensors list dnames = [attr for attr in dir(self) if any([attr.startswith(s) for s in ['eps','sigma']]) ] # filter out all data that begins with 'eps' or 'sigma' but are not strain or stress tensors - dnames = [d for d in dnames if any([d.endswith(s) for s in ['Lab','Ref','d']]) ] + dnames = [d for d in dnames if any([d.endswith(s) for s in ['Lab','Ref','_d']]) ] # stop if no strain / stress tensor list fond if len(dnames) == 0: @@ -581,29 +585,8 @@ def invariants(T): I3 = np.linalg.det(T) return I1, I2, I3 - - -def invariants_quantities(T): - """ - compute relevant invariant parameters of strain / stress tensor T - - Returns - -------- - P (float) : -I1/3 : hydrostatic pressure / volumetric strain - τ_oct (float) : √(2*J2/3) : octahedral shear stress / strain - """ - T_dev = deviatoric(T) - I1, I2, I3 = invariants(T) - J1, J2, J3 = invariants(T_dev) - - return -I1/3, np.sqrt(2*J2/3) - - - - - -# Cij_symmetry for stiffness tensor (copied from matscipy.elasticity : https://github.com/libAtoms/matscipy/blob/master/matscipy/elasticity.py) +# Cij_symmetry for stiffness tensor (from Mouhat & Coudert 2014) ######################## Cij_symmetry = { @@ -614,12 +597,12 @@ def invariants_quantities(T): [0, 0, 0, 0, 4, 0], [0, 0, 0, 0, 0, 4]]), - 'trigonal_high': np.array([[1, 7, 8, 9, 10, 0], - [7, 1, 8, 0,-9, 0], - [8, 8, 3, 0, 0, 0], + 'trigonal_high': np.array([[1, 7, 8, 9, 0, 0], + [7, 1, 8, -9, 0, 0], + [8, 8, 3, 0, 0, 0], [9, -9, 0, 4, 0, 0], - [10, 0, 0, 0, 4, 0], - [0, 0, 0, 0, 0, 6]]), + [0, 0, 0, 0, 4, 9], + [0, 0, 0, 0, 9, 6]]), 'trigonal_low': np.array([[1, 7, 8, 9, 10, 0 ], [7, 1, 8, -9, -10, 0 ], @@ -665,5 +648,6 @@ def invariants_quantities(T): } -Cij_symmetry['hexagonal'] = Cij_symmetry['trigonal_high'] +Cij_symmetry['hexagonal'] = Cij_symmetry['tetragonal_high'] +Cij_symmetry['tetragonal'] = Cij_symmetry['tetragonal_high'] Cij_symmetry[None] = Cij_symmetry['triclinic']