diff --git a/PRIVATE b/PRIVATE index 92b7b1314..f529a16d1 160000 --- a/PRIVATE +++ b/PRIVATE @@ -1 +1 @@ -Subproject commit 92b7b1314a9c576a20f073a230e2aaf811cb932a +Subproject commit f529a16d100434c736476944a3a696d8f95ac770 diff --git a/python/damask/__init__.py b/python/damask/__init__.py index 600f64138..9d0350337 100644 --- a/python/damask/__init__.py +++ b/python/damask/__init__.py @@ -15,12 +15,12 @@ from . import seeds # noqa from . import mechanics # noqa from . import solver # noqa from . import grid_filters # noqa -from ._lattice import Symmetry, Lattice# noqa -from ._table import Table # noqa +from . import lattice # noqa from ._rotation import Rotation # noqa +from ._orientation import Orientation # noqa +from ._table import Table # noqa from ._vtk import VTK # noqa from ._colormap import Colormap # noqa -from ._orientation import Orientation # noqa from ._config import Config # noqa from ._configmaterial import ConfigMaterial # noqa from ._geom import Geom # noqa diff --git a/python/damask/_configmaterial.py b/python/damask/_configmaterial.py index cdeab3a53..4f3eba816 100644 --- a/python/damask/_configmaterial.py +++ b/python/damask/_configmaterial.py @@ -3,8 +3,8 @@ import copy import numpy as np from . import Config -from . import Lattice from . import Rotation +from . import Orientation class ConfigMaterial(Config): """Material configuration.""" @@ -152,7 +152,7 @@ class ConfigMaterial(Config): for k,v in self['phase'].items(): if 'lattice' in v: try: - Lattice(v['lattice']) + Orientation(lattice=v['lattice']) except KeyError: s = v['lattice'] print(f"Invalid lattice: '{s}' in phase '{k}'") diff --git a/python/damask/_lattice.py b/python/damask/_lattice.py deleted file mode 100644 index 143fa50f1..000000000 --- a/python/damask/_lattice.py +++ /dev/null @@ -1,646 +0,0 @@ -import numpy as np - - -class Symmetry: - """ - Symmetry-related operations for crystal systems. - - References - ---------- - https://en.wikipedia.org/wiki/Crystal_system - - """ - - crystal_systems = [None,'orthorhombic','tetragonal','hexagonal','cubic'] - - def __init__(self, system = None): - """ - Symmetry Definition. - - Parameters - ---------- - system : {None,'orthorhombic','tetragonal','hexagonal','cubic'}, optional - Name of the crystal system. Defaults to 'None'. - - """ - if system is not None and system.lower() not in self.crystal_systems: - raise KeyError(f'Crystal system "{system}" is unknown') - - self.system = system.lower() if isinstance(system,str) else system - - - def __copy__(self): - """Copy.""" - return self.__class__(self.system) - - copy = __copy__ - - - def __repr__(self): - """Readable string.""" - return f'{self.system}' - - - def __eq__(self, other): - """ - Equal to other. - - Parameters - ---------- - other : Symmetry - Symmetry to check for equality. - - """ - return self.system == other.system - - def __neq__(self, other): - """ - Not Equal to other. - - Parameters - ---------- - other : Symmetry - Symmetry to check for inequality. - - """ - return not self.__eq__(other) - - def __cmp__(self,other): - """ - Linear ordering. - - Parameters - ---------- - other : Symmetry - Symmetry to check for for order. - - """ - myOrder = self.crystal_systems.index(self.system) - otherOrder = self.crystal_systems.index(other.system) - return (myOrder > otherOrder) - (myOrder < otherOrder) - - - @property - def symmetry_operations(self): - """Symmetry operations as quaternions.""" - if self.system == 'cubic': - sym_quats = [ - [ 1.0, 0.0, 0.0, 0.0 ], - [ 0.0, 1.0, 0.0, 0.0 ], - [ 0.0, 0.0, 1.0, 0.0 ], - [ 0.0, 0.0, 0.0, 1.0 ], - [ 0.0, 0.0, 0.5*np.sqrt(2), 0.5*np.sqrt(2) ], - [ 0.0, 0.0, 0.5*np.sqrt(2),-0.5*np.sqrt(2) ], - [ 0.0, 0.5*np.sqrt(2), 0.0, 0.5*np.sqrt(2) ], - [ 0.0, 0.5*np.sqrt(2), 0.0, -0.5*np.sqrt(2) ], - [ 0.0, 0.5*np.sqrt(2),-0.5*np.sqrt(2), 0.0 ], - [ 0.0, -0.5*np.sqrt(2),-0.5*np.sqrt(2), 0.0 ], - [ 0.5, 0.5, 0.5, 0.5 ], - [-0.5, 0.5, 0.5, 0.5 ], - [-0.5, 0.5, 0.5, -0.5 ], - [-0.5, 0.5, -0.5, 0.5 ], - [-0.5, -0.5, 0.5, 0.5 ], - [-0.5, -0.5, 0.5, -0.5 ], - [-0.5, -0.5, -0.5, 0.5 ], - [-0.5, 0.5, -0.5, -0.5 ], - [-0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], - [ 0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], - [-0.5*np.sqrt(2), 0.0, 0.5*np.sqrt(2), 0.0 ], - [-0.5*np.sqrt(2), 0.0, -0.5*np.sqrt(2), 0.0 ], - [-0.5*np.sqrt(2), 0.5*np.sqrt(2), 0.0, 0.0 ], - [-0.5*np.sqrt(2),-0.5*np.sqrt(2), 0.0, 0.0 ], - ] - elif self.system == 'hexagonal': - sym_quats = [ - [ 1.0, 0.0, 0.0, 0.0 ], - [-0.5*np.sqrt(3), 0.0, 0.0, -0.5 ], - [ 0.5, 0.0, 0.0, 0.5*np.sqrt(3) ], - [ 0.0, 0.0, 0.0, 1.0 ], - [-0.5, 0.0, 0.0, 0.5*np.sqrt(3) ], - [-0.5*np.sqrt(3), 0.0, 0.0, 0.5 ], - [ 0.0, 1.0, 0.0, 0.0 ], - [ 0.0, -0.5*np.sqrt(3), 0.5, 0.0 ], - [ 0.0, 0.5, -0.5*np.sqrt(3), 0.0 ], - [ 0.0, 0.0, 1.0, 0.0 ], - [ 0.0, -0.5, -0.5*np.sqrt(3), 0.0 ], - [ 0.0, 0.5*np.sqrt(3), 0.5, 0.0 ], - ] - elif self.system == 'tetragonal': - sym_quats = [ - [ 1.0, 0.0, 0.0, 0.0 ], - [ 0.0, 1.0, 0.0, 0.0 ], - [ 0.0, 0.0, 1.0, 0.0 ], - [ 0.0, 0.0, 0.0, 1.0 ], - [ 0.0, 0.5*np.sqrt(2), 0.5*np.sqrt(2), 0.0 ], - [ 0.0, -0.5*np.sqrt(2), 0.5*np.sqrt(2), 0.0 ], - [ 0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], - [-0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], - ] - elif self.system == 'orthorhombic': - sym_quats = [ - [ 1.0,0.0,0.0,0.0 ], - [ 0.0,1.0,0.0,0.0 ], - [ 0.0,0.0,1.0,0.0 ], - [ 0.0,0.0,0.0,1.0 ], - ] - else: - sym_quats = [ - [ 1.0,0.0,0.0,0.0 ], - ] - return np.array(sym_quats) - - - def in_FZ(self,rho): - """ - Check whether given Rodrigues-Frank vector falls into fundamental zone. - - Fundamental zone in Rodrigues space is point symmetric around origin. - """ - if(rho.shape[-1] != 3): - raise ValueError('Input is not a Rodrigues-Frank vector field.') - - rho_abs = np.abs(rho) - - with np.errstate(invalid='ignore'): - # using '*'/prod for 'and' - if self.system == 'cubic': - return np.where(np.prod(np.sqrt(2)-1. >= rho_abs,axis=-1) * - (1. >= np.sum(rho_abs,axis=-1)),True,False) - elif self.system == 'hexagonal': - return np.where(np.prod(1. >= rho_abs,axis=-1) * - (2. >= np.sqrt(3)*rho_abs[...,0] + rho_abs[...,1]) * - (2. >= np.sqrt(3)*rho_abs[...,1] + rho_abs[...,0]) * - (2. >= np.sqrt(3) + rho_abs[...,2]),True,False) - elif self.system == 'tetragonal': - return np.where(np.prod(1. >= rho_abs[...,:2],axis=-1) * - (np.sqrt(2) >= rho_abs[...,0] + rho_abs[...,1]) * - (np.sqrt(2) >= rho_abs[...,2] + 1.),True,False) - elif self.system == 'orthorhombic': - return np.where(np.prod(1. >= rho_abs,axis=-1),True,False) - else: - return np.where(np.all(np.isfinite(rho_abs),axis=-1),True,False) - - - def in_disorientation_SST(self,rho): - """ - Check whether given Rodrigues-Frank vector (of misorientation) falls into standard stereographic triangle. - - References - ---------- - A. Heinz and P. Neumann, Acta Crystallographica Section A 47:780-789, 1991 - https://doi.org/10.1107/S0108767391006864 - - """ - if(rho.shape[-1] != 3): - raise ValueError('Input is not a Rodrigues-Frank vector field.') - - with np.errstate(invalid='ignore'): - # using '*' for 'and' - if self.system == 'cubic': - return np.where((rho[...,0] >= rho[...,1]) * \ - (rho[...,1] >= rho[...,2]) * \ - (rho[...,2] >= 0),True,False) - elif self.system == 'hexagonal': - return np.where((rho[...,0] >= rho[...,1]*np.sqrt(3)) * \ - (rho[...,1] >= 0) * \ - (rho[...,2] >= 0),True,False) - elif self.system == 'tetragonal': - return np.where((rho[...,0] >= rho[...,1]) * \ - (rho[...,1] >= 0) * \ - (rho[...,2] >= 0),True,False) - elif self.system == 'orthorhombic': - return np.where((rho[...,0] >= 0) * \ - (rho[...,1] >= 0) * \ - (rho[...,2] >= 0),True,False) - else: - return np.ones_like(rho[...,0],dtype=bool) - - - #ToDo: IPF color in separate function - def in_SST(self,vector,proper=False,color=False): - """ - Check whether given vector falls into standard stereographic triangle of own symmetry. - - proper considers only vectors with z >= 0, hence uses two neighboring SSTs. - Return inverse pole figure color if requested. - Bases are computed from - - >>> basis = {'cubic' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red - ... [1.,0.,1.]/np.sqrt(2.), # direction of green - ... [1.,1.,1.]/np.sqrt(3.)]).T), # direction of blue - ... 'hexagonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red - ... [1.,0.,0.], # direction of green - ... [np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T), # direction of blue - ... 'tetragonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red - ... [1.,0.,0.], # direction of green - ... [1.,1.,0.]/np.sqrt(2.)]).T), # direction of blue - ... 'orthorhombic': np.linalg.inv(np.array([[0.,0.,1.], # direction of red - ... [1.,0.,0.], # direction of green - ... [0.,1.,0.]]).T), # direction of blue - ... } - - """ - if(vector.shape[-1] != 3): - raise ValueError('Input is not a 3D vector field.') - - if self.system == 'cubic': - basis = {'improper':np.array([ [-1. , 0. , 1. ], - [ np.sqrt(2.) , -np.sqrt(2.) , 0. ], - [ 0. , np.sqrt(3.) , 0. ] ]), - 'proper':np.array([ [ 0. , -1. , 1. ], - [-np.sqrt(2.) , np.sqrt(2.) , 0. ], - [ np.sqrt(3.) , 0. , 0. ] ]), - } - elif self.system == 'hexagonal': - basis = {'improper':np.array([ [ 0. , 0. , 1. ], - [ 1. , -np.sqrt(3.) , 0. ], - [ 0. , 2. , 0. ] ]), - 'proper':np.array([ [ 0. , 0. , 1. ], - [-1. , np.sqrt(3.) , 0. ], - [ np.sqrt(3.) , -1. , 0. ] ]), - } - elif self.system == 'tetragonal': - basis = {'improper':np.array([ [ 0. , 0. , 1. ], - [ 1. , -1. , 0. ], - [ 0. , np.sqrt(2.) , 0. ] ]), - 'proper':np.array([ [ 0. , 0. , 1. ], - [-1. , 1. , 0. ], - [ np.sqrt(2.) , 0. , 0. ] ]), - } - elif self.system == 'orthorhombic': - basis = {'improper':np.array([ [ 0., 0., 1.], - [ 1., 0., 0.], - [ 0., 1., 0.] ]), - 'proper':np.array([ [ 0., 0., 1.], - [-1., 0., 0.], - [ 0., 1., 0.] ]), - } - else: # direct exit for unspecified symmetry - if color: - return (np.ones_like(vector[...,0],bool),np.zeros_like(vector)) - else: - return np.ones_like(vector[...,0],bool) - - - b_i = np.broadcast_to(basis['improper'],vector.shape+(3,)) - if proper: - b_p = np.broadcast_to(basis['proper'], vector.shape+(3,)) - improper = np.all(np.around(np.einsum('...ji,...i',b_i,vector),12)>=0.0,axis=-1,keepdims=True) - theComponents = np.where(np.broadcast_to(improper,vector.shape), - np.around(np.einsum('...ji,...i',b_i,vector),12), - np.around(np.einsum('...ji,...i',b_p,vector),12)) - else: - vector_ = np.block([vector[...,0:2],np.abs(vector[...,2:3])]) # z component projects identical - theComponents = np.around(np.einsum('...ji,...i',b_i,vector_),12) - - in_SST = np.all(theComponents >= 0.0,axis=-1) - - if color: # have to return color array - with np.errstate(invalid='ignore',divide='ignore'): - rgb = (theComponents/np.linalg.norm(theComponents,axis=-1,keepdims=True))**0.5 # smoothen color ramps - rgb = np.minimum(1.,rgb) # limit to maximum intensity - rgb /= np.max(rgb,axis=-1,keepdims=True) # normalize to (HS)V = 1 - rgb[np.broadcast_to(~in_SST.reshape(vector[...,0].shape+(1,)),vector.shape)] = 0.0 - return (in_SST,rgb) - else: - return in_SST - - -# ****************************************************************************************** -class Lattice: # ToDo: Make a subclass of Symmetry! - """ - Bravais lattice. - - This contains only a mapping from Bravais lattice to symmetry - and orientation relationships. It could include twin and slip systems. - - References - ---------- - https://en.wikipedia.org/wiki/Bravais_lattice - - """ - - lattices = { - 'iso': {'system':None}, - 'triclinic':{'system':None}, - 'bct': {'system':'tetragonal'}, - 'hex': {'system':'hexagonal'}, - 'fcc': {'system':'cubic','c/a':1.0}, - 'bcc': {'system':'cubic','c/a':1.0}, - } - - - def __init__(self,lattice,c_over_a=None): - """ - New lattice of given type. - - Parameters - ---------- - lattice : str - Bravais lattice. - - """ - self.lattice = lattice - self.symmetry = Symmetry(self.lattices[lattice]['system']) - - # transition to subclass - self.system = self.symmetry.system - self.in_SST = self.symmetry.in_SST - self.in_FZ = self.symmetry.in_FZ - self.in_disorientation_SST = self.symmetry.in_disorientation_SST - - def __repr__(self): - """Report basic lattice information.""" - return f'Bravais lattice {self.lattice} ({self.symmetry} crystal system)' - - - # Kurdjomov--Sachs orientation relationship for fcc <-> bcc transformation - # from S. Morito et al., Journal of Alloys and Compounds 577:s587-s592, 2013 - # also see K. Kitahara et al., Acta Materialia 54:1279-1288, 2006 - _KS = {'mapping':{'fcc':0,'bcc':1}, - 'planes': np.array([ - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, -1],[ 0, 1, 1]], - [[ 1, 1, -1],[ 0, 1, 1]], - [[ 1, 1, -1],[ 0, 1, 1]], - [[ 1, 1, -1],[ 0, 1, 1]], - [[ 1, 1, -1],[ 0, 1, 1]], - [[ 1, 1, -1],[ 0, 1, 1]]],dtype='float'), - 'directions': np.array([ - [[ -1, 0, 1],[ -1, -1, 1]], - [[ -1, 0, 1],[ -1, 1, -1]], - [[ 0, 1, -1],[ -1, -1, 1]], - [[ 0, 1, -1],[ -1, 1, -1]], - [[ 1, -1, 0],[ -1, -1, 1]], - [[ 1, -1, 0],[ -1, 1, -1]], - [[ 1, 0, -1],[ -1, -1, 1]], - [[ 1, 0, -1],[ -1, 1, -1]], - [[ -1, -1, 0],[ -1, -1, 1]], - [[ -1, -1, 0],[ -1, 1, -1]], - [[ 0, 1, 1],[ -1, -1, 1]], - [[ 0, 1, 1],[ -1, 1, -1]], - [[ 0, -1, 1],[ -1, -1, 1]], - [[ 0, -1, 1],[ -1, 1, -1]], - [[ -1, 0, -1],[ -1, -1, 1]], - [[ -1, 0, -1],[ -1, 1, -1]], - [[ 1, 1, 0],[ -1, -1, 1]], - [[ 1, 1, 0],[ -1, 1, -1]], - [[ -1, 1, 0],[ -1, -1, 1]], - [[ -1, 1, 0],[ -1, 1, -1]], - [[ 0, -1, -1],[ -1, -1, 1]], - [[ 0, -1, -1],[ -1, 1, -1]], - [[ 1, 0, 1],[ -1, -1, 1]], - [[ 1, 0, 1],[ -1, 1, -1]]],dtype='float')} - - # Greninger--Troiano orientation relationship for fcc <-> bcc transformation - # from Y. He et al., Journal of Applied Crystallography 39:72-81, 2006 - _GT = {'mapping':{'fcc':0,'bcc':1}, - 'planes': np.array([ - [[ 1, 1, 1],[ 1, 0, 1]], - [[ 1, 1, 1],[ 1, 1, 0]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ -1, -1, 1],[ -1, 0, 1]], - [[ -1, -1, 1],[ -1, -1, 0]], - [[ -1, -1, 1],[ 0, -1, 1]], - [[ -1, 1, 1],[ -1, 0, 1]], - [[ -1, 1, 1],[ -1, 1, 0]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 1, 0, 1]], - [[ 1, -1, 1],[ 1, -1, 0]], - [[ 1, -1, 1],[ 0, -1, 1]], - [[ 1, 1, 1],[ 1, 1, 0]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 1, 0, 1]], - [[ -1, -1, 1],[ -1, -1, 0]], - [[ -1, -1, 1],[ 0, -1, 1]], - [[ -1, -1, 1],[ -1, 0, 1]], - [[ -1, 1, 1],[ -1, 1, 0]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ -1, 0, 1]], - [[ 1, -1, 1],[ 1, -1, 0]], - [[ 1, -1, 1],[ 0, -1, 1]], - [[ 1, -1, 1],[ 1, 0, 1]]],dtype='float'), - 'directions': np.array([ - [[ -5,-12, 17],[-17, -7, 17]], - [[ 17, -5,-12],[ 17,-17, -7]], - [[-12, 17, -5],[ -7, 17,-17]], - [[ 5, 12, 17],[ 17, 7, 17]], - [[-17, 5,-12],[-17, 17, -7]], - [[ 12,-17, -5],[ 7,-17,-17]], - [[ -5, 12,-17],[-17, 7,-17]], - [[ 17, 5, 12],[ 17, 17, 7]], - [[-12,-17, 5],[ -7,-17, 17]], - [[ 5,-12,-17],[ 17, -7,-17]], - [[-17, -5, 12],[-17,-17, 7]], - [[ 12, 17, 5],[ 7, 17, 17]], - [[ -5, 17,-12],[-17, 17, -7]], - [[-12, -5, 17],[ -7,-17, 17]], - [[ 17,-12, -5],[ 17, -7,-17]], - [[ 5,-17,-12],[ 17,-17, -7]], - [[ 12, 5, 17],[ 7, 17, 17]], - [[-17, 12, -5],[-17, 7,-17]], - [[ -5,-17, 12],[-17,-17, 7]], - [[-12, 5,-17],[ -7, 17,-17]], - [[ 17, 12, 5],[ 17, 7, 17]], - [[ 5, 17, 12],[ 17, 17, 7]], - [[ 12, -5,-17],[ 7,-17,-17]], - [[-17,-12, 5],[-17,-7, 17]]],dtype='float')} - - # Greninger--Troiano' orientation relationship for fcc <-> bcc transformation - # from Y. He et al., Journal of Applied Crystallography 39:72-81, 2006 - _GTprime = {'mapping':{'fcc':0,'bcc':1}, - 'planes': np.array([ - [[ 7, 17, 17],[ 12, 5, 17]], - [[ 17, 7, 17],[ 17, 12, 5]], - [[ 17, 17, 7],[ 5, 17, 12]], - [[ -7,-17, 17],[-12, -5, 17]], - [[-17, -7, 17],[-17,-12, 5]], - [[-17,-17, 7],[ -5,-17, 12]], - [[ 7,-17,-17],[ 12, -5,-17]], - [[ 17, -7,-17],[ 17,-12, -5]], - [[ 17,-17, -7],[ 5,-17,-12]], - [[ -7, 17,-17],[-12, 5,-17]], - [[-17, 7,-17],[-17, 12, -5]], - [[-17, 17, -7],[ -5, 17,-12]], - [[ 7, 17, 17],[ 12, 17, 5]], - [[ 17, 7, 17],[ 5, 12, 17]], - [[ 17, 17, 7],[ 17, 5, 12]], - [[ -7,-17, 17],[-12,-17, 5]], - [[-17, -7, 17],[ -5,-12, 17]], - [[-17,-17, 7],[-17, -5, 12]], - [[ 7,-17,-17],[ 12,-17, -5]], - [[ 17, -7,-17],[ 5, -12,-17]], - [[ 17,-17, -7],[ 17, -5,-12]], - [[ -7, 17,-17],[-12, 17, -5]], - [[-17, 7,-17],[ -5, 12,-17]], - [[-17, 17, -7],[-17, 5,-12]]],dtype='float'), - 'directions': np.array([ - [[ 0, 1, -1],[ 1, 1, -1]], - [[ -1, 0, 1],[ -1, 1, 1]], - [[ 1, -1, 0],[ 1, -1, 1]], - [[ 0, -1, -1],[ -1, -1, -1]], - [[ 1, 0, 1],[ 1, -1, 1]], - [[ 1, -1, 0],[ 1, -1, -1]], - [[ 0, 1, -1],[ -1, 1, -1]], - [[ 1, 0, 1],[ 1, 1, 1]], - [[ -1, -1, 0],[ -1, -1, 1]], - [[ 0, -1, -1],[ 1, -1, -1]], - [[ -1, 0, 1],[ -1, -1, 1]], - [[ -1, -1, 0],[ -1, -1, -1]], - [[ 0, -1, 1],[ 1, -1, 1]], - [[ 1, 0, -1],[ 1, 1, -1]], - [[ -1, 1, 0],[ -1, 1, 1]], - [[ 0, 1, 1],[ -1, 1, 1]], - [[ -1, 0, -1],[ -1, -1, -1]], - [[ -1, 1, 0],[ -1, 1, -1]], - [[ 0, -1, 1],[ -1, -1, 1]], - [[ -1, 0, -1],[ -1, 1, -1]], - [[ 1, 1, 0],[ 1, 1, 1]], - [[ 0, 1, 1],[ 1, 1, 1]], - [[ 1, 0, -1],[ 1, -1, -1]], - [[ 1, 1, 0],[ 1, 1, -1]]],dtype='float')} - - # Nishiyama--Wassermann orientation relationship for fcc <-> bcc transformation - # from H. Kitahara et al., Materials Characterization 54:378-386, 2005 - _NW = {'mapping':{'fcc':0,'bcc':1}, - 'planes': np.array([ - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ 1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ -1, 1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ 1, -1, 1],[ 0, 1, 1]], - [[ -1, -1, 1],[ 0, 1, 1]], - [[ -1, -1, 1],[ 0, 1, 1]], - [[ -1, -1, 1],[ 0, 1, 1]]],dtype='float'), - 'directions': np.array([ - [[ 2, -1, -1],[ 0, -1, 1]], - [[ -1, 2, -1],[ 0, -1, 1]], - [[ -1, -1, 2],[ 0, -1, 1]], - [[ -2, -1, -1],[ 0, -1, 1]], - [[ 1, 2, -1],[ 0, -1, 1]], - [[ 1, -1, 2],[ 0, -1, 1]], - [[ 2, 1, -1],[ 0, -1, 1]], - [[ -1, -2, -1],[ 0, -1, 1]], - [[ -1, 1, 2],[ 0, -1, 1]], - [[ 2, -1, 1],[ 0, -1, 1]], #It is wrong in the paper, but matrix is correct - [[ -1, 2, 1],[ 0, -1, 1]], - [[ -1, -1, -2],[ 0, -1, 1]]],dtype='float')} - - # Pitsch orientation relationship for fcc <-> bcc transformation - # from Y. He et al., Acta Materialia 53:1179-1190, 2005 - _Pitsch = {'mapping':{'fcc':0,'bcc':1}, - 'planes': np.array([ - [[ 0, 1, 0],[ -1, 0, 1]], - [[ 0, 0, 1],[ 1, -1, 0]], - [[ 1, 0, 0],[ 0, 1, -1]], - [[ 1, 0, 0],[ 0, -1, -1]], - [[ 0, 1, 0],[ -1, 0, -1]], - [[ 0, 0, 1],[ -1, -1, 0]], - [[ 0, 1, 0],[ -1, 0, -1]], - [[ 0, 0, 1],[ -1, -1, 0]], - [[ 1, 0, 0],[ 0, -1, -1]], - [[ 1, 0, 0],[ 0, -1, 1]], - [[ 0, 1, 0],[ 1, 0, -1]], - [[ 0, 0, 1],[ -1, 1, 0]]],dtype='float'), - 'directions': np.array([ - [[ 1, 0, 1],[ 1, -1, 1]], - [[ 1, 1, 0],[ 1, 1, -1]], - [[ 0, 1, 1],[ -1, 1, 1]], - [[ 0, 1, -1],[ -1, 1, -1]], - [[ -1, 0, 1],[ -1, -1, 1]], - [[ 1, -1, 0],[ 1, -1, -1]], - [[ 1, 0, -1],[ 1, -1, -1]], - [[ -1, 1, 0],[ -1, 1, -1]], - [[ 0, -1, 1],[ -1, -1, 1]], - [[ 0, 1, 1],[ -1, 1, 1]], - [[ 1, 0, 1],[ 1, -1, 1]], - [[ 1, 1, 0],[ 1, 1, -1]]],dtype='float')} - - # Bain orientation relationship for fcc <-> bcc transformation - # from Y. He et al., Journal of Applied Crystallography 39:72-81, 2006 - _Bain = {'mapping':{'fcc':0,'bcc':1}, - 'planes': np.array([ - [[ 1, 0, 0],[ 1, 0, 0]], - [[ 0, 1, 0],[ 0, 1, 0]], - [[ 0, 0, 1],[ 0, 0, 1]]],dtype='float'), - 'directions': np.array([ - [[ 0, 1, 0],[ 0, 1, 1]], - [[ 0, 0, 1],[ 1, 0, 1]], - [[ 1, 0, 0],[ 1, 1, 0]]],dtype='float')} - - - def relation_operations(self,model): - """ - Crystallographic orientation relationships for phase transformations. - - References - ---------- - S. Morito et al., Journal of Alloys and Compounds 577:s587-s592, 2013 - https://doi.org/10.1016/j.jallcom.2012.02.004 - - K. Kitahara et al., Acta Materialia 54(5):1279-1288, 2006 - https://doi.org/10.1016/j.actamat.2005.11.001 - - Y. He et al., Journal of Applied Crystallography 39:72-81, 2006 - https://doi.org/10.1107/S0021889805038276 - - H. Kitahara et al., Materials Characterization 54(4-5):378-386, 2005 - https://doi.org/10.1016/j.matchar.2004.12.015 - - Y. He et al., Acta Materialia 53(4):1179-1190, 2005 - https://doi.org/10.1016/j.actamat.2004.11.021 - - """ - models={'KS':self._KS, 'GT':self._GT, 'GT_prime':self._GTprime, - 'NW':self._NW, 'Pitsch': self._Pitsch, 'Bain':self._Bain} - try: - relationship = models[model] - except KeyError : - raise KeyError(f'Orientation relationship "{model}" is unknown') - - if self.lattice not in relationship['mapping']: - raise ValueError(f'Relationship "{model}" not supported for lattice "{self.lattice}"') - - r = {'lattice':Lattice((set(relationship['mapping'])-{self.lattice}).pop()), # target lattice - 'rotations':[] } - - myPlane_id = relationship['mapping'][self.lattice] - otherPlane_id = (myPlane_id+1)%2 - myDir_id = myPlane_id +2 - otherDir_id = otherPlane_id +2 - - for miller in np.hstack((relationship['planes'],relationship['directions'])): - myPlane = miller[myPlane_id]/ np.linalg.norm(miller[myPlane_id]) - myDir = miller[myDir_id]/ np.linalg.norm(miller[myDir_id]) - myMatrix = np.array([myDir,np.cross(myPlane,myDir),myPlane]) - - otherPlane = miller[otherPlane_id]/ np.linalg.norm(miller[otherPlane_id]) - otherDir = miller[otherDir_id]/ np.linalg.norm(miller[otherDir_id]) - otherMatrix = np.array([otherDir,np.cross(otherPlane,otherDir),otherPlane]) - - r['rotations'].append(np.dot(otherMatrix.T,myMatrix)) - - r['rotations'] = np.array(r['rotations']) - - return r diff --git a/python/damask/_orientation.py b/python/damask/_orientation.py index 0bb0e1bc8..e0a00e0a0 100644 --- a/python/damask/_orientation.py +++ b/python/damask/_orientation.py @@ -1,202 +1,1234 @@ import numpy as np -from . import Lattice from . import Rotation +from . import util +from . import mechanics -class Orientation: # ToDo: make subclass of lattice and Rotation? +__parameter_doc__ = \ + """lattice : str + Either a crystal family out of [triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic] + or a Bravais lattice out of [aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF]. + When specifying a Bravais lattice, additional lattice parameters might be required: + a : float, optional + Length of lattice parameter "a". + b : float, optional + Length of lattice parameter "b". + c : float, optional + Length of lattice parameter "c". + alpha : float, optional + Angle between b and c lattice basis. + beta : float, optional + Angle between c and a lattice basis. + gamma : float, optional + Angle between a and b lattice basis. + degrees : bool, optional + Angles are given in degrees. Defaults to False. + + """ + + +def extend_docstring(): + """Decorator: Append Orientation parameter documentation to function's docstring.""" + def _decorator(func): + func.__doc__ += __parameter_doc__ + return func + return _decorator + + +def extended_docstring(f): + """Decorator: Combine Orientation parameter documentation with another function's docstring.""" + def _decorator(func): + func.__doc__ = f.__doc__ + __parameter_doc__ + return func + return _decorator + + +class Orientation(Rotation): """ - Crystallographic orientation. + Representation of crystallographic orientation as combination of rotation and either crystal family or Bravais lattice. - A crystallographic orientation contains a rotation and a lattice. + The crystal family is one of Orientation.crystal_families: + + - triclinic + - monoclinic + - orthorhombic + - tetragonal + - hexagonal + - cubic + + and enables symmetry-related operations such as + "equivalent", "reduced", "disorientation", "IPF_color", or "to_SST". + + The Bravais lattice is one of Orientation.lattice_symmetries: + + - aP : triclinic primitive + - mP : monoclinic primitive + - mS : ... base-centered + - oP : orthorhombic primitive + - oS : ... base-centered + - oI : ... body-centered + - oF : ... face-centered + - tP : tetragonal primitive + - tI : ... body-centered + - hP : hexagonal primitive + - cP : cubic primitive + - cI : ... body-centered + - cF : ... face-centered + + and inherits the corresponding crystal family. + Specifying a Bravais lattice, compared to just the crystal family, + extends the functionality of Orientation objects to include operations such as + "Schmid", "related", or "to_pole" that require a lattice type and its parameters. + + Examples + -------- + An array of 3 x 5 random orientations reduced to the fundamental zone of tetragonal symmetry: + >>> damask.Orientation.from_random(shape=(3,5),lattice='tetragonal').reduced + + Disorientation between two specific orientations of hexagonal symmetry: + >>> a = damask.Orientation.from_Eulers(phi=[123,32,21],degrees=True,lattice='hexagonal') + >>> b = damask.Orientation.from_Eulers(phi=[104,11,87],degrees=True,lattice='hexagonal') + >>> a.disorientation(b) + + Inverse pole figure color of the e_3 direction for a crystal in "Cube" orientation with cubic symmetry: + >>> o = damask.Orientation(lattice='cubic') + >>> o.IPF_color(o.to_SST(np.array([0,0,1]))) + + Schmid matrix (in lab frame) of slip systems of a face-centered cubic crystal in "Goss" orientation: + >>> damask.Orientation.from_Eulers(phi=[0,45,0],degrees=True,lattice='cF').Schmid('slip') """ - __slots__ = ['rotation','lattice'] + crystal_families = ['triclinic', + 'monoclinic', + 'orthorhombic', + 'tetragonal', + 'hexagonal', + 'cubic'] - def __repr__(self): - """Report lattice type and orientation.""" - return self.lattice.__repr__()+'\n'+self.rotation.__repr__() + lattice_symmetries = { + 'aP': 'triclinic', - def __init__(self, rotation, lattice): + 'mP': 'monoclinic', + 'mS': 'monoclinic', + + 'oP': 'orthorhombic', + 'oS': 'orthorhombic', + 'oI': 'orthorhombic', + 'oF': 'orthorhombic', + + 'tP': 'tetragonal', + 'tI': 'tetragonal', + + 'hP': 'hexagonal', + + 'cP': 'cubic', + 'cI': 'cubic', + 'cF': 'cubic', + } + + + @extend_docstring() + def __init__(self, + rotation = None, + lattice = None, + a = None,b = None,c = None, + alpha = None,beta = None,gamma = None, + degrees = False, + **kwargs): """ - New orientation from rotation and lattice. + Initialize orientation object. Parameters ---------- - rotation : Rotation - Rotation specifying the lattice orientation. - lattice : Lattice - Lattice type of the crystal. + rotation : list, numpy.ndarray, Rotation, optional + Unit quaternion in positive real hemisphere. + Use .from_quaternion to perform a sanity check. + Defaults to no rotation. """ - if isinstance(lattice, Lattice): + from damask.lattice import kinematics + + Rotation.__init__(self) if rotation is None else Rotation.__init__(self,rotation=rotation) + + if ( lattice is not None + and lattice not in self.lattice_symmetries + and lattice not in self.crystal_families): + raise KeyError(f'Lattice "{lattice}" is unknown') + + self.family = None + self.lattice = None + self.a = None + self.b = None + self.c = None + self.alpha = None + self.beta = None + self.gamma = None + self.kinematics = None + + if lattice in self.lattice_symmetries: + self.family = self.lattice_symmetries[lattice] self.lattice = lattice - else: - self.lattice = Lattice(lattice) # assume string + self.a = 1 if a is None else a + self.b = b + self.c = c + self.alpha = (np.radians(alpha) if degrees else alpha) if alpha is not None else None + self.beta = (np.radians(beta) if degrees else beta) if beta is not None else None + self.gamma = (np.radians(gamma) if degrees else gamma) if gamma is not None else None - if isinstance(rotation, Rotation): - self.rotation = rotation - else: - self.rotation = Rotation.from_quaternion(rotation) # assume quaternion + self.a = float(self.a) if self.a is not None else \ + (self.b / self.ratio['b'] if self.b is not None and self.ratio['b'] is not None else + self.c / self.ratio['c'] if self.c is not None and self.ratio['c'] is not None else None) + self.b = float(self.b) if self.b is not None else \ + (self.a * self.ratio['b'] if self.a is not None and self.ratio['b'] is not None else + self.c / self.ratio['c'] * self.ratio['b'] + if self.c is not None and self.ratio['b'] is not None and self.ratio['c'] is not None else None) + self.c = float(self.c) if self.c is not None else \ + (self.a * self.ratio['c'] if self.a is not None and self.ratio['c'] is not None else + self.b / self.ratio['b'] * self.ratio['c'] + if self.c is not None and self.ratio['b'] is not None and self.ratio['c'] is not None else None) + self.alpha = self.alpha if self.alpha is not None else self.immutable['alpha'] if 'alpha' in self.immutable else None + self.beta = self.beta if self.beta is not None else self.immutable['beta'] if 'beta' in self.immutable else None + self.gamma = self.gamma if self.gamma is not None else self.immutable['gamma'] if 'gamma' in self.immutable else None - def __getitem__(self,item): - """Iterate over leading/leftmost dimension of Orientation array.""" - return self.__class__(self.rotation[item],self.lattice) + if \ + (self.a is None) \ + or (self.b is None or ('b' in self.immutable and self.b != self.immutable['b'] * self.a)) \ + or (self.c is None or ('c' in self.immutable and self.c != self.immutable['c'] * self.b)) \ + or (self.alpha is None or ('alpha' in self.immutable and self.alpha != self.immutable['alpha'])) \ + or (self.beta is None or ( 'beta' in self.immutable and self.beta != self.immutable['beta'])) \ + or (self.gamma is None or ('gamma' in self.immutable and self.gamma != self.immutable['gamma'])): + raise ValueError (f'Incompatible parameters {self.parameters} for crystal family {self.family}') + + if np.any(np.array([self.alpha,self.beta,self.gamma]) <= 0): + raise ValueError ('Lattice angles must be positive') + if np.any([np.roll([self.alpha,self.beta,self.gamma],r)[0] + > np.sum(np.roll([self.alpha,self.beta,self.gamma],r)[1:]) for r in range(3)]): + raise ValueError ('Each lattice angle must be less than sum of others') + + if self.lattice in kinematics: + master = kinematics[self.lattice] + self.kinematics = {} + for m in master: + self.kinematics[m] = {'direction':master[m][:,0:3],'plane':master[m][:,3:6]} \ + if master[m].shape[-1] == 6 else \ + {'direction':self.Bravais_to_Miller(uvtw=master[m][:,0:4]), + 'plane': self.Bravais_to_Miller(hkil=master[m][:,4:8])} + elif lattice in self.crystal_families: + self.family = lattice - # ToDo: Discuss vectorization/calling signature - def disorientation(self, - other, - SST = True, - symmetries = False): + def __repr__(self): + """Represent.""" + return '\n'.join(([] if self.lattice is None else [f'Bravais lattice {self.lattice}']) + + ([] if self.family is None else [f'Crystal family {self.family}']) + + [super().__repr__()]) + + + def __copy__(self,**kwargs): + """Copy.""" + return self.__class__(rotation=kwargs['rotation'] if 'rotation' in kwargs else self.quaternion, + lattice =kwargs['lattice'] if 'lattice' in kwargs else self.lattice + if self.lattice is not None else self.family, + a =kwargs['a'] if 'a' in kwargs else self.a, + b =kwargs['b'] if 'b' in kwargs else self.b, + c =kwargs['c'] if 'c' in kwargs else self.c, + alpha =kwargs['alpha'] if 'alpha' in kwargs else self.alpha, + beta =kwargs['beta'] if 'beta' in kwargs else self.beta, + gamma =kwargs['gamma'] if 'gamma' in kwargs else self.gamma, + degrees =kwargs['degrees'] if 'degrees' in kwargs else None, + ) + + copy = __copy__ + + + def __eq__(self,other): """ - Disorientation between myself and given other orientation. + Equal to other. - Rotation axis falls into SST if SST == True. - - Currently requires same symmetry for both orientations. - Look into A. Heinz and P. Neumann 1991 for cases with differing sym. + Parameters + ---------- + other : Orientation + Orientation to check for equality. """ - if self.lattice.symmetry != other.lattice.symmetry: - raise NotImplementedError('disorientation between different symmetry classes not supported yet.') + return super().__eq__(other) \ + and self.family == other.family \ + and self.lattice == other.lattice \ + and self.parameters == other.parameters - mySymEqs = self.equivalent if SST else self.equivalent[0] #ToDo: This is just me! # take all or only first sym operation - otherSymEqs = other.equivalent - for i,sA in enumerate(mySymEqs): - aInv = sA.rotation.inversed() - for j,sB in enumerate(otherSymEqs): - b = sB.rotation - r = b*aInv - for k in range(2): - r.inverse() - breaker = self.lattice.in_FZ(r.as_Rodrigues(vector=True)) \ - and (not SST or other.lattice.in_disorientation_SST(r.as_Rodrigues(vector=True))) - if breaker: break - if breaker: break - if breaker: break + def __matmul__(self,other): + """ + Rotation of vector, second or fourth order tensor, or rotation object. + + Parameters + ---------- + other : numpy.ndarray, Rotation, or Orientation + Vector, second or fourth order tensor, or rotation object that is rotated. + + Returns + ------- + other_rot : numpy.ndarray or Rotation + Rotated vector, second or fourth order tensor, or rotation object. + + """ + return self.copy(rotation=Rotation.__matmul__(self,Rotation(other.quaternion))) \ + if isinstance(other,self.__class__) else \ + Rotation.__matmul__(self,other) + + + @classmethod + @extended_docstring(Rotation.from_random) + def from_random(cls,**kwargs): + return cls(rotation=Rotation.from_random(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_quaternion) + def from_quaternion(cls,**kwargs): + return cls(rotation=Rotation.from_quaternion(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_Eulers) + def from_Eulers(cls,**kwargs): + return cls(rotation=Rotation.from_Eulers(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_axis_angle) + def from_axis_angle(cls,**kwargs): + return cls(rotation=Rotation.from_axis_angle(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_basis) + def from_basis(cls,**kwargs): + return cls(rotation=Rotation.from_basis(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_matrix) + def from_matrix(cls,**kwargs): + return cls(rotation=Rotation.from_matrix(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_Rodrigues) + def from_Rodrigues(cls,**kwargs): + return cls(rotation=Rotation.from_Rodrigues(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_homochoric) + def from_homochoric(cls,**kwargs): + return cls(rotation=Rotation.from_homochoric(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_cubochoric) + def from_cubochoric(cls,**kwargs): + return cls(rotation=Rotation.from_cubochoric(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_spherical_component) + def from_spherical_component(cls,**kwargs): + return cls(rotation=Rotation.from_spherical_component(**kwargs),**kwargs) + + + @classmethod + @extended_docstring(Rotation.from_fiber_component) + def from_fiber_component(cls,**kwargs): + return cls(rotation=Rotation.from_fiber_component(**kwargs),**kwargs) + + + @classmethod + @extend_docstring() + def from_directions(cls,uvw,hkl,**kwargs): + """ + Initialize orientation object from two crystallographic directions. + + Parameters + ---------- + uvw : list, numpy.ndarray of shape (...,3) + lattice direction aligned with lab frame x-direction. + hkl : list, numpy.ndarray of shape (...,3) + lattice plane normal aligned with lab frame z-direction. + + """ + o = cls(**kwargs) + x = o.to_frame(uvw=uvw) + z = o.to_frame(hkl=hkl) + om = np.stack([x,np.cross(z,x),z],axis=-2) + return o.copy(rotation=Rotation.from_matrix(mechanics.transpose(om/np.linalg.norm(om,axis=-1,keepdims=True)))) - return (Orientation(r,self.lattice), i,j, k == 1) if symmetries else r # disorientation ... - # ... own sym, other sym, - # self-->other: True, self<--other: False @property - def in_FZ(self): - """Check if orientations fall into Fundamental Zone.""" - return self.lattice.in_FZ(self.rotation.as_Rodrigues(vector=True)) + def symmetry_operations(self): + """Symmetry operations as Rotations.""" + if self.family == 'cubic': + sym_quats = [ + [ 1.0, 0.0, 0.0, 0.0 ], + [ 0.0, 1.0, 0.0, 0.0 ], + [ 0.0, 0.0, 1.0, 0.0 ], + [ 0.0, 0.0, 0.0, 1.0 ], + [ 0.0, 0.0, 0.5*np.sqrt(2), 0.5*np.sqrt(2) ], + [ 0.0, 0.0, 0.5*np.sqrt(2),-0.5*np.sqrt(2) ], + [ 0.0, 0.5*np.sqrt(2), 0.0, 0.5*np.sqrt(2) ], + [ 0.0, 0.5*np.sqrt(2), 0.0, -0.5*np.sqrt(2) ], + [ 0.0, 0.5*np.sqrt(2),-0.5*np.sqrt(2), 0.0 ], + [ 0.0, -0.5*np.sqrt(2),-0.5*np.sqrt(2), 0.0 ], + [ 0.5, 0.5, 0.5, 0.5 ], + [-0.5, 0.5, 0.5, 0.5 ], + [-0.5, 0.5, 0.5, -0.5 ], + [-0.5, 0.5, -0.5, 0.5 ], + [-0.5, -0.5, 0.5, 0.5 ], + [-0.5, -0.5, 0.5, -0.5 ], + [-0.5, -0.5, -0.5, 0.5 ], + [-0.5, 0.5, -0.5, -0.5 ], + [-0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], + [ 0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], + [-0.5*np.sqrt(2), 0.0, 0.5*np.sqrt(2), 0.0 ], + [-0.5*np.sqrt(2), 0.0, -0.5*np.sqrt(2), 0.0 ], + [-0.5*np.sqrt(2), 0.5*np.sqrt(2), 0.0, 0.0 ], + [-0.5*np.sqrt(2),-0.5*np.sqrt(2), 0.0, 0.0 ], + ] + elif self.family == 'hexagonal': + sym_quats = [ + [ 1.0, 0.0, 0.0, 0.0 ], + [-0.5*np.sqrt(3), 0.0, 0.0, -0.5 ], + [ 0.5, 0.0, 0.0, 0.5*np.sqrt(3) ], + [ 0.0, 0.0, 0.0, 1.0 ], + [-0.5, 0.0, 0.0, 0.5*np.sqrt(3) ], + [-0.5*np.sqrt(3), 0.0, 0.0, 0.5 ], + [ 0.0, 1.0, 0.0, 0.0 ], + [ 0.0, -0.5*np.sqrt(3), 0.5, 0.0 ], + [ 0.0, 0.5, -0.5*np.sqrt(3), 0.0 ], + [ 0.0, 0.0, 1.0, 0.0 ], + [ 0.0, -0.5, -0.5*np.sqrt(3), 0.0 ], + [ 0.0, 0.5*np.sqrt(3), 0.5, 0.0 ], + ] + elif self.family == 'tetragonal': + sym_quats = [ + [ 1.0, 0.0, 0.0, 0.0 ], + [ 0.0, 1.0, 0.0, 0.0 ], + [ 0.0, 0.0, 1.0, 0.0 ], + [ 0.0, 0.0, 0.0, 1.0 ], + [ 0.0, 0.5*np.sqrt(2), 0.5*np.sqrt(2), 0.0 ], + [ 0.0, -0.5*np.sqrt(2), 0.5*np.sqrt(2), 0.0 ], + [ 0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], + [-0.5*np.sqrt(2), 0.0, 0.0, 0.5*np.sqrt(2) ], + ] + elif self.family == 'orthorhombic': + sym_quats = [ + [ 1.0,0.0,0.0,0.0 ], + [ 0.0,1.0,0.0,0.0 ], + [ 0.0,0.0,1.0,0.0 ], + [ 0.0,0.0,0.0,1.0 ], + ] + elif self.family == 'monoclinic': + sym_quats = [ + [ 1.0,0.0,0.0,0.0 ], + [ 0.0,0.0,1.0,0.0 ], + ] + elif self.family == 'triclinic': + sym_quats = [ + [ 1.0,0.0,0.0,0.0 ], + ] + else: + raise KeyError(f'Crystal family "{self.family}" is unknown') + + return Rotation.from_quaternion(sym_quats,accept_homomorph=True) @property def equivalent(self): """ - Orientations which are symmetrically equivalent. + Orientations that are symmetrically equivalent. - One dimension (length according to number of symmetrically equivalent orientations) + One dimension (length corresponds to number of symmetrically equivalent orientations) is added to the left of the Rotation array. """ - o = self.lattice.symmetry.symmetry_operations - o = o.reshape(o.shape[:1]+(1,)*len(self.rotation.shape)+(4,)) - o = Rotation(np.broadcast_to(o,o.shape[:1]+self.rotation.quaternion.shape)) + if self.family is None: + raise ValueError('Missing crystal symmetry') - s = np.broadcast_to(self.rotation.quaternion,o.shape[:1]+self.rotation.quaternion.shape) + o = self.symmetry_operations.broadcast_to(self.symmetry_operations.shape+self.shape,mode='right') + return self.copy(rotation=o@Rotation(self.quaternion).broadcast_to(o.shape,mode='left')) - return self.__class__(o@Rotation(s),self.lattice) + + @property + def reduced(self): + """Select symmetrically equivalent orientation that falls into fundamental zone according to symmetry.""" + if self.family is None: + raise ValueError('Missing crystal symmetry') + + eq = self.equivalent + ok = eq.in_FZ + ok &= np.cumsum(ok,axis=0) == 1 + loc = np.where(ok) + sort = 0 if len(loc) == 1 else np.lexsort(loc[:0:-1]) + return eq[ok][sort].reshape(self.shape) + + + @property + def in_FZ(self): + """ + Check whether orientation falls into fundamental zone of own symmetry. + + Returns + ------- + in : numpy.ndarray of quaternion.shape + Boolean array indicating whether Rodrigues-Frank vector falls into fundamental zone. + + Notes + ----- + Fundamental zones in Rodrigues space are point-symmetric around origin. + + References + ---------- + A. Heinz and P. Neumann, Acta Crystallographica Section A 47:780-789, 1991 + https://doi.org/10.1107/S0108767391006864 + + """ + if self.family is None: + raise ValueError('Missing crystal symmetry') + + rho_abs = np.abs(self.as_Rodrigues(vector=True)) + + with np.errstate(invalid='ignore'): + # using '*'/prod for 'and' + if self.family == 'cubic': + return (np.prod(np.sqrt(2)-1. >= rho_abs,axis=-1) * + (1. >= np.sum(rho_abs,axis=-1))).astype(np.bool) + elif self.family == 'hexagonal': + return (np.prod(1. >= rho_abs,axis=-1) * + (2. >= np.sqrt(3)*rho_abs[...,0] + rho_abs[...,1]) * + (2. >= np.sqrt(3)*rho_abs[...,1] + rho_abs[...,0]) * + (2. >= np.sqrt(3) + rho_abs[...,2])).astype(np.bool) + elif self.family == 'tetragonal': + return (np.prod(1. >= rho_abs[...,:2],axis=-1) * + (np.sqrt(2) >= rho_abs[...,0] + rho_abs[...,1]) * + (np.sqrt(2) >= rho_abs[...,2] + 1.)).astype(np.bool) + elif self.family == 'orthorhombic': + return (np.prod(1. >= rho_abs,axis=-1)).astype(np.bool) + elif self.family == 'monoclinic': + return (1. >= rho_abs[...,1]).astype(np.bool) + else: + return np.all(np.isfinite(rho_abs),axis=-1) + + + @property + def in_disorientation_FZ(self): + """ + Check whether orientation falls into fundamental zone of disorientations. + + Returns + ------- + in : numpy.ndarray of quaternion.shape + Boolean array indicating whether Rodrigues-Frank vector falls into disorientation FZ. + + References + ---------- + A. Heinz and P. Neumann, Acta Crystallographica Section A 47:780-789, 1991 + https://doi.org/10.1107/S0108767391006864 + + """ + if self.family is None: + raise ValueError('Missing crystal symmetry') + + rho = self.as_Rodrigues(vector=True) + + with np.errstate(invalid='ignore'): + if self.family == 'cubic': + return ((rho[...,0] >= rho[...,1]) & + (rho[...,1] >= rho[...,2]) & + (rho[...,2] >= 0)).astype(np.bool) + elif self.family == 'hexagonal': + return ((rho[...,0] >= rho[...,1]*np.sqrt(3)) & + (rho[...,1] >= 0) & + (rho[...,2] >= 0)).astype(np.bool) + elif self.family == 'tetragonal': + return ((rho[...,0] >= rho[...,1]) & + (rho[...,1] >= 0) & + (rho[...,2] >= 0)).astype(np.bool) + elif self.family == 'orthorhombic': + return ((rho[...,0] >= 0) & + (rho[...,1] >= 0) & + (rho[...,2] >= 0)).astype(np.bool) + elif self.family == 'monoclinic': + return ((rho[...,1] >= 0) & + (rho[...,2] >= 0)).astype(np.bool) + else: + return np.ones_like(rho[...,0],dtype=bool) + + + def relation_operations(self,model,return_lattice=False): + """ + Crystallographic orientation relationships for phase transformations. + + Parameters + ---------- + model : str + Name of orientation relationship. + return_lattice : bool, optional + Return the target lattice in addition. + + Returns + ------- + operations : Rotations + Rotations characterizing the orientation relationship. + + References + ---------- + S. Morito et al., Journal of Alloys and Compounds 577:s587-s592, 2013 + https://doi.org/10.1016/j.jallcom.2012.02.004 + + K. Kitahara et al., Acta Materialia 54(5):1279-1288, 2006 + https://doi.org/10.1016/j.actamat.2005.11.001 + + Y. He et al., Journal of Applied Crystallography 39:72-81, 2006 + https://doi.org/10.1107/S0021889805038276 + + H. Kitahara et al., Materials Characterization 54(4-5):378-386, 2005 + https://doi.org/10.1016/j.matchar.2004.12.015 + + Y. He et al., Acta Materialia 53(4):1179-1190, 2005 + https://doi.org/10.1016/j.actamat.2004.11.021 + + """ + from damask.lattice import relations + + if model not in relations: + raise KeyError(f'Orientation relationship "{model}" is unknown') + r = relations[model] + + if self.lattice not in r: + raise KeyError(f'Relationship "{model}" not supported for lattice "{self.lattice}"') + + sl = self.lattice + ol = (set(r)-{sl}).pop() + m = r[sl] + o = r[ol] + + p_,_p = np.zeros(m.shape[:-1]+(3,)),np.zeros(o.shape[:-1]+(3,)) + p_[...,0,:] = m[...,0,:] if m.shape[-1] == 3 else self.Bravais_to_Miller(uvtw=m[...,0,0:4]) + p_[...,1,:] = m[...,1,:] if m.shape[-1] == 3 else self.Bravais_to_Miller(hkil=m[...,1,0:4]) + _p[...,0,:] = o[...,0,:] if o.shape[-1] == 3 else self.Bravais_to_Miller(uvtw=o[...,0,0:4]) + _p[...,1,:] = o[...,1,:] if o.shape[-1] == 3 else self.Bravais_to_Miller(hkil=o[...,1,0:4]) + + return (Rotation.from_parallel(p_,_p),ol) \ + if return_lattice else \ + Rotation.from_parallel(p_,_p) def related(self,model): """ - Orientations related by the given orientation relationship. + Orientations derived from the given relationship. One dimension (length according to number of related orientations) is added to the left of the Rotation array. """ - o = Rotation.from_matrix(self.lattice.relation_operations(model)['rotations']).as_quaternion() - o = o.reshape(o.shape[:1]+(1,)*len(self.rotation.shape)+(4,)) - o = Rotation(np.broadcast_to(o,o.shape[:1]+self.rotation.quaternion.shape)) - - s = np.broadcast_to(self.rotation.quaternion,o.shape[:1]+self.rotation.quaternion.shape) - - return self.__class__(o@Rotation(s),self.lattice.relation_operations(model)['lattice']) + o,lattice = self.relation_operations(model,return_lattice=True) + target = Orientation(lattice=lattice) + o = o.broadcast_to(o.shape+self.shape,mode='right') + return self.copy(rotation=o@Rotation(self.quaternion).broadcast_to(o.shape,mode='left'), + lattice=lattice, + b = self.b if target.ratio['b'] is None else self.a*target.ratio['b'], + c = self.c if target.ratio['c'] is None else self.a*target.ratio['c'], + alpha = None if 'alpha' in target.immutable else self.alpha, + beta = None if 'beta' in target.immutable else self.beta, + gamma = None if 'gamma' in target.immutable else self.gamma, + ) @property - def reduced(self): - """Transform orientation to fall into fundamental zone according to symmetry.""" - eq = self.equivalent - in_FZ = eq.in_FZ - - # remove duplicates (occur for highly symmetric orientations) - found = np.zeros_like(in_FZ[0],dtype=bool) - q = self.rotation.quaternion[0] - for s in range(in_FZ.shape[0]): - #something fishy... why does q needs to be initialized? - q = np.where(np.expand_dims(np.logical_and(in_FZ[s],~found),-1),eq.rotation.quaternion[s],q) - found = np.logical_or(in_FZ[s],found) - - return self.__class__(q,self.lattice) + def parameters(self): + """Return lattice parameters a, b, c, alpha, beta, gamma.""" + return (self.a,self.b,self.c,self.alpha,self.beta,self.gamma) - def inverse_pole(self,axis,proper=False,SST=True): - """Axis rotated according to orientation (using crystal symmetry to ensure location falls into SST).""" - if SST: - eq = self.equivalent - pole = eq.rotation @ np.broadcast_to(axis/np.linalg.norm(axis),eq.rotation.shape+(3,)) - in_SST = self.lattice.in_SST(pole,proper=proper) - - # remove duplicates (occur for highly symmetric orientations) - found = np.zeros_like(in_SST[0],dtype=bool) - p = pole[0] - for s in range(in_SST.shape[0]): - p = np.where(np.expand_dims(np.logical_and(in_SST[s],~found),-1),pole[s],p) - found = np.logical_or(in_SST[s],found) - - return p + @property + def immutable(self): + """Return immutable parameters of own lattice.""" + if self.family == 'triclinic': + return {} + if self.family == 'monoclinic': + return { + 'alpha': np.pi/2., + 'gamma': np.pi/2., + } + if self.family == 'orthorhombic': + return { + 'alpha': np.pi/2., + 'beta': np.pi/2., + 'gamma': np.pi/2., + } + if self.family == 'tetragonal': + return { + 'b': 1.0, + 'alpha': np.pi/2., + 'beta': np.pi/2., + 'gamma': np.pi/2., + } + if self.family == 'hexagonal': + return { + 'b': 1.0, + 'alpha': np.pi/2., + 'beta': np.pi/2., + 'gamma': 2.*np.pi/3., + } + if self.family == 'cubic': + return { + 'b': 1.0, + 'c': 1.0, + 'alpha': np.pi/2., + 'beta': np.pi/2., + 'gamma': np.pi/2., + } else: - return self.rotation @ np.broadcast_to(axis/np.linalg.norm(axis),self.rotation.shape+(3,)) + raise KeyError(f'Crystal family "{self.family}" is unknown') + @property + def ratio(self): + """Return axes ratios of own lattice.""" + _ratio = { 'hexagonal': {'c': np.sqrt(8./3.)}} + + return dict(b = self.immutable['b'] + if 'b' in self.immutable else + _ratio[self.family]['b'] if self.family in _ratio and 'b' in _ratio[self.family] else None, + c = self.immutable['c'] + if 'c' in self.immutable else + _ratio[self.family]['c'] if self.family in _ratio and 'c' in _ratio[self.family] else None, + ) + + + @property + def basis_real(self): + """ + Calculate orthogonal real space crystal basis. + + References + ---------- + C.T. Young and J.L. Lytton, J. Appl. Phys. 43:1408–1417, 1972 + "Computer Generation and Identification of Kikuchi Projections" + https://doi.org/10.1063/1.1661333 + + """ + if None in self.parameters: + raise KeyError('Missing crystal lattice parameters') + return np.array([ + [1,0,0], + [np.cos(self.gamma),np.sin(self.gamma),0], + [np.cos(self.beta), + (np.cos(self.alpha)-np.cos(self.beta)*np.cos(self.gamma)) /np.sin(self.gamma), + np.sqrt(1 - np.cos(self.alpha)**2 - np.cos(self.beta)**2 - np.cos(self.gamma)**2 + + 2 * np.cos(self.alpha) * np.cos(self.beta) * np.cos(self.gamma))/np.sin(self.gamma)], + ],dtype=float).T \ + * np.array([self.a,self.b,self.c]) + + + @property + def basis_reciprocal(self): + """Calculate reciprocal (dual) crystal basis.""" + return np.linalg.inv(self.basis_real.T) + + + def in_SST(self,vector,proper=False): + """ + Check whether given crystal frame vector falls into standard stereographic triangle of own symmetry. + + Parameters + ---------- + vector : numpy.ndarray of shape (...,3) + Vector to check. + proper : bool, optional + Consider only vectors with z >= 0, hence combine two neighboring SSTs. + Defaults to False. + + Returns + ------- + in : numpy.ndarray of shape (...) + Boolean array indicating whether vector falls into SST. + + References + ---------- + Bases are computed from + + >>> basis = { + ... 'cubic' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,1.]/np.sqrt(2.), # green + ... [1.,1.,1.]/np.sqrt(3.)]).T), # blue + ... 'hexagonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,0.], # green + ... [np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T), # blue + ... 'tetragonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,0.], # green + ... [1.,1.,0.]/np.sqrt(2.)]).T), # blue + ... 'orthorhombic': np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,0.], # green + ... [0.,1.,0.]]).T), # blue + ... } + + """ + if not isinstance(vector,np.ndarray) or vector.shape[-1] != 3: + raise ValueError('Input is not a field of three-dimensional vectors.') + + if self.family == 'cubic': + basis = {'improper':np.array([ [-1. , 0. , 1. ], + [ np.sqrt(2.) , -np.sqrt(2.) , 0. ], + [ 0. , np.sqrt(3.) , 0. ] ]), + 'proper':np.array([ [ 0. , -1. , 1. ], + [-np.sqrt(2.) , np.sqrt(2.) , 0. ], + [ np.sqrt(3.) , 0. , 0. ] ]), + } + elif self.family == 'hexagonal': + basis = {'improper':np.array([ [ 0. , 0. , 1. ], + [ 1. , -np.sqrt(3.) , 0. ], + [ 0. , 2. , 0. ] ]), + 'proper':np.array([ [ 0. , 0. , 1. ], + [-1. , np.sqrt(3.) , 0. ], + [ np.sqrt(3.) , -1. , 0. ] ]), + } + elif self.family == 'tetragonal': + basis = {'improper':np.array([ [ 0. , 0. , 1. ], + [ 1. , -1. , 0. ], + [ 0. , np.sqrt(2.) , 0. ] ]), + 'proper':np.array([ [ 0. , 0. , 1. ], + [-1. , 1. , 0. ], + [ np.sqrt(2.) , 0. , 0. ] ]), + } + elif self.family == 'orthorhombic': + basis = {'improper':np.array([ [ 0., 0., 1.], + [ 1., 0., 0.], + [ 0., 1., 0.] ]), + 'proper':np.array([ [ 0., 0., 1.], + [-1., 0., 0.], + [ 0., 1., 0.] ]), + } + else: # direct exit for unspecified symmetry + return np.ones_like(vector[...,0],bool) + + if proper: + components_proper = np.around(np.einsum('...ji,...i', + np.broadcast_to(basis['proper'], vector.shape+(3,)), + vector), 12) + components_improper = np.around(np.einsum('...ji,...i', + np.broadcast_to(basis['improper'], vector.shape+(3,)), + vector), 12) + return np.all(components_proper >= 0.0,axis=-1) \ + | np.all(components_improper >= 0.0,axis=-1) + else: + components = np.around(np.einsum('...ji,...i', + np.broadcast_to(basis['improper'], vector.shape+(3,)), + np.block([vector[...,:2],np.abs(vector[...,2:3])])), 12) + + return np.all(components >= 0.0,axis=-1) + + + def IPF_color(self,vector,proper=False): + """ + Map vector to RGB color within standard stereographic triangle of own symmetry. + + Parameters + ---------- + vector : numpy.ndarray of shape (...,3) + Vector to colorize. + proper : bool, optional + Consider only vectors with z >= 0, hence combine two neighboring SSTs (with mirrored colors). + Defaults to False. + + Returns + ------- + rgb : numpy.ndarray of shape (...,3) + RGB array of IPF colors. + + References + ---------- + Bases are computed from + + >>> basis = { + ... 'cubic' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,1.]/np.sqrt(2.), # green + ... [1.,1.,1.]/np.sqrt(3.)]).T), # blue + ... 'hexagonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,0.], # green + ... [np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T), # blue + ... 'tetragonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,0.], # green + ... [1.,1.,0.]/np.sqrt(2.)]).T), # blue + ... 'orthorhombic': np.linalg.inv(np.array([[0.,0.,1.], # direction of red + ... [1.,0.,0.], # green + ... [0.,1.,0.]]).T), # blue + ... } + + """ + if vector.shape[-1] != 3: + raise ValueError('Input is not a field of three-dimensional vectors.') + + if self.family == 'cubic': + basis = {'improper':np.array([ [-1. , 0. , 1. ], + [ np.sqrt(2.) , -np.sqrt(2.) , 0. ], + [ 0. , np.sqrt(3.) , 0. ] ]), + 'proper':np.array([ [ 0. , -1. , 1. ], + [-np.sqrt(2.) , np.sqrt(2.) , 0. ], + [ np.sqrt(3.) , 0. , 0. ] ]), + } + elif self.family == 'hexagonal': + basis = {'improper':np.array([ [ 0. , 0. , 1. ], + [ 1. , -np.sqrt(3.) , 0. ], + [ 0. , 2. , 0. ] ]), + 'proper':np.array([ [ 0. , 0. , 1. ], + [-1. , np.sqrt(3.) , 0. ], + [ np.sqrt(3.) , -1. , 0. ] ]), + } + elif self.family == 'tetragonal': + basis = {'improper':np.array([ [ 0. , 0. , 1. ], + [ 1. , -1. , 0. ], + [ 0. , np.sqrt(2.) , 0. ] ]), + 'proper':np.array([ [ 0. , 0. , 1. ], + [-1. , 1. , 0. ], + [ np.sqrt(2.) , 0. , 0. ] ]), + } + elif self.family == 'orthorhombic': + basis = {'improper':np.array([ [ 0., 0., 1.], + [ 1., 0., 0.], + [ 0., 1., 0.] ]), + 'proper':np.array([ [ 0., 0., 1.], + [-1., 0., 0.], + [ 0., 1., 0.] ]), + } + else: # direct exit for unspecified symmetry + return np.zeros_like(vector) + + if proper: + components_proper = np.around(np.einsum('...ji,...i', + np.broadcast_to(basis['proper'], vector.shape+(3,)), + vector), 12) + components_improper = np.around(np.einsum('...ji,...i', + np.broadcast_to(basis['improper'], vector.shape+(3,)), + vector), 12) + in_SST = np.all(components_proper >= 0.0,axis=-1) \ + | np.all(components_improper >= 0.0,axis=-1) + components = np.where((in_SST & np.all(components_proper >= 0.0,axis=-1))[...,np.newaxis], + components_proper,components_improper) + else: + components = np.around(np.einsum('...ji,...i', + np.broadcast_to(basis['improper'], vector.shape+(3,)), + np.block([vector[...,:2],np.abs(vector[...,2:3])])), 12) + + in_SST = np.all(components >= 0.0,axis=-1) + + with np.errstate(invalid='ignore',divide='ignore'): + rgb = (components/np.linalg.norm(components,axis=-1,keepdims=True))**0.5 # smoothen color ramps + rgb = np.clip(rgb,0.,1.) # clip intensity + rgb /= np.max(rgb,axis=-1,keepdims=True) # normalize to (HS)V = 1 + rgb[np.broadcast_to(~in_SST[...,np.newaxis],rgb.shape)] = 0.0 + return rgb + + + def disorientation(self,other,return_operators=False): + """ + Calculate disorientation between myself and given other orientation. + + Parameters + ---------- + other : Orientation + Orientation to calculate disorientation for. + Shape of other blends with shape of own rotation array. + For example, shapes of (2,3) for own rotations and (3,2) for other's result in (2,3,2) disorientations. + return_operators : bool, optional + Return index pair of symmetrically equivalent orientations that result in disorientation axis falling into FZ. + Defaults to False. + + Returns + ------- + disorientation : Orientation + Disorientation between self and other. + operators : numpy.ndarray int of shape (...,2), conditional + Index of symmetrically equivalent orientation that rotated vector to the SST. + + Notes + ----- + Currently requires same crystal family for both orientations. + For extension to cases with differing symmetry see A. Heinz and P. Neumann 1991 and 10.1107/S0021889808016373. + + """ + if self.family is None or other.family is None: + raise ValueError('Missing crystal symmetry') + if self.family != other.family: + raise NotImplementedError('Disorientation between different crystal families not supported yet.') + + blend = util.shapeblender(self.shape,other.shape) + s = self.equivalent + o = other.equivalent + + s_ = s.reshape((s.shape[0],1)+ self.shape).broadcast_to((s.shape[0],o.shape[0])+blend,mode='right') + o_ = o.reshape((1,o.shape[0])+other.shape).broadcast_to((s.shape[0],o.shape[0])+blend,mode='right') + r_ = s_.misorientation(o_) + _r = ~r_ + + forward = r_.in_FZ & r_.in_disorientation_FZ + reverse = _r.in_FZ & _r.in_disorientation_FZ + ok = forward | reverse + ok &= (np.cumsum(ok.reshape((-1,)+ok.shape[2:]),axis=0) == 1).reshape(ok.shape) + r = np.where(np.any(forward[...,np.newaxis],axis=(0,1),keepdims=True), + r_.quaternion, + _r.quaternion) + loc = np.where(ok) + sort = 0 if len(loc) == 2 else np.lexsort(loc[:1:-1]) + quat = r[ok][sort].reshape(blend+(4,)) + + return ( + (self.copy(rotation=quat), + (np.vstack(loc[:2]).T)[sort].reshape(blend+(2,))) + if return_operators else + self.copy(rotation=quat) + ) + + + def average(self,weights=None,return_cloud=False): + """ + Return orientation average over last dimension. + + Parameters + ---------- + weights : numpy.ndarray, optional + Relative weights of orientations. + return_cloud : bool, optional + Return the set of symmetrically equivalent orientations that was used in averaging. + Defaults to False. + + Returns + ------- + average : Orientation + Weighted average of original Orientation field. + cloud : Orientations, conditional + Set of symmetrically equivalent orientations that were used in averaging. + + References + ---------- + J.C. Glez and J. Driver, J. Appl. Cryst. 34:280-288, 2001 + "Orientation distribution analysis in deformed grains" + https://doi.org/10.1107/S0021889801003077 + + """ + if self.family is None: + raise ValueError('Missing crystal symmetry') - def IPF_color(self,axis): #ToDo axis or direction? - """TSL color of inverse pole figure for given axis.""" eq = self.equivalent - pole = eq.rotation @ np.broadcast_to(axis/np.linalg.norm(axis),eq.rotation.shape+(3,)) - in_SST, color = self.lattice.in_SST(pole,color=True) - - # remove duplicates (occur for highly symmetric orientations) - found = np.zeros_like(in_SST[0],dtype=bool) - c = color[0] - for s in range(in_SST.shape[0]): - c = np.where(np.expand_dims(np.logical_and(in_SST[s],~found),-1),color[s],c) - found = np.logical_or(in_SST[s],found) - - return c + m = eq.misorientation(self[...,0].reshape((1,)+self.shape[:-1]+(1,)) + .broadcast_to(eq.shape))\ + .as_axis_angle()[...,3] + r = Rotation(np.squeeze(np.take_along_axis(eq.quaternion, + np.argmin(m,axis=0)[np.newaxis,...,np.newaxis], + axis=0), + axis=0)) + return ( + (self.copy(rotation=Rotation(r).average(weights)), + self.copy(rotation=Rotation(r))) + if return_cloud else + self.copy(rotation=Rotation(r).average(weights)) + ) - # ToDo: Discuss vectorization/calling signature - @staticmethod - def from_average(orientations, - weights = []): - """Create orientation from average of list of orientations.""" - # further read: Orientation distribution analysis in deformed grains - # https://doi.org/10.1107/S0021889801003077 - if not all(isinstance(item, Orientation) for item in orientations): - raise TypeError("Only instances of Orientation can be averaged.") + def to_SST(self,vector,proper=False,return_operators=False): + """ + Rotate vector to ensure it falls into (improper or proper) standard stereographic triangle of crystal symmetry. - closest = [] - ref = orientations[0] - for o in orientations: - closest.append(o.equivalent[ - ref.disorientation(o, - SST = False, # select (o[ther]'s) sym orientation - symmetries = True)[2]].rotation) # with lowest misorientation + Parameters + ---------- + vector : numpy.ndarray of shape (...,3) + Lab frame vector to align with crystal frame direction. + Shape of other blends with shape of own rotation array. + For example, a rotation array of shape (3,2) and a (2,4) vector array result in (3,2,4) outputs. + proper : bool, optional + Consider only vectors with z >= 0, hence combine two neighboring SSTs. + Defaults to False. + return_operators : bool, optional + Return the symmetrically equivalent orientation that rotated vector to SST. + Defaults to False. - return Orientation(Rotation.from_average(closest,weights),ref.lattice) + Returns + ------- + vector_SST : numpy.ndarray of shape (...,3) + Rotated vector falling into SST. + operators : numpy.ndarray int of shape (...), conditional + Index of symmetrically equivalent orientation that rotated vector to SST. + + """ + if self.family is None: + raise ValueError('Missing crystal symmetry') + + eq = self.equivalent + blend = util.shapeblender(eq.shape,vector.shape[:-1]) + poles = eq.broadcast_to(blend,mode='right') @ np.broadcast_to(vector,blend+(3,)) + ok = self.in_SST(poles,proper=proper) + ok &= np.cumsum(ok,axis=0) == 1 + loc = np.where(ok) + sort = 0 if len(loc) == 1 else np.lexsort(loc[:0:-1]) + return ( + (poles[ok][sort].reshape(blend[1:]+(3,)), (np.vstack(loc[:1]).T)[sort].reshape(blend[1:])) + if return_operators else + poles[ok][sort].reshape(blend[1:]+(3,)) + ) - # ToDo: Discuss vectorization/calling signature - def average(self,other): - """Calculate the average rotation.""" - return Orientation.from_average([self,other]) + @classmethod + def Bravais_to_Miller(cls,uvtw=None,hkil=None): + """ + Transform 4 Miller–Bravais indices to 3 Miller indices of crystal direction [uvw] or plane normal (hkl). + + Parameters + ---------- + uvtw | hkil : numpy.ndarray of shape (...,4) + Miller–Bravais indices of crystallographic direction [uvtw] or plane normal (hkil). + + Returns + ------- + uvw | hkl : numpy.ndarray of shape (...,3) + Miller indices of [uvw] direction or (hkl) plane normal. + + """ + if (uvtw is not None) ^ (hkil is None): + raise KeyError('Specify either "uvtw" or "hkil"') + axis,basis = (np.array(uvtw),np.array([[1,0,-1,0], + [0,1,-1,0], + [0,0, 0,1]])) \ + if hkil is None else \ + (np.array(hkil),np.array([[1,0,0,0], + [0,1,0,0], + [0,0,0,1]])) + return np.einsum('il,...l->...i',basis,axis) + + + @classmethod + def Miller_to_Bravais(cls,uvw=None,hkl=None): + """ + Transform 3 Miller indices to 4 Miller–Bravais indices of crystal direction [uvtw] or plane normal (hkil). + + Parameters + ---------- + uvw | hkl : numpy.ndarray of shape (...,3) + Miller indices of crystallographic direction [uvw] or plane normal (hkl). + + Returns + ------- + uvtw | hkil : numpy.ndarray of shape (...,4) + Miller–Bravais indices of [uvtw] direction or (hkil) plane normal. + + """ + if (uvw is not None) ^ (hkl is None): + raise KeyError('Specify either "uvw" or "hkl"') + axis,basis = (np.array(uvw),np.array([[ 2,-1, 0], + [-1, 2, 0], + [-1,-1, 0], + [ 0, 0, 3]])/3) \ + if hkl is None else \ + (np.array(hkl),np.array([[ 1, 0, 0], + [ 0, 1, 0], + [-1,-1, 0], + [ 0, 0, 1]])) + return np.einsum('il,...l->...i',basis,axis) + + + def to_lattice(self,direction=None,plane=None): + """ + Calculate lattice vector corresponding to crystal frame direction or plane normal. + + Parameters + ---------- + direction | normal : numpy.ndarray of shape (...,3) + Vector along direction or plane normal. + + Returns + ------- + Miller : numpy.ndarray of shape (...,3) + lattice vector of direction or plane. + Use util.scale_to_coprime to convert to (integer) Miller indices. + + """ + if (direction is not None) ^ (plane is None): + raise KeyError('Specify either "direction" or "plane"') + axis,basis = (np.array(direction),self.basis_reciprocal.T) \ + if plane is None else \ + (np.array(plane),self.basis_real.T) + return np.einsum('il,...l->...i',basis,axis) + + + def to_frame(self,uvw=None,hkl=None,with_symmetry=False): + """ + Calculate crystal frame vector along lattice direction [uvw] or plane normal (hkl). + + Parameters + ---------- + uvw | hkl : numpy.ndarray of shape (...,3) + Miller indices of crystallographic direction or plane normal. + with_symmetry : bool, optional + Calculate all N symmetrically equivalent vectors. + + Returns + ------- + vector : numpy.ndarray of shape (...,3) or (N,...,3) + Crystal frame vector (or vectors if with_symmetry) along [uvw] direction or (hkl) plane normal. + + """ + if (uvw is not None) ^ (hkl is None): + raise KeyError('Specify either "uvw" or "hkl"') + axis,basis = (np.array(uvw),self.basis_real) \ + if hkl is None else \ + (np.array(hkl),self.basis_reciprocal) + return (self.symmetry_operations.broadcast_to(self.symmetry_operations.shape+axis.shape[:-1],mode='right') + @ np.broadcast_to(np.einsum('il,...l->...i',basis,axis),self.symmetry_operations.shape+axis.shape) + if with_symmetry else + np.einsum('il,...l->...i',basis,axis)) + + + def to_pole(self,uvw=None,hkl=None,with_symmetry=False): + """ + Calculate lab frame vector along lattice direction [uvw] or plane normal (hkl). + + Parameters + ---------- + uvw | hkl : numpy.ndarray of shape (...,3) + Miller indices of crystallographic direction or plane normal. + with_symmetry : bool, optional + Calculate all N symmetrically equivalent vectors. + + Returns + ------- + vector : numpy.ndarray of shape (...,3) or (N,...,3) + Lab frame vector (or vectors if with_symmetry) along [uvw] direction or (hkl) plane normal. + + """ + v = self.to_frame(uvw=uvw,hkl=hkl,with_symmetry=with_symmetry) + return ~(self if self.shape+v.shape[:-1] == () else self.broadcast_to(self.shape+v.shape[:-1],mode='right')) \ + @ np.broadcast_to(v,self.shape+v.shape) + + + def Schmid(self,mode): + u""" + Calculate Schmid matrix P = d ⨂ n in the lab frame for given lattice shear kinematics. + + Parameters + ---------- + mode : str + Type of kinematics, e.g. 'slip' or 'twin'. + + Returns + ------- + P : numpy.ndarray of shape (...,N,3,3) + Schmid matrix for each of the N deformation systems. + + """ + d = self.to_frame(uvw=self.kinematics[mode]['direction'],with_symmetry=False) + p = self.to_frame(hkl=self.kinematics[mode]['plane'] ,with_symmetry=False) + P = np.einsum('...i,...j->...ij',d/np.linalg.norm(d,axis=-1,keepdims=True), + p/np.linalg.norm(p,axis=-1,keepdims=True)) + + return ~self.broadcast_to( self.shape+P.shape[:-2],mode='right') \ + @ np.broadcast_to(P,self.shape+P.shape) diff --git a/python/damask/_result.py b/python/damask/_result.py index f0406ecad..5bed4347f 100644 --- a/python/damask/_result.py +++ b/python/damask/_result.py @@ -15,7 +15,6 @@ from numpy.lib import recfunctions as rfn import damask from . import VTK from . import Table -from . import Rotation from . import Orientation from . import grid_filters from . import mechanics @@ -743,11 +742,13 @@ class Result: def _add_IPF_color(q,l): m = util.scale_to_coprime(np.array(l)) - o = Orientation(Rotation(rfn.structured_to_unstructured(q['data'])), - lattice = q['meta']['Lattice']) + o = Orientation(rotation = (rfn.structured_to_unstructured(q['data'])), + lattice = {'fcc':'cF', + 'bcc':'cI', + 'hex':'hP'}[q['meta']['Lattice']]) return { - 'data': np.uint8(o.IPF_color(l)*255), + 'data': np.uint8(o.IPF_color(o.to_SST(l))*255), 'label': 'IPFcolor_[{} {} {}]'.format(*m), 'meta' : { 'Unit': '8-bit RGB', @@ -897,42 +898,47 @@ class Result: self._add_generic_pointwise(self._add_PK2,{'P':P,'F':F}) - @staticmethod - def _add_pole(q,p,polar): - pole = np.array(p) - unit_pole = pole/np.linalg.norm(pole) - m = util.scale_to_coprime(pole) - rot = Rotation(q['data'].view(np.double).reshape(-1,4)) +# The add_pole functionality needs discussion. +# The new Crystal object can perform such a calculation but the outcome depends on the lattice parameters +# as well as on whether a direction or plane is concerned (see the DAMASK_examples/pole_figure notebook). +# Below code appears to be too simplistic. - rotatedPole = rot @ np.broadcast_to(unit_pole,rot.shape+(3,)) # rotate pole according to crystal orientation - xy = rotatedPole[:,0:2]/(1.+abs(unit_pole[2])) # stereographic projection - coords = xy if not polar else \ - np.block([np.sqrt(xy[:,0:1]*xy[:,0:1]+xy[:,1:2]*xy[:,1:2]),np.arctan2(xy[:,1:2],xy[:,0:1])]) - return { - 'data': coords, - 'label': 'p^{}_[{} {} {})'.format(u'rφ' if polar else 'xy',*m), - 'meta' : { - 'Unit': '1', - 'Description': '{} coordinates of stereographic projection of pole (direction/plane) in crystal frame'\ - .format('Polar' if polar else 'Cartesian'), - 'Creator': 'add_pole' - } - } - def add_pole(self,q,p,polar=False): - """ - Add coordinates of stereographic projection of given pole in crystal frame. - - Parameters - ---------- - q : str - Label of the dataset containing the crystallographic orientation as quaternions. - p : numpy.array of shape (3) - Crystallographic direction or plane. - polar : bool, optional - Give pole in polar coordinates. Defaults to False. - - """ - self._add_generic_pointwise(self._add_pole,{'q':q},{'p':p,'polar':polar}) + # @staticmethod + # def _add_pole(q,p,polar): + # pole = np.array(p) + # unit_pole = pole/np.linalg.norm(pole) + # m = util.scale_to_coprime(pole) + # rot = Rotation(q['data'].view(np.double).reshape(-1,4)) + # + # rotatedPole = rot @ np.broadcast_to(unit_pole,rot.shape+(3,)) # rotate pole according to crystal orientation + # xy = rotatedPole[:,0:2]/(1.+abs(unit_pole[2])) # stereographic projection + # coords = xy if not polar else \ + # np.block([np.sqrt(xy[:,0:1]*xy[:,0:1]+xy[:,1:2]*xy[:,1:2]),np.arctan2(xy[:,1:2],xy[:,0:1])]) + # return { + # 'data': coords, + # 'label': 'p^{}_[{} {} {})'.format(u'rφ' if polar else 'xy',*m), + # 'meta' : { + # 'Unit': '1', + # 'Description': '{} coordinates of stereographic projection of pole (direction/plane) in crystal frame'\ + # .format('Polar' if polar else 'Cartesian'), + # 'Creator': 'add_pole' + # } + # } + # def add_pole(self,q,p,polar=False): + # """ + # Add coordinates of stereographic projection of given pole in crystal frame. + # + # Parameters + # ---------- + # q : str + # Label of the dataset containing the crystallographic orientation as quaternions. + # p : numpy.array of shape (3) + # Crystallographic direction or plane. + # polar : bool, optional + # Give pole in polar coordinates. Defaults to False. + # + # """ + # self._add_generic_pointwise(self._add_pole,{'q':q},{'p':p,'polar':polar}) @staticmethod diff --git a/python/damask/_rotation.py b/python/damask/_rotation.py index b0c9dab5b..deae8b6ac 100644 --- a/python/damask/_rotation.py +++ b/python/damask/_rotation.py @@ -13,18 +13,18 @@ _R1 = (3.*np.pi/4.)**(1./3.) class Rotation: u""" - Orientation stored with functionality for conversion to different representations. + Rotation with functionality for conversion between different representations. The following conventions apply: - - coordinate frames are right-handed. - - a rotation angle ω is taken to be positive for a counterclockwise rotation + - Coordinate frames are right-handed. + - A rotation angle ω is taken to be positive for a counterclockwise rotation when viewing from the end point of the rotation axis towards the origin. - - rotations will be interpreted in the passive sense. + - Rotations will be interpreted in the passive sense. - Euler angle triplets are implemented using the Bunge convention, - with the angular ranges as [0,2π], [0,π], [0,2π]. - - the rotation angle ω is limited to the interval [0,π]. - - the real part of a quaternion is positive, Re(q) > 0 + with angular ranges of [0,2π], [0,π], [0,2π]. + - The rotation angle ω is limited to the interval [0,π]. + - The real part of a quaternion is positive, Re(q) > 0 - P = -1 (as default). Examples @@ -33,7 +33,7 @@ class Rotation: coordinates "b" expressed in system "B": - b = Q @ a - - b = np.dot(Q.asMatrix(),a) + - b = np.dot(Q.as_matrix(),a) References ---------- @@ -44,20 +44,83 @@ class Rotation: __slots__ = ['quaternion'] - def __init__(self,quaternion = np.array([1.0,0.0,0.0,0.0])): + def __init__(self,rotation = np.array([1.0,0.0,0.0,0.0])): """ - Initializes to identity unless specified. + Initialize rotation object. Parameters ---------- - quaternion : numpy.ndarray, optional + rotation : list, numpy.ndarray, Rotation, optional Unit quaternion in positive real hemisphere. Use .from_quaternion to perform a sanity check. + Defaults to no rotation. """ - if quaternion.shape[-1] != 4: - raise ValueError('Not a quaternion') - self.quaternion = quaternion.copy() + if isinstance(rotation,Rotation): + self.quaternion = rotation.quaternion.copy() + elif np.array(rotation).shape[-1] == 4: + self.quaternion = np.array(rotation) + else: + raise ValueError('"rotation" is neither a Rotation nor a quaternion') + + + def __repr__(self): + """Represent rotation as unit quaternion, rotation matrix, and Bunge-Euler angles.""" + return 'Quaternions:\n'+str(self.quaternion) \ + if self.quaternion.shape != (4,) else \ + '\n'.join([ + 'Quaternion: (real={:.3f}, imag=<{:+.3f}, {:+.3f}, {:+.3f}>)'.format(*(self.quaternion)), + 'Matrix:\n{}'.format(np.round(self.as_matrix(),8)), + 'Bunge Eulers / deg: ({:3.2f}, {:3.2f}, {:3.2f})'.format(*self.as_Eulers(degrees=True)), + ]) + + + # ToDo: Check difference __copy__ vs __deepcopy__ + def __copy__(self,**kwargs): + """Copy.""" + return self.__class__(rotation=kwargs['rotation'] if 'rotation' in kwargs else self.quaternion) + + copy = __copy__ + + + def __getitem__(self,item): + """Return slice according to item.""" + return self.copy() \ + if self.shape == () else \ + self.copy(rotation=self.quaternion[item+(slice(None),)] if isinstance(item,tuple) else self.quaternion[item]) + + + def __eq__(self,other): + """ + Equal to other. + + Equality is determined taking limited floating point precision into + account. See numpy.allclose for details. + + Parameters + ---------- + other : Rotation + Rotation to check for equality. + + """ + return np.prod(self.shape,dtype=int) == np.prod(other.shape,dtype=int) \ + and np.allclose(self.quaternion,other.quaternion) + + + def __neq__(self,other): + """ + Not Equal to other. + + Equality is determined taking limited floating point precision into + account. See numpy.allclose for details. + + Parameters + ---------- + other : Rotation + Rotation to check for inequality. + + """ + return not self.__eq__(other) @property @@ -65,39 +128,36 @@ class Rotation: return self.quaternion.shape[:-1] - # ToDo: Check difference __copy__ vs __deepcopy__ - def __copy__(self): - """Copy.""" - return self.__class__(self.quaternion) - - copy = __copy__ - - - def __repr__(self): - """Orientation displayed as unit quaternion, rotation matrix, and Bunge-Euler angles.""" - if self.quaternion.shape != (4,): - return 'Quaternions:\n'+str(self.quaternion) # ToDo: could be nicer ... - return '\n'.join([ - 'Quaternion: (real={:.3f}, imag=<{:+.3f}, {:+.3f}, {:+.3f}>)'.format(*(self.quaternion)), - 'Matrix:\n{}'.format(np.round(self.as_matrix(),8)), - 'Bunge Eulers / deg: ({:3.2f}, {:3.2f}, {:3.2f})'.format(*self.as_Eulers(degrees=True)), - ]) - - - def __getitem__(self,item): - """Iterate over leading/leftmost dimension of Rotation array.""" - if self.shape == (): return self.copy() - if isinstance(item,tuple) and len(item) >= len(self): - raise IndexError('Too many indices') - return self.__class__(self.quaternion[item]) - - def __len__(self): """Length of leading/leftmost dimension of Rotation array.""" return 0 if self.shape == () else self.shape[0] - def __matmul__(self, other): + def __invert__(self): + """Inverse rotation (backward rotation).""" + dup = self.copy() + dup.quaternion[...,1:] *= -1 + return dup + + + def __pow__(self,pwr): + """ + Raise quaternion to power. + + Equivalent to performing the rotation 'pwr' times. + + Parameters + ---------- + pwr : float + Power to raise quaternion to. + + """ + phi = np.arccos(self.quaternion[...,0:1]) + p = self.quaternion[...,1:]/np.linalg.norm(self.quaternion[...,1:],axis=-1,keepdims=True) + return self.copy(rotation=Rotation(np.block([np.cos(pwr*phi),np.sin(pwr*phi)*p]))._standardize()) + + + def __matmul__(self,other): """ Rotation of vector, second or fourth order tensor, or rotation object. @@ -112,14 +172,14 @@ class Rotation: Rotated vector, second or fourth order tensor, or rotation object. """ - if isinstance(other, Rotation): + if isinstance(other,Rotation): q_m = self.quaternion[...,0:1] p_m = self.quaternion[...,1:] q_o = other.quaternion[...,0:1] p_o = other.quaternion[...,1:] q = (q_m*q_o - np.einsum('...i,...i',p_m,p_o).reshape(self.shape+(1,))) p = q_m*p_o + q_o*p_m + _P * np.cross(p_m,p_o) - return self.__class__(np.block([q,p]))._standardize() + return Rotation(np.block([q,p]))._standardize() elif isinstance(other,np.ndarray): if self.shape + (3,) == other.shape: @@ -146,27 +206,89 @@ class Rotation: def _standardize(self): - """Standardize (ensure positive real hemisphere).""" + """Standardize quaternion (ensure positive real hemisphere).""" self.quaternion[self.quaternion[...,0] < 0.0] *= -1 return self - def inverse(self): - """In-place inverse rotation (backward rotation).""" - self.quaternion[...,1:] *= -1 - return self - def __invert__(self): - """Inverse rotation (backward rotation).""" - return self.copy().inverse() + def append(self,other): + """Extend rotation array along first dimension with other array.""" + return self.copy(rotation=np.vstack((self.quaternion,other.quaternion))) - def inversed(self): - """Inverse rotation (backward rotation).""" - return ~ self + + def flatten(self,order = 'C'): + """Flatten quaternion array.""" + return self.copy(rotation=self.quaternion.reshape((-1,4),order=order)) + + + def reshape(self,shape,order = 'C'): + """Reshape quaternion array.""" + if isinstance(shape,(int,np.integer)): shape = (shape,) + return self.copy(rotation=self.quaternion.reshape(tuple(shape)+(4,),order=order)) + + + def broadcast_to(self,shape,mode = 'right'): + """ + Broadcast quaternion array to shape. + + Parameters + ---------- + shape : tuple + Shape of broadcasted array. + mode : str, optional + Where to preferentially locate missing dimensions. + Either 'left' or 'right' (default). + + """ + if isinstance(shape,(int,np.integer)): shape = (shape,) + return self.copy(rotation=np.broadcast_to(self.quaternion.reshape(util.shapeshifter(self.shape,shape,mode)+(4,)), + shape+(4,))) + + + def average(self,weights = None): + """ + Average rotations along last dimension. + + Parameters + ---------- + weights : list of floats, optional + Relative weight of each rotation. + + Returns + ------- + average : Rotation + Weighted average of original Rotation field. + + References + ---------- + Quaternion averaging + F. Landis Markley, Yang Cheng, John L. Crassidis, Yaakov Oshman + Journal of Guidance, Control, and Dynamics 30(4):1193-1197, 2007 + 10.2514/1.28949 + + """ + def _M(quat): + """Intermediate representation supporting quaternion averaging.""" + return np.einsum('...i,...j',quat,quat) + + if not weights: + weights = np.ones(self.shape,dtype=float) + + eig, vec = np.linalg.eig(np.sum(_M(self.quaternion) * weights[...,np.newaxis,np.newaxis],axis=-3) \ + /np.sum( weights[...,np.newaxis,np.newaxis],axis=-3)) + + return Rotation.from_quaternion(np.real( + np.squeeze( + np.take_along_axis(vec, + eig.argmax(axis=-1)[...,np.newaxis,np.newaxis], + axis=-1), + axis=-1)), + accept_homomorph = True) def misorientation(self,other): """ - Get Misorientation. + Calculate misorientation from self to other Rotation. Parameters ---------- @@ -177,33 +299,6 @@ class Rotation: return other@~self - def broadcast_to(self,shape): - if isinstance(shape,(int,np.integer)): shape = (shape,) - if self.shape == (): - q = np.broadcast_to(self.quaternion,shape+(4,)) - else: - q = np.block([np.broadcast_to(self.quaternion[...,0:1],shape).reshape(shape+(1,)), - np.broadcast_to(self.quaternion[...,1:2],shape).reshape(shape+(1,)), - np.broadcast_to(self.quaternion[...,2:3],shape).reshape(shape+(1,)), - np.broadcast_to(self.quaternion[...,3:4],shape).reshape(shape+(1,))]) - return self.__class__(q) - - - def average(self,other): #ToDo: discuss calling for vectors - """ - Calculate the average rotation. - - Parameters - ---------- - other : Rotation - Rotation from which the average is rotated. - - """ - if self.quaternion.shape != (4,) or other.quaternion.shape != (4,): - raise NotImplementedError('Support for multiple rotations missing') - return Rotation.from_average([self,other]) - - ################################################################################################ # convert to different orientation representations (numpy arrays) @@ -326,20 +421,6 @@ class Rotation: """ return Rotation._qu2cu(self.quaternion) - @property - def M(self): # ToDo not sure about the name: as_M or M? we do not have a from_M - """ - Intermediate representation supporting quaternion averaging. - - References - ---------- - F. Landis Markley et al., Journal of Guidance, Control, and Dynamics 30(4):1193-1197, 2007 - https://doi.org/10.2514/1.28949 - - """ - return np.einsum('...i,...j',self.quaternion,self.quaternion) - - ################################################################################################ # Static constructors. The input data needs to follow the conventions, options allow to # relax the conventions. @@ -347,7 +428,7 @@ class Rotation: def from_quaternion(q, accept_homomorph = False, P = -1, - acceptHomomorph = None): # old name (for compatibility) + **kwargs): """ Initialize from quaternion. @@ -363,15 +444,13 @@ class Rotation: Convention used. Defaults to -1. """ - if acceptHomomorph is not None: - accept_homomorph = acceptHomomorph # for compatibility qu = np.array(q,dtype=float) if qu.shape[:-2:-1] != (4,): raise ValueError('Invalid shape.') if abs(P) != 1: raise ValueError('P ∉ {-1,1}') - if P == 1: qu[...,1:4] *= -1 + qu[...,1:4] *= -P if accept_homomorph: qu[qu[...,0] < 0.0] *= -1 else: @@ -384,7 +463,8 @@ class Rotation: @staticmethod def from_Eulers(phi, - degrees = False): + degrees = False, + **kwargs): """ Initialize from Bunge-Euler angles. @@ -411,7 +491,8 @@ class Rotation: def from_axis_angle(axis_angle, degrees = False, normalize = False, - P = -1): + P = -1, + **kwargs): """ Initialize from Axis angle pair. @@ -434,7 +515,7 @@ class Rotation: if abs(P) != 1: raise ValueError('P ∉ {-1,1}') - if P == 1: ax[...,0:3] *= -1 + ax[...,0:3] *= -P if degrees: ax[..., 3] = np.radians(ax[...,3]) if normalize: ax[...,0:3] /= np.linalg.norm(ax[...,0:3],axis=-1,keepdims=True) if np.any(ax[...,3] < 0.0) or np.any(ax[...,3] > np.pi): @@ -448,14 +529,15 @@ class Rotation: @staticmethod def from_basis(basis, orthonormal = True, - reciprocal = False): + reciprocal = False, + **kwargs): """ Initialize from lattice basis vectors. Parameters ---------- basis : numpy.ndarray of shape (...,3,3) - Three lattice basis vectors in three dimensions. + Three three-dimensional lattice basis vectors. orthonormal : boolean, optional Basis is strictly orthonormal, i.e. is free of stretch components. Defaults to True. reciprocal : boolean, optional @@ -463,7 +545,7 @@ class Rotation: """ om = np.array(basis,dtype=float) - if om.shape[:-3:-1] != (3,3): + if om.shape[-2:] != (3,3): raise ValueError('Invalid shape.') if reciprocal: @@ -482,7 +564,7 @@ class Rotation: return Rotation(Rotation._om2qu(om)) @staticmethod - def from_matrix(R): + def from_matrix(R,**kwargs): """ Initialize from rotation matrix. @@ -494,10 +576,40 @@ class Rotation: """ return Rotation.from_basis(R) + @staticmethod + def from_parallel(a,b, + **kwargs): + """ + Initialize from pairs of two orthogonal lattice basis vectors. + + Parameters + ---------- + a : numpy.ndarray of shape (...,2,3) + Two three-dimensional lattice vectors of first orthogonal basis. + b : numpy.ndarray of shape (...,2,3) + Corresponding three-dimensional lattice vectors of second basis. + + """ + a_ = np.array(a) + b_ = np.array(b) + if a_.shape[-2:] != (2,3) or b_.shape[-2:] != (2,3) or a_.shape != b_.shape: + raise ValueError('Invalid shape.') + am = np.stack([ a_[...,0,:], + a_[...,1,:], + np.cross(a_[...,0,:],a_[...,1,:]) ],axis=-2) + bm = np.stack([ b_[...,0,:], + b_[...,1,:], + np.cross(b_[...,0,:],b_[...,1,:]) ],axis=-2) + + return Rotation.from_basis(np.swapaxes(am/np.linalg.norm(am,axis=-1,keepdims=True),-1,-2))\ + .misorientation(Rotation.from_basis(np.swapaxes(bm/np.linalg.norm(bm,axis=-1,keepdims=True),-1,-2))) + + @staticmethod def from_Rodrigues(rho, normalize = False, - P = -1): + P = -1, + **kwargs): """ Initialize from Rodrigues-Frank vector. @@ -518,7 +630,7 @@ class Rotation: if abs(P) != 1: raise ValueError('P ∉ {-1,1}') - if P == 1: ro[...,0:3] *= -1 + ro[...,0:3] *= -P if normalize: ro[...,0:3] /= np.linalg.norm(ro[...,0:3],axis=-1,keepdims=True) if np.any(ro[...,3] < 0.0): raise ValueError('Rodrigues vector rotation angle not positive.') @@ -529,7 +641,8 @@ class Rotation: @staticmethod def from_homochoric(h, - P = -1): + P = -1, + **kwargs): """ Initialize from homochoric vector. @@ -547,7 +660,7 @@ class Rotation: if abs(P) != 1: raise ValueError('P ∉ {-1,1}') - if P == 1: ho *= -1 + ho *= -P if np.any(np.linalg.norm(ho,axis=-1) >_R1+1e-9): raise ValueError('Homochoric coordinate outside of the sphere.') @@ -556,7 +669,8 @@ class Rotation: @staticmethod def from_cubochoric(c, - P = -1): + P = -1, + **kwargs): """ Initialize from cubochoric vector. @@ -577,46 +691,15 @@ class Rotation: if np.abs(np.max(cu)) > np.pi**(2./3.) * 0.5+1e-9: raise ValueError('Cubochoric coordinate outside of the cube.') - ho = Rotation._cu2ho(cu) - if P == 1: ho *= -1 + ho = -P * Rotation._cu2ho(cu) return Rotation(Rotation._ho2qu(ho)) @staticmethod - def from_average(rotations,weights = None): - """ - Average rotation. - - References - ---------- - F. Landis Markley et al., Journal of Guidance, Control, and Dynamics 30(4):1193-1197, 2007 - https://doi.org/10.2514/1.28949 - - Parameters - ---------- - rotations : list of Rotations - Rotations to average from - weights : list of floats, optional - Weights for each rotation used for averaging - - """ - if not all(isinstance(item, Rotation) for item in rotations): - raise TypeError('Only instances of Rotation can be averaged.') - - N = len(rotations) - if not weights: - weights = np.ones(N,dtype='i') - - for i,(r,n) in enumerate(zip(rotations,weights)): - M = r.M * n if i == 0 \ - else M + r.M * n # noqa add (multiples) of this rotation to average noqa - eig, vec = np.linalg.eig(M/N) - - return Rotation.from_quaternion(np.real(vec.T[eig.argmax()]),accept_homomorph = True) - - @staticmethod - def from_random(shape=None,seed=None): + def from_random(shape = None, + seed = None, + **kwargs): """ Draw random rotation. @@ -633,12 +716,7 @@ class Rotation: """ rng = np.random.default_rng(seed) - if shape is None: - r = rng.random(3) - elif hasattr(shape, '__iter__'): - r = rng.random(tuple(shape)+(3,)) - else: - r = rng.random((shape,3)) + r = rng.random(3 if shape is None else tuple(shape)+(3,) if hasattr(shape, '__iter__') else (shape,3)) A = np.sqrt(r[...,2]) B = np.sqrt(1.0-r[...,2]) @@ -647,14 +725,17 @@ class Rotation: np.cos(2.0*np.pi*r[...,1])*B, np.sin(2.0*np.pi*r[...,0])*A],axis=-1) - return Rotation(q.reshape(r.shape[:-1]+(4,)) if shape is not None else q)._standardize() - - # for compatibility - __mul__ = __matmul__ + return Rotation(q if shape is None else q.reshape(r.shape[:-1]+(4,)))._standardize() @staticmethod - def from_ODF(weights,Eulers,N=500,degrees=True,fractions=True,seed=None): + def from_ODF(weights, + Eulers, + N = 500, + degrees = True, + fractions = True, + seed = None, + **kwargs): """ Sample discrete values from a binned ODF. @@ -707,7 +788,12 @@ class Rotation: @staticmethod - def from_spherical_component(center,sigma,N=500,degrees=True,seed=None): + def from_spherical_component(center, + sigma, + N = 500, + degrees = True, + seed = None, + **kwargs): """ Calculate set of rotations with Gaussian distribution around center. @@ -738,7 +824,13 @@ class Rotation: @staticmethod - def from_fiber_component(alpha,beta,sigma=0.0,N=500,degrees=True,seed=None): + def from_fiber_component(alpha, + beta, + sigma = 0.0, + N = 500, + degrees = True, + seed = None, + **kwargs): """ Calculate set of rotations with Gaussian distribution around direction. diff --git a/python/damask/_table.py b/python/damask/_table.py index dd8df97b0..e90650eea 100644 --- a/python/damask/_table.py +++ b/python/damask/_table.py @@ -175,7 +175,7 @@ class Table: @property def labels(self): - return list(self.shapes.keys()) + return list(self.shapes) def get(self,label): diff --git a/python/damask/lattice.py b/python/damask/lattice.py new file mode 100644 index 000000000..768dc5ace --- /dev/null +++ b/python/damask/lattice.py @@ -0,0 +1,420 @@ +import numpy as _np + +kinematics = { + 'cF': { + 'slip' : _np.array([ + [+0,+1,-1 , +1,+1,+1], + [-1,+0,+1 , +1,+1,+1], + [+1,-1,+0 , +1,+1,+1], + [+0,-1,-1 , -1,-1,+1], + [+1,+0,+1 , -1,-1,+1], + [-1,+1,+0 , -1,-1,+1], + [+0,-1,+1 , +1,-1,-1], + [-1,+0,-1 , +1,-1,-1], + [+1,+1,+0 , +1,-1,-1], + [+0,+1,+1 , -1,+1,-1], + [+1,+0,-1 , -1,+1,-1], + [-1,-1,+0 , -1,+1,-1], + [+1,+1,+0 , +1,-1,+0], + [+1,-1,+0 , +1,+1,+0], + [+1,+0,+1 , +1,+0,-1], + [+1,+0,-1 , +1,+0,+1], + [+0,+1,+1 , +0,+1,-1], + [+0,+1,-1 , +0,+1,+1], + ],'d'), + 'twin' : _np.array([ + [-2, 1, 1, 1, 1, 1], + [ 1,-2, 1, 1, 1, 1], + [ 1, 1,-2, 1, 1, 1], + [ 2,-1, 1, -1,-1, 1], + [-1, 2, 1, -1,-1, 1], + [-1,-1,-2, -1,-1, 1], + [-2,-1,-1, 1,-1,-1], + [ 1, 2,-1, 1,-1,-1], + [ 1,-1, 2, 1,-1,-1], + [ 2, 1,-1, -1, 1,-1], + [-1,-2,-1, -1, 1,-1], + [-1, 1, 2, -1, 1,-1], + ],dtype=float), + }, + 'cI': { + 'slip' : _np.array([ + [+1,-1,+1 , +0,+1,+1], + [-1,-1,+1 , +0,+1,+1], + [+1,+1,+1 , +0,-1,+1], + [-1,+1,+1 , +0,-1,+1], + [-1,+1,+1 , +1,+0,+1], + [-1,-1,+1 , +1,+0,+1], + [+1,+1,+1 , -1,+0,+1], + [+1,-1,+1 , -1,+0,+1], + [-1,+1,+1 , +1,+1,+0], + [-1,+1,-1 , +1,+1,+0], + [+1,+1,+1 , -1,+1,+0], + [+1,+1,-1 , -1,+1,+0], + [-1,+1,+1 , +2,+1,+1], + [+1,+1,+1 , -2,+1,+1], + [+1,+1,-1 , +2,-1,+1], + [+1,-1,+1 , +2,+1,-1], + [+1,-1,+1 , +1,+2,+1], + [+1,+1,-1 , -1,+2,+1], + [+1,+1,+1 , +1,-2,+1], + [-1,+1,+1 , +1,+2,-1], + [+1,+1,-1 , +1,+1,+2], + [+1,-1,+1 , -1,+1,+2], + [-1,+1,+1 , +1,-1,+2], + [+1,+1,+1 , +1,+1,-2], + ],'d'), + 'twin' : _np.array([ + [-1, 1, 1, 2, 1, 1], + [ 1, 1, 1, -2, 1, 1], + [ 1, 1,-1, 2,-1, 1], + [ 1,-1, 1, 2, 1,-1], + [ 1,-1, 1, 1, 2, 1], + [ 1, 1,-1, -1, 2, 1], + [ 1, 1, 1, 1,-2, 1], + [-1, 1, 1, 1, 2,-1], + [ 1, 1,-1, 1, 1, 2], + [ 1,-1, 1, -1, 1, 2], + [-1, 1, 1, 1,-1, 2], + [ 1, 1, 1, 1, 1,-2], + ],dtype=float), + }, + 'hP': { + 'slip' : _np.array([ + [+2,-1,-1,+0 , +0,+0,+0,+1], + [-1,+2,-1,+0 , +0,+0,+0,+1], + [-1,-1,+2,+0 , +0,+0,+0,+1], + [+2,-1,-1,+0 , +0,+1,-1,+0], + [-1,+2,-1,+0 , -1,+0,+1,+0], + [-1,-1,+2,+0 , +1,-1,+0,+0], + [-1,+1,+0,+0 , +1,+1,-2,+0], + [+0,-1,+1,+0 , -2,+1,+1,+0], + [+1,+0,-1,+0 , +1,-2,+1,+0], + [-1,+2,-1,+0 , +1,+0,-1,+1], + [-2,+1,+1,+0 , +0,+1,-1,+1], + [-1,-1,+2,+0 , -1,+1,+0,+1], + [+1,-2,+1,+0 , -1,+0,+1,+1], + [+2,-1,-1,+0 , +0,-1,+1,+1], + [+1,+1,-2,+0 , +1,-1,+0,+1], + [-2,+1,+1,+3 , +1,+0,-1,+1], + [-1,-1,+2,+3 , +1,+0,-1,+1], + [-1,-1,+2,+3 , +0,+1,-1,+1], + [+1,-2,+1,+3 , +0,+1,-1,+1], + [+1,-2,+1,+3 , -1,+1,+0,+1], + [+2,-1,-1,+3 , -1,+1,+0,+1], + [+2,-1,-1,+3 , -1,+0,+1,+1], + [+1,+1,-2,+3 , -1,+0,+1,+1], + [+1,+1,-2,+3 , +0,-1,+1,+1], + [-1,+2,-1,+3 , +0,-1,+1,+1], + [-1,+2,-1,+3 , +1,-1,+0,+1], + [-2,+1,+1,+3 , +1,-1,+0,+1], + [-1,-1,+2,+3 , +1,+1,-2,+2], + [+1,-2,+1,+3 , -1,+2,-1,+2], + [+2,-1,-1,+3 , -2,+1,+1,+2], + [+1,+1,-2,+3 , -1,-1,+2,+2], + [-1,+2,-1,+3 , +1,-2,+1,+2], + [-2,+1,+1,+3 , +2,-1,-1,+2], + ],'d'), + 'twin' : _np.array([ + [-1, 0, 1, 1, 1, 0, -1, 2], # shear = (3-(c/a)^2)/(sqrt(3) c/a) <-10.1>{10.2} + [ 0, -1, 1, 1, 0, 1, -1, 2], + [ 1, -1, 0, 1, -1, 1, 0, 2], + [ 1, 0, -1, 1, -1, 0, 1, 2], + [ 0, 1, -1, 1, 0, -1, 1, 2], + [-1, 1, 0, 1, 1, -1, 0, 2], + [-1, -1, 2, 6, 1, 1, -2, 1], # shear = 1/(c/a) <11.6>{-1-1.1} + [ 1, -2, 1, 6, -1, 2, -1, 1], + [ 2, -1, -1, 6, -2, 1, 1, 1], + [ 1, 1, -2, 6, -1, -1, 2, 1], + [-1, 2, -1, 6, 1, -2, 1, 1], + [-2, 1, 1, 6, 2, -1, -1, 1], + [ 1, 0, -1, -2, 1, 0, -1, 1], # shear = (4(c/a)^2-9)/(4 sqrt(3) c/a) <10.-2>{10.1} + [ 0, 1, -1, -2, 0, 1, -1, 1], + [-1, 1, 0, -2, -1, 1, 0, 1], + [-1, 0, 1, -2, -1, 0, 1, 1], + [ 0, -1, 1, -2, 0, -1, 1, 1], + [ 1, -1, 0, -2, 1, -1, 0, 1], + [ 1, 1, -2, -3, 1, 1, -2, 2], # shear = 2((c/a)^2-2)/(3 c/a) <11.-3>{11.2} + [-1, 2, -1, -3, -1, 2, -1, 2], + [-2, 1, 1, -3, -2, 1, 1, 2], + [-1, -1, 2, -3, -1, -1, 2, 2], + [ 1, -2, 1, -3, 1, -2, 1, 2], + [ 2, -1, -1, -3, 2, -1, -1, 2], + ],dtype=float), + }, +} + +# Kurdjomov--Sachs orientation relationship for fcc <-> bcc transformation +# from S. Morito et al., Journal of Alloys and Compounds 577:s587-s592, 2013 +# also see K. Kitahara et al., Acta Materialia 54:1279-1288, 2006 + +relations = { + 'KS': { + 'cF' : _np.array([ + [[ -1, 0, 1],[ 1, 1, 1]], + [[ -1, 0, 1],[ 1, 1, 1]], + [[ 0, 1, -1],[ 1, 1, 1]], + [[ 0, 1, -1],[ 1, 1, 1]], + [[ 1, -1, 0],[ 1, 1, 1]], + [[ 1, -1, 0],[ 1, 1, 1]], + [[ 1, 0, -1],[ 1, -1, 1]], + [[ 1, 0, -1],[ 1, -1, 1]], + [[ -1, -1, 0],[ 1, -1, 1]], + [[ -1, -1, 0],[ 1, -1, 1]], + [[ 0, 1, 1],[ 1, -1, 1]], + [[ 0, 1, 1],[ 1, -1, 1]], + [[ 0, -1, 1],[ -1, 1, 1]], + [[ 0, -1, 1],[ -1, 1, 1]], + [[ -1, 0, -1],[ -1, 1, 1]], + [[ -1, 0, -1],[ -1, 1, 1]], + [[ 1, 1, 0],[ -1, 1, 1]], + [[ 1, 1, 0],[ -1, 1, 1]], + [[ -1, 1, 0],[ 1, 1, -1]], + [[ -1, 1, 0],[ 1, 1, -1]], + [[ 0, -1, -1],[ 1, 1, -1]], + [[ 0, -1, -1],[ 1, 1, -1]], + [[ 1, 0, 1],[ 1, 1, -1]], + [[ 1, 0, 1],[ 1, 1, -1]], + ],dtype=float), + 'cI' : _np.array([ + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + [[ -1, -1, 1],[ 0, 1, 1]], + [[ -1, 1, -1],[ 0, 1, 1]], + ],dtype=float), + }, + 'GT': { + 'cF' : _np.array([ + [[ -5,-12, 17],[ 1, 1, 1]], + [[ 17, -5,-12],[ 1, 1, 1]], + [[-12, 17, -5],[ 1, 1, 1]], + [[ 5, 12, 17],[ -1, -1, 1]], + [[-17, 5,-12],[ -1, -1, 1]], + [[ 12,-17, -5],[ -1, -1, 1]], + [[ -5, 12,-17],[ -1, 1, 1]], + [[ 17, 5, 12],[ -1, 1, 1]], + [[-12,-17, 5],[ -1, 1, 1]], + [[ 5,-12,-17],[ 1, -1, 1]], + [[-17, -5, 12],[ 1, -1, 1]], + [[ 12, 17, 5],[ 1, -1, 1]], + [[ -5, 17,-12],[ 1, 1, 1]], + [[-12, -5, 17],[ 1, 1, 1]], + [[ 17,-12, -5],[ 1, 1, 1]], + [[ 5,-17,-12],[ -1, -1, 1]], + [[ 12, 5, 17],[ -1, -1, 1]], + [[-17, 12, -5],[ -1, -1, 1]], + [[ -5,-17, 12],[ -1, 1, 1]], + [[-12, 5,-17],[ -1, 1, 1]], + [[ 17, 12, 5],[ -1, 1, 1]], + [[ 5, 17, 12],[ 1, -1, 1]], + [[ 12, -5,-17],[ 1, -1, 1]], + [[-17,-12, 5],[ 1, -1, 1]], + ],dtype=float), + 'cI' : _np.array([ + [[-17, -7, 17],[ 1, 0, 1]], + [[ 17,-17, -7],[ 1, 1, 0]], + [[ -7, 17,-17],[ 0, 1, 1]], + [[ 17, 7, 17],[ -1, 0, 1]], + [[-17, 17, -7],[ -1, -1, 0]], + [[ 7,-17,-17],[ 0, -1, 1]], + [[-17, 7,-17],[ -1, 0, 1]], + [[ 17, 17, 7],[ -1, 1, 0]], + [[ -7,-17, 17],[ 0, 1, 1]], + [[ 17, -7,-17],[ 1, 0, 1]], + [[-17,-17, 7],[ 1, -1, 0]], + [[ 7, 17, 17],[ 0, -1, 1]], + [[-17, 17, -7],[ 1, 1, 0]], + [[ -7,-17, 17],[ 0, 1, 1]], + [[ 17, -7,-17],[ 1, 0, 1]], + [[ 17,-17, -7],[ -1, -1, 0]], + [[ 7, 17, 17],[ 0, -1, 1]], + [[-17, 7,-17],[ -1, 0, 1]], + [[-17,-17, 7],[ -1, 1, 0]], + [[ -7, 17,-17],[ 0, 1, 1]], + [[ 17, 7, 17],[ -1, 0, 1]], + [[ 17, 17, 7],[ 1, -1, 0]], + [[ 7,-17,-17],[ 0, -1, 1]], + [[-17, -7, 17],[ 1, 0, 1]], + ],dtype=float), + }, + 'GT_prime': { + 'cF' : _np.array([ + [[ 0, 1, -1],[ 7, 17, 17]], + [[ -1, 0, 1],[ 17, 7, 17]], + [[ 1, -1, 0],[ 17, 17, 7]], + [[ 0, -1, -1],[ -7,-17, 17]], + [[ 1, 0, 1],[-17, -7, 17]], + [[ 1, -1, 0],[-17,-17, 7]], + [[ 0, 1, -1],[ 7,-17,-17]], + [[ 1, 0, 1],[ 17, -7,-17]], + [[ -1, -1, 0],[ 17,-17, -7]], + [[ 0, -1, -1],[ -7, 17,-17]], + [[ -1, 0, 1],[-17, 7,-17]], + [[ -1, -1, 0],[-17, 17, -7]], + [[ 0, -1, 1],[ 7, 17, 17]], + [[ 1, 0, -1],[ 17, 7, 17]], + [[ -1, 1, 0],[ 17, 17, 7]], + [[ 0, 1, 1],[ -7,-17, 17]], + [[ -1, 0, -1],[-17, -7, 17]], + [[ -1, 1, 0],[-17,-17, 7]], + [[ 0, -1, 1],[ 7,-17,-17]], + [[ -1, 0, -1],[ 17, -7,-17]], + [[ 1, 1, 0],[ 17,-17, -7]], + [[ 0, 1, 1],[ -7, 17,-17]], + [[ 1, 0, -1],[-17, 7,-17]], + [[ 1, 1, 0],[-17, 17, -7]], + ],dtype=float), + 'cI' : _np.array([ + [[ 1, 1, -1],[ 12, 5, 17]], + [[ -1, 1, 1],[ 17, 12, 5]], + [[ 1, -1, 1],[ 5, 17, 12]], + [[ -1, -1, -1],[-12, -5, 17]], + [[ 1, -1, 1],[-17,-12, 5]], + [[ 1, -1, -1],[ -5,-17, 12]], + [[ -1, 1, -1],[ 12, -5,-17]], + [[ 1, 1, 1],[ 17,-12, -5]], + [[ -1, -1, 1],[ 5,-17,-12]], + [[ 1, -1, -1],[-12, 5,-17]], + [[ -1, -1, 1],[-17, 12, -5]], + [[ -1, -1, -1],[ -5, 17,-12]], + [[ 1, -1, 1],[ 12, 17, 5]], + [[ 1, 1, -1],[ 5, 12, 17]], + [[ -1, 1, 1],[ 17, 5, 12]], + [[ -1, 1, 1],[-12,-17, 5]], + [[ -1, -1, -1],[ -5,-12, 17]], + [[ -1, 1, -1],[-17, -5, 12]], + [[ -1, -1, 1],[ 12,-17, -5]], + [[ -1, 1, -1],[ 5,-12,-17]], + [[ 1, 1, 1],[ 17, -5,-12]], + [[ 1, 1, 1],[-12, 17, -5]], + [[ 1, -1, -1],[ -5, 12,-17]], + [[ 1, 1, -1],[-17, 5,-12]], + ],dtype=float), + }, + 'NW': { + 'cF' : _np.array([ + [[ 2, -1, -1],[ 1, 1, 1]], + [[ -1, 2, -1],[ 1, 1, 1]], + [[ -1, -1, 2],[ 1, 1, 1]], + [[ -2, -1, -1],[ -1, 1, 1]], + [[ 1, 2, -1],[ -1, 1, 1]], + [[ 1, -1, 2],[ -1, 1, 1]], + [[ 2, 1, -1],[ 1, -1, 1]], + [[ -1, -2, -1],[ 1, -1, 1]], + [[ -1, 1, 2],[ 1, -1, 1]], + [[ 2, -1, 1],[ -1, -1, 1]], + [[ -1, 2, 1],[ -1, -1, 1]], + [[ -1, -1, -2],[ -1, -1, 1]], + ],dtype=float), + 'cI' : _np.array([ + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + [[ 0, -1, 1],[ 0, 1, 1]], + ],dtype=float), + }, + 'Pitsch': { + 'cF' : _np.array([ + [[ 1, 0, 1],[ 0, 1, 0]], + [[ 1, 1, 0],[ 0, 0, 1]], + [[ 0, 1, 1],[ 1, 0, 0]], + [[ 0, 1, -1],[ 1, 0, 0]], + [[ -1, 0, 1],[ 0, 1, 0]], + [[ 1, -1, 0],[ 0, 0, 1]], + [[ 1, 0, -1],[ 0, 1, 0]], + [[ -1, 1, 0],[ 0, 0, 1]], + [[ 0, -1, 1],[ 1, 0, 0]], + [[ 0, 1, 1],[ 1, 0, 0]], + [[ 1, 0, 1],[ 0, 1, 0]], + [[ 1, 1, 0],[ 0, 0, 1]], + ],dtype=float), + 'cI' : _np.array([ + [[ 1, -1, 1],[ -1, 0, 1]], + [[ 1, 1, -1],[ 1, -1, 0]], + [[ -1, 1, 1],[ 0, 1, -1]], + [[ -1, 1, -1],[ 0, -1, -1]], + [[ -1, -1, 1],[ -1, 0, -1]], + [[ 1, -1, -1],[ -1, -1, 0]], + [[ 1, -1, -1],[ -1, 0, -1]], + [[ -1, 1, -1],[ -1, -1, 0]], + [[ -1, -1, 1],[ 0, -1, -1]], + [[ -1, 1, 1],[ 0, -1, 1]], + [[ 1, -1, 1],[ 1, 0, -1]], + [[ 1, 1, -1],[ -1, 1, 0]], + ],dtype=float), + }, + 'Bain': { + 'cF' : _np.array([ + [[ 0, 1, 0],[ 1, 0, 0]], + [[ 0, 0, 1],[ 0, 1, 0]], + [[ 1, 0, 0],[ 0, 0, 1]], + ],dtype=float), + 'cI' : _np.array([ + [[ 0, 1, 1],[ 1, 0, 0]], + [[ 1, 0, 1],[ 0, 1, 0]], + [[ 1, 1, 0],[ 0, 0, 1]], + ],dtype=float), + }, + 'Burgers' : { + 'cI' : _np.array([ + [[ -1, 1, 1],[ 1, 1, 0]], + [[ -1, 1, -1],[ 1, 1, 0]], + [[ 1, 1, 1],[ 1, -1, 0]], + [[ 1, 1, -1],[ 1, -1, 0]], + + [[ 1, 1, -1],[ 1, 0, 1]], + [[ -1, 1, 1],[ 1, 0, 1]], + [[ 1, 1, 1],[ -1, 0, 1]], + [[ 1, -1, 1],[ -1, 0, 1]], + + [[ -1, 1, -1],[ 0, 1, 1]], + [[ 1, 1, -1],[ 0, 1, 1]], + [[ -1, 1, 1],[ 0, -1, 1]], + [[ 1, 1, 1],[ 0, -1, 1]], + ],dtype=float), + 'hP' : _np.array([ + [[ -1, 2, -1, 0],[ 0, 0, 0, 1]], + [[ -1, -1, 2, 0],[ 0, 0, 0, 1]], + [[ -1, 2, -1, 0],[ 0, 0, 0, 1]], + [[ -1, -1, 2, 0],[ 0, 0, 0, 1]], + + [[ -1, 2, -1, 0],[ 0, 0, 0, 1]], + [[ -1, -1, 2, 0],[ 0, 0, 0, 1]], + [[ -1, 2, -1, 0],[ 0, 0, 0, 1]], + [[ -1, -1, 2, 0],[ 0, 0, 0, 1]], + + [[ -1, 2, -1, 0],[ 0, 0, 0, 1]], + [[ -1, -1, 2, 0],[ 0, 0, 0, 1]], + [[ -1, 2, -1, 0],[ 0, 0, 0, 1]], + [[ -1, -1, 2, 0],[ 0, 0, 0, 1]], + ],dtype=float), + }, +} diff --git a/python/damask/util.py b/python/damask/util.py index 914c3951c..d02877f8f 100644 --- a/python/damask/util.py +++ b/python/damask/util.py @@ -3,6 +3,7 @@ import datetime import os import subprocess import shlex +import re import fractions from functools import reduce from optparse import Option @@ -20,10 +21,13 @@ __all__=[ 'execute', 'show_progress', 'scale_to_coprime', + 'project_stereographic', 'hybrid_IA', 'return_message', 'extendableOption', - 'execution_stamp' + 'execution_stamp', + 'shapeshifter', + 'shapeblender', ] #################################################################################################### @@ -182,6 +186,28 @@ def scale_to_coprime(v): return m +def project_stereographic(vector,normalize=False): + """ + Apply stereographic projection to vector. + + Parameters + ---------- + vector : numpy.ndarray of shape (...,3) + Vector coordinates to be projected. + normalize : bool + Ensure unit length for vector. Defaults to False. + + Returns + ------- + coordinates : numpy.ndarray of shape (...,2) + Projected coordinates. + + """ + v_ = vector/np.linalg.norm(vector,axis=-1,keepdims=True) if normalize else vector + return np.block([v_[...,:2]/(1+np.abs(v_[...,2:3])), + np.zeros_like(v_[...,2:3])]) + + def execution_stamp(class_name,function_name=None): """Timestamp the execution of a (function within a) class.""" now = datetime.datetime.now().astimezone().strftime('%Y-%m-%d %H:%M:%S%z') @@ -203,6 +229,77 @@ def hybrid_IA(dist,N,seed=None): return np.repeat(np.arange(len(dist)),repeats)[np.random.default_rng(seed).permutation(N_inv_samples)[:N]] +def shapeshifter(fro,to,mode='left',keep_ones=False): + """ + Return a tuple that reshapes 'fro' to become broadcastable to 'to'. + + Parameters + ---------- + fro : tuple + Original shape of array. + to : tuple + Target shape of array after broadcasting. + len(to) cannot be less than len(fro). + mode : str, optional + Indicates whether new axes are preferably added to + either 'left' or 'right' of the original shape. + Defaults to 'left'. + keep_ones : bool, optional + Treat '1' in fro as literal value instead of dimensional placeholder. + Defaults to False. + + """ + beg = dict(left ='(^.*\\b)', + right='(^.*?\\b)') + sep = dict(left ='(.*\\b)', + right='(.*?\\b)') + end = dict(left ='(.*?$)', + right='(.*$)') + fro = (1,) if not len(fro) else fro + to = (1,) if not len(to) else to + try: + grp = re.match(beg[mode] + +f',{sep[mode]}'.join(map(lambda x: f'{x}' + if x>1 or (keep_ones and len(fro)>1) else + '\\d+',fro)) + +f',{end[mode]}', + ','.join(map(str,to))+',').groups() + except AttributeError: + raise ValueError(f'Shapes can not be shifted {fro} --> {to}') + fill = () + for g,d in zip(grp,fro+(None,)): + fill += (1,)*g.count(',')+(d,) + return fill[:-1] + + +def shapeblender(a,b): + """ + Return a shape that overlaps the rightmost entries of 'a' with the leftmost of 'b'. + + Parameters + ---------- + a : tuple + Shape of first array. + b : tuple + Shape of second array. + + Examples + -------- + >>> shapeblender((4,4,3),(3,2,1)) + (4,4,3,2,1) + >>> shapeblender((1,2),(1,2,3)) + (1,2,3) + >>> shapeblender((1,),(2,2,1)) + (1,2,2,1) + >>> shapeblender((3,2),(3,2)) + (3,2) + + """ + i = min(len(a),len(b)) + while i > 0 and a[-i:] != b[:i]: i -= 1 + return a + b[i:] + + #################################################################################################### # Classes #################################################################################################### diff --git a/python/tests/reference/ConfigMaterial/material.yaml b/python/tests/reference/ConfigMaterial/material.yaml index 97c6504bb..b616e1ef8 100644 --- a/python/tests/reference/ConfigMaterial/material.yaml +++ b/python/tests/reference/ConfigMaterial/material.yaml @@ -36,9 +36,9 @@ phase: elasticity: {C_11: 106.75e9, C_12: 60.41e9, C_44: 28.34e9, type: hooke} generic: output: [F, P, Fe, Fp, Lp] - lattice: fcc + lattice: cF Steel: elasticity: {C_11: 233.3e9, C_12: 135.5e9, C_44: 118.0e9, type: hooke} generic: output: [F, P, Fe, Fp, Lp] - lattice: bcc + lattice: cI diff --git a/python/tests/reference/Orientation/cF_Bain.txt b/python/tests/reference/Orientation/cF_Bain.txt new file mode 100644 index 000000000..876cf3888 --- /dev/null +++ b/python/tests/reference/Orientation/cF_Bain.txt @@ -0,0 +1,5 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +180.0 45.00000000000001 180.0 1 1 +270.0 45.00000000000001 90.0 1 2 +315.0 0.0 0.0 1 3 diff --git a/python/tests/reference/Orientation/cF_GT.txt b/python/tests/reference/Orientation/cF_GT.txt new file mode 100644 index 000000000..cefae431a --- /dev/null +++ b/python/tests/reference/Orientation/cF_GT.txt @@ -0,0 +1,26 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +146.75362934444064 9.976439066337804 256.395594327347 1 1 +356.59977719102034 43.39784965440254 12.173896584899929 1 2 +75.92521636876346 43.82007387041961 277.8843642946069 1 3 +326.75362934444064 9.976439066337806 76.39559432734703 1 4 +176.59977719102034 43.397849654402556 192.17389658489986 1 5 +255.92521636876344 43.82007387041961 97.88436429460687 1 6 +213.24637065555936 9.976439066337804 103.604405672653 1 7 +3.400222808979685 43.39784965440255 347.8261034151001 1 8 +284.0747836312365 43.82007387041961 82.11563570539313 1 9 +33.24637065555936 9.976439066337804 283.60440567265294 1 10 +183.40022280897963 43.397849654402556 167.8261034151001 1 11 +104.07478363123654 43.82007387041961 262.1156357053931 1 12 +273.4002228089796 43.397849654402556 77.82610341510008 1 13 +123.24637065555939 9.976439066337806 193.60440567265297 1 14 +194.07478363123653 43.82007387041961 172.11563570539317 1 15 +93.40022280897969 43.39784965440255 257.8261034151001 1 16 +303.24637065555936 9.976439066337804 13.604405672652977 1 17 +14.074783631236542 43.82007387041961 352.1156357053931 1 18 +86.59977719102032 43.39784965440254 282.17389658489986 1 19 +236.75362934444058 9.976439066337804 166.39559432734703 1 20 +165.92521636876344 43.82007387041961 187.88436429460683 1 21 +266.59977719102034 43.39784965440254 102.17389658489992 1 22 +56.75362934444064 9.976439066337804 346.395594327347 1 23 +345.9252163687635 43.82007387041961 7.884364294606862 1 24 diff --git a/python/tests/reference/Orientation/cF_GT_prime.txt b/python/tests/reference/Orientation/cF_GT_prime.txt new file mode 100644 index 000000000..44a9b25ec --- /dev/null +++ b/python/tests/reference/Orientation/cF_GT_prime.txt @@ -0,0 +1,26 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +166.39559432734697 9.976439066337804 236.75362934444058 1 1 +352.1156357053931 43.82007387041961 14.074783631236542 1 2 +77.82610341510008 43.397849654402556 273.4002228089796 1 3 +346.395594327347 9.976439066337804 56.75362934444064 1 4 +172.11563570539317 43.82007387041961 194.07478363123653 1 5 +257.8261034151001 43.39784965440255 93.40022280897969 1 6 +193.604405672653 9.976439066337804 123.24637065555939 1 7 +7.884364294606862 43.82007387041961 345.9252163687635 1 8 +282.17389658489986 43.39784965440254 86.59977719102032 1 9 +13.604405672652977 9.976439066337804 303.24637065555936 1 10 +187.88436429460683 43.82007387041961 165.92521636876344 1 11 +102.17389658489992 43.39784965440254 266.59977719102034 1 12 +277.8843642946069 43.82007387041961 75.92521636876346 1 13 +103.604405672653 9.976439066337804 213.24637065555936 1 14 +192.17389658489986 43.397849654402556 176.59977719102034 1 15 +97.88436429460687 43.82007387041961 255.92521636876344 1 16 +283.60440567265294 9.976439066337804 33.24637065555936 1 17 +12.173896584899929 43.39784965440254 356.59977719102034 1 18 +82.11563570539313 43.82007387041961 284.0747836312365 1 19 +256.395594327347 9.976439066337804 146.75362934444064 1 20 +167.8261034151001 43.397849654402556 183.40022280897963 1 21 +262.1156357053931 43.82007387041961 104.07478363123654 1 22 +76.39559432734703 9.976439066337806 326.75362934444064 1 23 +347.8261034151001 43.39784965440255 3.400222808979685 1 24 diff --git a/python/tests/reference/Orientation/cF_KS.txt b/python/tests/reference/Orientation/cF_KS.txt new file mode 100644 index 000000000..93fdcf07e --- /dev/null +++ b/python/tests/reference/Orientation/cF_KS.txt @@ -0,0 +1,26 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +114.20342833932975 10.52877936550932 204.20342833932972 1 1 +94.3573968784815 80.40593177313954 311.22729452432543 1 2 +175.6426031215185 80.40593177313954 48.77270547567447 1 3 +155.79657166067025 10.52877936550932 155.79657166067025 1 4 +99.62136089109411 85.70366403943004 318.04510841542015 1 5 +170.37863910890587 85.70366403943002 41.954891584579855 1 6 +85.64260312151852 80.40593177313954 48.77270547567448 1 7 +65.79657166067024 10.52877936550932 155.79657166067025 1 8 +9.621360891094124 85.70366403943004 318.04510841542015 1 9 +80.37863910890587 85.70366403943004 41.95489158457987 1 10 +24.203428339329758 10.52877936550932 204.20342833932975 1 11 +4.357396878481486 80.40593177313954 311.2272945243255 1 12 +204.20342833932972 10.52877936550932 204.20342833932972 1 13 +184.35739687848147 80.40593177313954 311.2272945243255 1 14 +265.64260312151845 80.40593177313953 48.77270547567449 1 15 +245.79657166067025 10.528779365509317 155.79657166067025 1 16 +189.62136089109413 85.70366403943004 318.04510841542015 1 17 +260.3786391089059 85.70366403943002 41.954891584579855 1 18 +170.37863910890587 94.29633596056996 138.04510841542015 1 19 +99.62136089109411 94.29633596056998 221.95489158457983 1 20 +155.79657166067025 169.4712206344907 24.203428339329754 1 21 +175.64260312151848 99.59406822686046 131.22729452432552 1 22 +94.35739687848151 99.59406822686046 228.77270547567446 1 23 +114.20342833932975 169.4712206344907 335.7965716606702 1 24 diff --git a/python/tests/reference/Orientation/cF_NW.txt b/python/tests/reference/Orientation/cF_NW.txt new file mode 100644 index 000000000..cc9c95a05 --- /dev/null +++ b/python/tests/reference/Orientation/cF_NW.txt @@ -0,0 +1,14 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +96.91733794010702 83.13253115922213 314.5844440567886 1 1 +173.082662059893 83.13253115922211 45.41555594321143 1 2 +135.0 9.735610317245317 180.0 1 3 +263.082662059893 83.13253115922213 45.415555943211444 1 4 +186.91733794010702 83.13253115922211 314.5844440567886 1 5 +224.99999999999997 9.735610317245317 180.0 1 6 +83.082662059893 83.13253115922213 45.415555943211444 1 7 +6.917337940106983 83.13253115922211 314.5844440567886 1 8 +45.0 9.73561031724532 180.0 1 9 +13.638707279476469 45.81931182053557 80.40196970123216 1 10 +256.36129272052347 45.81931182053556 279.59803029876775 1 11 +315.0 99.73561031724536 0.0 1 12 diff --git a/python/tests/reference/Orientation/cF_Pitsch.txt b/python/tests/reference/Orientation/cF_Pitsch.txt new file mode 100644 index 000000000..aa0c32365 --- /dev/null +++ b/python/tests/reference/Orientation/cF_Pitsch.txt @@ -0,0 +1,14 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +135.41555594321144 83.13253115922213 173.082662059893 1 1 +260.26438968275465 90.0 135.0 1 2 +260.40196970123213 45.81931182053557 13.638707279476478 1 3 +314.5844440567886 83.13253115922213 96.91733794010702 1 4 +350.40196970123213 45.81931182053557 283.6387072794765 1 5 +170.26438968275465 90.0 224.99999999999997 1 6 +315.4155559432114 83.13253115922213 353.08266205989304 1 7 +99.73561031724536 90.0 225.0 1 8 +279.59803029876787 45.819311820535574 166.36129272052352 1 9 +134.58444405678856 83.13253115922213 276.91733794010696 1 10 +9.598030298767851 45.819311820535574 76.36129272052355 1 11 +9.735610317245369 90.0 315.0 1 12 diff --git a/python/tests/reference/Orientation/cF_slip.txt b/python/tests/reference/Orientation/cF_slip.txt new file mode 100644 index 000000000..18aa03f24 --- /dev/null +++ b/python/tests/reference/Orientation/cF_slip.txt @@ -0,0 +1,19 @@ +3x3:1_Schmid 3x3:2_Schmid 3x3:3_Schmid 3x3:4_Schmid 3x3:5_Schmid 3x3:6_Schmid 3x3:7_Schmid 3x3:8_Schmid 3x3:9_Schmid +0.0 0.0 0.0 0.408248290463863 0.408248290463863 0.40824829046386296 -0.4082482904638631 -0.4082482904638631 -0.408248290463863 +-0.408248290463863 -0.408248290463863 -0.40824829046386296 2.4997998108697446e-17 2.4997998108697446e-17 2.499799810869744e-17 0.4082482904638631 0.4082482904638631 0.408248290463863 +0.408248290463863 0.408248290463863 0.40824829046386296 -0.4082482904638631 -0.4082482904638631 -0.408248290463863 0.0 0.0 0.0 +4.999599621739488e-17 4.9995996217394874e-17 -4.999599621739488e-17 0.408248290463863 0.40824829046386296 -0.408248290463863 0.408248290463863 0.40824829046386296 -0.408248290463863 +-0.408248290463863 -0.40824829046386296 0.408248290463863 -2.499799810869744e-17 -2.4997998108697437e-17 2.499799810869744e-17 -0.408248290463863 -0.40824829046386296 0.408248290463863 +0.408248290463863 0.40824829046386296 -0.408248290463863 -0.4082482904638631 -0.408248290463863 0.4082482904638631 0.0 0.0 0.0 +0.0 0.0 0.0 -0.408248290463863 0.408248290463863 0.408248290463863 0.4082482904638631 -0.4082482904638631 -0.4082482904638631 +-0.408248290463863 0.408248290463863 0.408248290463863 -2.499799810869744e-17 2.499799810869744e-17 2.499799810869744e-17 -0.408248290463863 0.408248290463863 0.408248290463863 +0.408248290463863 -0.408248290463863 -0.408248290463863 0.408248290463863 -0.408248290463863 -0.408248290463863 0.0 0.0 0.0 +-4.999599621739488e-17 4.999599621739488e-17 -4.9995996217394874e-17 -0.408248290463863 0.408248290463863 -0.40824829046386296 -0.408248290463863 0.408248290463863 -0.40824829046386296 +-0.408248290463863 0.408248290463863 -0.40824829046386296 2.4997998108697446e-17 -2.4997998108697446e-17 2.499799810869744e-17 0.4082482904638631 -0.4082482904638631 0.408248290463863 +0.408248290463863 -0.408248290463863 0.40824829046386296 0.408248290463863 -0.408248290463863 0.40824829046386296 0.0 0.0 0.0 +0.4999999999999999 -0.4999999999999999 0.0 0.4999999999999999 -0.4999999999999999 0.0 0.0 0.0 0.0 +0.5 0.4999999999999999 -6.123233995736766e-17 -0.5000000000000001 -0.5 6.123233995736767e-17 0.0 0.0 0.0 +0.4999999999999999 -3.0616169978683824e-17 -0.4999999999999999 3.0616169978683824e-17 -1.874699728327322e-33 -3.0616169978683824e-17 0.4999999999999999 -3.0616169978683824e-17 -0.4999999999999999 +0.5 -3.061616997868383e-17 0.4999999999999999 -3.0616169978683836e-17 1.8746997283273227e-33 -3.061616997868383e-17 -0.5000000000000001 3.0616169978683836e-17 -0.5 +0.0 6.123233995736765e-17 -6.123233995736765e-17 0.0 0.4999999999999999 -0.4999999999999999 0.0 0.4999999999999999 -0.4999999999999999 +0.0 0.0 0.0 0.0 0.5 0.4999999999999999 0.0 -0.5000000000000001 -0.5 diff --git a/python/tests/reference/Orientation/cF_twin.txt b/python/tests/reference/Orientation/cF_twin.txt new file mode 100644 index 000000000..8b0961823 --- /dev/null +++ b/python/tests/reference/Orientation/cF_twin.txt @@ -0,0 +1,13 @@ +3x3:1_Schmid 3x3:2_Schmid 3x3:3_Schmid 3x3:4_Schmid 3x3:5_Schmid 3x3:6_Schmid 3x3:7_Schmid 3x3:8_Schmid 3x3:9_Schmid +-0.4714045207910318 -0.4714045207910318 -0.47140452079103173 0.2357022603955159 0.2357022603955159 0.23570226039551587 0.2357022603955159 0.2357022603955159 0.23570226039551587 +0.2357022603955159 0.2357022603955159 0.23570226039551587 -0.4714045207910318 -0.4714045207910318 -0.47140452079103173 0.2357022603955159 0.2357022603955159 0.23570226039551587 +0.23570226039551587 0.23570226039551587 0.23570226039551584 0.23570226039551587 0.23570226039551587 0.23570226039551584 -0.4714045207910318 -0.4714045207910318 -0.47140452079103173 +-0.4714045207910318 -0.47140452079103173 0.4714045207910318 0.23570226039551587 0.23570226039551584 -0.23570226039551587 -0.2357022603955159 -0.23570226039551587 0.2357022603955159 +0.23570226039551584 0.23570226039551578 -0.23570226039551584 -0.4714045207910318 -0.47140452079103173 0.4714045207910318 -0.2357022603955159 -0.23570226039551587 0.2357022603955159 +0.2357022603955159 0.23570226039551587 -0.2357022603955159 0.2357022603955159 0.23570226039551587 -0.2357022603955159 0.4714045207910317 0.47140452079103157 -0.4714045207910317 +-0.4714045207910318 0.4714045207910318 0.4714045207910318 -0.2357022603955159 0.2357022603955159 0.2357022603955159 -0.2357022603955159 0.2357022603955159 0.2357022603955159 +0.23570226039551595 -0.23570226039551595 -0.23570226039551595 0.4714045207910318 -0.4714045207910318 -0.4714045207910318 -0.2357022603955159 0.2357022603955159 0.2357022603955159 +0.2357022603955159 -0.2357022603955159 -0.2357022603955159 -0.23570226039551587 0.23570226039551587 0.23570226039551587 0.4714045207910318 -0.4714045207910318 -0.4714045207910318 +-0.4714045207910318 0.4714045207910318 -0.47140452079103173 -0.23570226039551587 0.23570226039551587 -0.23570226039551584 0.2357022603955159 -0.2357022603955159 0.23570226039551587 +0.23570226039551595 -0.23570226039551595 0.2357022603955159 0.4714045207910318 -0.4714045207910318 0.47140452079103173 0.2357022603955159 -0.2357022603955159 0.23570226039551587 +0.23570226039551584 -0.23570226039551584 0.23570226039551578 -0.23570226039551595 0.23570226039551595 -0.2357022603955159 -0.4714045207910318 0.4714045207910318 -0.47140452079103173 diff --git a/python/tests/reference/Orientation/cI_Bain.txt b/python/tests/reference/Orientation/cI_Bain.txt new file mode 100644 index 000000000..e0bc4f6c7 --- /dev/null +++ b/python/tests/reference/Orientation/cI_Bain.txt @@ -0,0 +1,5 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +0.0 45.00000000000001 0.0 1 1 +90.0 45.00000000000001 270.0 1 2 +45.00000000000001 0.0 0.0 1 3 diff --git a/python/tests/reference/Orientation/cI_GT.txt b/python/tests/reference/Orientation/cI_GT.txt new file mode 100644 index 000000000..5d5102698 --- /dev/null +++ b/python/tests/reference/Orientation/cI_GT.txt @@ -0,0 +1,26 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +283.60440567265294 9.976439066337804 33.24637065555936 1 1 +167.8261034151001 43.397849654402556 183.40022280897963 1 2 +262.1156357053931 43.82007387041961 104.07478363123654 1 3 +103.604405672653 9.976439066337804 213.24637065555936 1 4 +347.8261034151001 43.39784965440255 3.400222808979685 1 5 +82.11563570539313 43.82007387041961 284.0747836312365 1 6 +76.39559432734703 9.976439066337806 326.75362934444064 1 7 +192.17389658489986 43.397849654402556 176.59977719102034 1 8 +97.88436429460687 43.82007387041961 255.92521636876344 1 9 +256.395594327347 9.976439066337804 146.75362934444064 1 10 +12.173896584899929 43.39784965440254 356.59977719102034 1 11 +277.8843642946069 43.82007387041961 75.92521636876346 1 12 +102.17389658489992 43.39784965440254 266.59977719102034 1 13 +346.395594327347 9.976439066337804 56.75362934444064 1 14 +7.884364294606862 43.82007387041961 345.9252163687635 1 15 +282.17389658489986 43.39784965440254 86.59977719102032 1 16 +166.39559432734703 9.976439066337804 236.75362934444058 1 17 +187.88436429460683 43.82007387041961 165.92521636876344 1 18 +257.8261034151001 43.39784965440255 93.40022280897969 1 19 +13.604405672652977 9.976439066337804 303.24637065555936 1 20 +352.1156357053931 43.82007387041961 14.074783631236542 1 21 +77.82610341510008 43.397849654402556 273.4002228089796 1 22 +193.60440567265297 9.976439066337806 123.24637065555939 1 23 +172.11563570539317 43.82007387041961 194.07478363123653 1 24 diff --git a/python/tests/reference/Orientation/cI_GT_prime.txt b/python/tests/reference/Orientation/cI_GT_prime.txt new file mode 100644 index 000000000..e398d3139 --- /dev/null +++ b/python/tests/reference/Orientation/cI_GT_prime.txt @@ -0,0 +1,26 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +303.24637065555936 9.976439066337804 13.604405672652977 1 1 +165.92521636876344 43.82007387041961 187.88436429460683 1 2 +266.59977719102034 43.39784965440254 102.17389658489992 1 3 +123.24637065555939 9.976439066337804 193.604405672653 1 4 +345.9252163687635 43.82007387041961 7.884364294606862 1 5 +86.59977719102032 43.39784965440254 282.17389658489986 1 6 +56.75362934444064 9.976439066337804 346.395594327347 1 7 +194.07478363123653 43.82007387041961 172.11563570539317 1 8 +93.40022280897969 43.39784965440255 257.8261034151001 1 9 +236.75362934444058 9.976439066337804 166.39559432734697 1 10 +14.074783631236542 43.82007387041961 352.1156357053931 1 11 +273.4002228089796 43.397849654402556 77.82610341510008 1 12 +104.07478363123654 43.82007387041961 262.1156357053931 1 13 +326.75362934444064 9.976439066337806 76.39559432734703 1 14 +3.400222808979685 43.39784965440255 347.8261034151001 1 15 +284.0747836312365 43.82007387041961 82.11563570539313 1 16 +146.75362934444064 9.976439066337804 256.395594327347 1 17 +183.40022280897963 43.397849654402556 167.8261034151001 1 18 +255.92521636876344 43.82007387041961 97.88436429460687 1 19 +33.24637065555936 9.976439066337804 283.60440567265294 1 20 +356.59977719102034 43.39784965440254 12.173896584899929 1 21 +75.92521636876346 43.82007387041961 277.8843642946069 1 22 +213.24637065555936 9.976439066337804 103.604405672653 1 23 +176.59977719102034 43.397849654402556 192.17389658489986 1 24 diff --git a/python/tests/reference/Orientation/cI_KS.txt b/python/tests/reference/Orientation/cI_KS.txt new file mode 100644 index 000000000..34b393358 --- /dev/null +++ b/python/tests/reference/Orientation/cI_KS.txt @@ -0,0 +1,26 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +335.7965716606702 10.528779365509317 65.79657166067024 1 1 +228.77270547567446 80.40593177313953 85.64260312151849 1 2 +131.22729452432552 80.40593177313954 4.357396878481506 1 3 +24.20342833932977 10.52877936550932 24.20342833932976 1 4 +221.95489158457983 85.70366403943002 80.37863910890589 1 5 +138.04510841542015 85.70366403943004 9.621360891094124 1 6 +131.22729452432552 80.40593177313953 94.35739687848151 1 7 +24.203428339329765 10.52877936550932 114.20342833932976 1 8 +221.95489158457983 85.70366403943004 170.37863910890587 1 9 +138.04510841542015 85.70366403943004 99.62136089109411 1 10 +335.7965716606702 10.52877936550932 155.79657166067025 1 11 +228.77270547567448 80.40593177313954 175.6426031215185 1 12 +335.7965716606702 10.52877936550932 335.7965716606702 1 13 +228.77270547567448 80.40593177313954 355.6426031215185 1 14 +131.2272945243255 80.40593177313954 274.35739687848144 1 15 +24.203428339329747 10.52877936550932 294.2034283393298 1 16 +221.95489158457985 85.70366403943004 350.3786391089059 1 17 +138.04510841542015 85.70366403943004 279.6213608910941 1 18 +41.95489158457986 94.29633596056998 9.621360891094133 1 19 +318.04510841542015 94.29633596056996 80.37863910890589 1 20 +155.79657166067025 169.4712206344907 24.203428339329754 1 21 +48.77270547567448 99.59406822686046 4.357396878481504 1 22 +311.2272945243255 99.59406822686046 85.64260312151852 1 23 +204.20342833932975 169.4712206344907 65.79657166067024 1 24 diff --git a/python/tests/reference/Orientation/cI_NW.txt b/python/tests/reference/Orientation/cI_NW.txt new file mode 100644 index 000000000..754c69bba --- /dev/null +++ b/python/tests/reference/Orientation/cI_NW.txt @@ -0,0 +1,14 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +225.41555594321144 83.13253115922213 83.08266205989301 1 1 +134.58444405678856 83.13253115922211 6.917337940107012 1 2 +4.702125169424418e-15 9.735610317245317 45.0 1 3 +134.58444405678856 83.13253115922213 276.91733794010696 1 4 +225.4155559432114 83.13253115922213 353.082662059893 1 5 +0.0 9.735610317245317 315.0 1 6 +134.58444405678858 83.13253115922213 96.91733794010702 1 7 +225.41555594321142 83.13253115922213 173.082662059893 1 8 +0.0 9.735610317245317 135.0 1 9 +99.59803029876785 45.81931182053557 166.36129272052355 1 10 +260.40196970123213 45.81931182053556 283.6387072794765 1 11 +180.0 99.73561031724535 225.0 1 12 diff --git a/python/tests/reference/Orientation/cI_Pitsch.txt b/python/tests/reference/Orientation/cI_Pitsch.txt new file mode 100644 index 000000000..ef28bbb4d --- /dev/null +++ b/python/tests/reference/Orientation/cI_Pitsch.txt @@ -0,0 +1,14 @@ +1 header +1_Eulers 2_Eulers 3_Eulers 1_pos 2_pos +6.9173379401070045 83.13253115922213 44.58444405678856 1 1 +45.0 89.99999999999999 279.7356103172453 1 2 +166.36129272052352 45.819311820535574 279.59803029876787 1 3 +83.08266205989301 83.13253115922213 225.41555594321144 1 4 +256.3612927205235 45.819311820535574 189.59803029876787 1 5 +315.0 90.0 9.735610317245369 1 6 +186.917337940107 83.13253115922213 224.58444405678856 1 7 +315.0 90.0 80.26438968275463 1 8 +13.638707279476478 45.81931182053557 260.40196970123213 1 9 +263.082662059893 83.13253115922213 45.415555943211444 1 10 +103.63870727947646 45.819311820535574 170.40196970123213 1 11 +224.99999999999997 90.0 170.26438968275465 1 12 diff --git a/python/tests/reference/Orientation/cI_slip.txt b/python/tests/reference/Orientation/cI_slip.txt new file mode 100644 index 000000000..5401392c8 --- /dev/null +++ b/python/tests/reference/Orientation/cI_slip.txt @@ -0,0 +1,25 @@ +3x3:1_Schmid 3x3:2_Schmid 3x3:3_Schmid 3x3:4_Schmid 3x3:5_Schmid 3x3:6_Schmid 3x3:7_Schmid 3x3:8_Schmid 3x3:9_Schmid +0.0 0.4082482904638631 0.408248290463863 0.0 -0.408248290463863 -0.40824829046386296 0.0 0.4082482904638631 0.408248290463863 +0.0 -0.408248290463863 -0.40824829046386296 0.0 -0.408248290463863 -0.40824829046386296 0.0 0.4082482904638631 0.408248290463863 +0.0 -0.408248290463863 0.408248290463863 0.0 -0.408248290463863 0.408248290463863 0.0 -0.408248290463863 0.408248290463863 +0.0 0.40824829046386285 -0.40824829046386285 0.0 -0.408248290463863 0.408248290463863 0.0 -0.408248290463863 0.408248290463863 +-0.40824829046386296 2.4997998108697434e-17 -0.40824829046386285 0.4082482904638631 -2.4997998108697446e-17 0.408248290463863 0.4082482904638631 -2.4997998108697446e-17 0.408248290463863 +-0.408248290463863 2.499799810869744e-17 -0.40824829046386296 -0.408248290463863 2.499799810869744e-17 -0.40824829046386296 0.4082482904638631 -2.4997998108697446e-17 0.408248290463863 +-0.408248290463863 2.499799810869744e-17 0.408248290463863 -0.408248290463863 2.499799810869744e-17 0.408248290463863 -0.408248290463863 2.499799810869744e-17 0.408248290463863 +-0.408248290463863 2.499799810869744e-17 0.408248290463863 0.40824829046386296 -2.4997998108697437e-17 -0.40824829046386296 -0.408248290463863 2.499799810869744e-17 0.408248290463863 +-0.40824829046386296 -0.40824829046386285 4.999599621739487e-17 0.4082482904638631 0.408248290463863 -4.999599621739489e-17 0.4082482904638631 0.408248290463863 -4.999599621739489e-17 +-0.4082482904638631 -0.408248290463863 4.999599621739489e-17 0.408248290463863 0.40824829046386296 -4.999599621739488e-17 -0.4082482904638631 -0.408248290463863 4.999599621739489e-17 +-0.408248290463863 0.408248290463863 0.0 -0.408248290463863 0.408248290463863 0.0 -0.408248290463863 0.408248290463863 0.0 +-0.40824829046386296 0.40824829046386296 0.0 -0.40824829046386296 0.40824829046386296 0.0 0.408248290463863 -0.408248290463863 0.0 +-0.4714045207910316 -0.23570226039551578 -0.23570226039551576 0.4714045207910318 0.23570226039551587 0.23570226039551584 0.4714045207910318 0.23570226039551587 0.23570226039551584 +-0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 +0.47140452079103173 -0.2357022603955159 0.23570226039551584 0.47140452079103173 -0.2357022603955159 0.23570226039551584 -0.4714045207910318 0.23570226039551595 -0.23570226039551587 +0.4714045207910318 0.23570226039551587 -0.23570226039551595 -0.47140452079103173 -0.23570226039551584 0.2357022603955159 0.4714045207910318 0.23570226039551587 -0.23570226039551595 +0.2357022603955159 0.4714045207910318 0.23570226039551584 -0.23570226039551587 -0.47140452079103173 -0.23570226039551578 0.2357022603955159 0.4714045207910318 0.23570226039551584 +-0.23570226039551587 0.47140452079103173 0.23570226039551587 -0.23570226039551587 0.47140452079103173 0.23570226039551587 0.2357022603955159 -0.4714045207910318 -0.2357022603955159 +0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 +-0.2357022603955158 -0.4714045207910316 0.23570226039551584 0.2357022603955159 0.4714045207910318 -0.23570226039551595 0.2357022603955159 0.4714045207910318 -0.23570226039551595 +0.23570226039551587 0.23570226039551584 0.47140452079103173 0.23570226039551587 0.23570226039551584 0.47140452079103173 -0.2357022603955159 -0.23570226039551587 -0.4714045207910318 +-0.2357022603955159 0.2357022603955159 0.4714045207910318 0.23570226039551587 -0.23570226039551587 -0.47140452079103173 -0.2357022603955159 0.2357022603955159 0.4714045207910318 +-0.2357022603955158 0.2357022603955158 -0.4714045207910316 0.2357022603955159 -0.2357022603955159 0.4714045207910318 0.2357022603955159 -0.2357022603955159 0.4714045207910318 +0.2357022603955159 0.23570226039551587 -0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 diff --git a/python/tests/reference/Orientation/cI_twin.txt b/python/tests/reference/Orientation/cI_twin.txt new file mode 100644 index 000000000..accbd02d8 --- /dev/null +++ b/python/tests/reference/Orientation/cI_twin.txt @@ -0,0 +1,13 @@ +3x3:1_Schmid 3x3:2_Schmid 3x3:3_Schmid 3x3:4_Schmid 3x3:5_Schmid 3x3:6_Schmid 3x3:7_Schmid 3x3:8_Schmid 3x3:9_Schmid +-0.4714045207910316 -0.23570226039551578 -0.23570226039551576 0.4714045207910318 0.23570226039551587 0.23570226039551584 0.4714045207910318 0.23570226039551587 0.23570226039551584 +-0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 +0.47140452079103173 -0.2357022603955159 0.23570226039551584 0.47140452079103173 -0.2357022603955159 0.23570226039551584 -0.4714045207910318 0.23570226039551595 -0.23570226039551587 +0.4714045207910318 0.23570226039551587 -0.23570226039551595 -0.47140452079103173 -0.23570226039551584 0.2357022603955159 0.4714045207910318 0.23570226039551587 -0.23570226039551595 +0.2357022603955159 0.4714045207910318 0.23570226039551584 -0.23570226039551587 -0.47140452079103173 -0.23570226039551578 0.2357022603955159 0.4714045207910318 0.23570226039551584 +-0.23570226039551587 0.47140452079103173 0.23570226039551587 -0.23570226039551587 0.47140452079103173 0.23570226039551587 0.2357022603955159 -0.4714045207910318 -0.2357022603955159 +0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 0.2357022603955159 -0.4714045207910318 0.23570226039551595 +-0.2357022603955158 -0.4714045207910316 0.23570226039551584 0.2357022603955159 0.4714045207910318 -0.23570226039551595 0.2357022603955159 0.4714045207910318 -0.23570226039551595 +0.23570226039551587 0.23570226039551584 0.47140452079103173 0.23570226039551587 0.23570226039551584 0.47140452079103173 -0.2357022603955159 -0.23570226039551587 -0.4714045207910318 +-0.2357022603955159 0.2357022603955159 0.4714045207910318 0.23570226039551587 -0.23570226039551587 -0.47140452079103173 -0.2357022603955159 0.2357022603955159 0.4714045207910318 +-0.2357022603955158 0.2357022603955158 -0.4714045207910316 0.2357022603955159 -0.2357022603955159 0.4714045207910318 0.2357022603955159 -0.2357022603955159 0.4714045207910318 +0.2357022603955159 0.23570226039551587 -0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 0.2357022603955159 0.23570226039551587 -0.4714045207910318 diff --git a/python/tests/reference/Orientation/hP_slip.txt b/python/tests/reference/Orientation/hP_slip.txt new file mode 100644 index 000000000..a67eea150 --- /dev/null +++ b/python/tests/reference/Orientation/hP_slip.txt @@ -0,0 +1,34 @@ +3x3:1_Schmid 3x3:2_Schmid 3x3:3_Schmid 3x3:4_Schmid 3x3:5_Schmid 3x3:6_Schmid 3x3:7_Schmid 3x3:8_Schmid 3x3:9_Schmid +0.0 0.0 1.0 0.0 0.0 0.0 0.0 0.0 0.0 +0.0 0.0 -0.4999999999999998 0.0 0.0 0.8660254037844388 0.0 0.0 0.0 +0.0 0.0 -0.5000000000000001 0.0 0.0 -0.8660254037844386 0.0 0.0 0.0 +0.0 1.0 -5.914589856893347e-17 0.0 0.0 0.0 0.0 0.0 0.0 +0.4330127018922192 0.24999999999999975 -3.102315069664884e-17 -0.7500000000000002 -0.4330127018922192 5.373367321746164e-17 0.0 0.0 0.0 +-0.43301270189221935 0.25000000000000006 1.4502014121821253e-18 -0.7499999999999998 0.4330127018922193 2.5118225271075755e-18 0.0 0.0 0.0 +-0.4330127018922194 -0.7499999999999999 6.059609998111558e-17 0.2500000000000001 0.4330127018922194 -3.498517463593857e-17 0.0 0.0 0.0 +2.563950248511418e-16 -5.693113199781536e-32 -9.614043519462407e-33 1.0 -2.220446049250313e-16 -3.7496997163046135e-17 0.0 0.0 0.0 +0.4330127018922194 -0.75 2.8122747872284606e-17 0.25000000000000006 -0.43301270189221935 1.6236676054415494e-17 0.0 0.0 0.0 +-0.38254602783800284 -0.22086305214969287 -0.23426064283290896 0.6625891564490795 0.38254602783800284 0.40575133560034454 0.0 0.0 0.0 +0.0 -0.8834522085987724 -0.4685212856658182 0.0 0.0 0.0 0.0 0.0 0.0 +0.38254602783800307 -0.22086305214969315 -0.23426064283290912 0.6625891564490792 -0.38254602783800296 -0.40575133560034443 0.0 0.0 0.0 +-0.38254602783800284 -0.22086305214969287 0.23426064283290904 0.6625891564490795 0.38254602783800284 -0.4057513356003447 0.0 0.0 0.0 +0.0 -0.8834522085987724 0.46852128566581835 0.0 0.0 0.0 0.0 0.0 0.0 +0.38254602783800307 -0.22086305214969315 0.23426064283290912 0.6625891564490792 -0.38254602783800296 0.40575133560034443 0.0 0.0 0.0 +-0.39955629492721617 -0.23068393443263763 -0.24467726152216654 3.8591083978971935e-17 2.228057272357889e-17 2.3632116092343994e-17 0.6524726973924442 0.3767052874784087 0.39955629492721606 +-0.19977814746360817 -0.11534196721631886 -0.12233863076108333 -0.34602590164895664 -0.19977814746360797 -0.21189672420660496 0.6524726973924442 0.3767052874784087 0.39955629492721606 +0.0 -0.23068393443263785 -0.12233863076108334 0.0 -0.39955629492721617 -0.211896724206605 0.0 0.7534105749568177 0.3995562949272161 +0.0 0.23068393443263768 0.12233863076108326 0.0 -0.3995562949272162 -0.21189672420660502 0.0 0.7534105749568178 0.39955629492721617 +-0.199778147463608 0.11534196721631884 0.12233863076108324 0.34602590164895664 -0.1997781474636081 -0.211896724206605 -0.6524726973924442 0.3767052874784089 0.3995562949272161 +-0.3995562949272161 0.23068393443263774 0.24467726152216654 -3.859108397897193e-17 2.22805727235789e-17 2.3632116092343994e-17 -0.6524726973924441 0.37670528747840887 0.39955629492721606 +-0.39955629492721617 -0.23068393443263763 0.24467726152216662 -3.8591083978971935e-17 -2.228057272357889e-17 2.3632116092344003e-17 -0.6524726973924442 -0.3767052874784087 0.39955629492721617 +-0.1997781474636082 -0.11534196721631888 0.1223386307610834 -0.3460259016489568 -0.19977814746360806 0.21189672420660513 -0.6524726973924442 -0.3767052874784087 0.39955629492721617 +0.0 -0.23068393443263788 0.12233863076108341 0.0 -0.3995562949272163 0.21189672420660516 0.0 -0.7534105749568177 0.3995562949272162 +0.0 0.2306839344326376 -0.12233863076108324 0.0 -0.3995562949272163 0.21189672420660516 0.0 -0.7534105749568177 0.3995562949272162 +-0.19977814746360792 0.1153419672163188 -0.12233863076108319 0.34602590164895675 -0.19977814746360814 0.21189672420660505 0.6524726973924441 -0.37670528747840887 0.39955629492721606 +-0.3995562949272161 0.23068393443263774 -0.24467726152216654 3.859108397897193e-17 -2.22805727235789e-17 2.3632116092343994e-17 0.6524726973924441 -0.37670528747840887 0.39955629492721606 +-0.11134044285378089 -0.19284730395996755 -0.1363636363636364 -0.19284730395996755 -0.3340213285613424 -0.23618874648666507 0.36363636363636365 0.6298366572977734 0.44536177141512323 +-0.11134044285378081 0.1928473039599675 0.13636363636363633 0.19284730395996758 -0.3340213285613426 -0.2361887464866651 -0.3636363636363637 0.6298366572977737 0.4453617714151233 +-0.44536177141512323 9.889017858258314e-17 0.2727272727272727 -4.301519895922435e-17 9.551292858672588e-33 2.634132215859942e-17 -0.7272727272727272 1.6148698540002275e-16 0.44536177141512323 +-0.1113404428537809 -0.1928473039599676 0.13636363636363644 -0.1928473039599676 -0.33402132856134253 0.23618874648666516 -0.36363636363636365 -0.6298366572977734 0.44536177141512323 +-0.11134044285378074 0.19284730395996735 -0.13636363636363627 0.19284730395996758 -0.33402132856134253 0.23618874648666516 0.3636363636363636 -0.6298366572977734 0.44536177141512323 +-0.44536177141512334 9.889017858258316e-17 -0.2727272727272727 4.3015198959224354e-17 -9.55129285867259e-33 2.634132215859942e-17 0.7272727272727273 -1.6148698540002277e-16 0.44536177141512323 diff --git a/python/tests/reference/Orientation/hP_twin.txt b/python/tests/reference/Orientation/hP_twin.txt new file mode 100644 index 000000000..62ce6fdef --- /dev/null +++ b/python/tests/reference/Orientation/hP_twin.txt @@ -0,0 +1,25 @@ +3x3:1_Schmid 3x3:2_Schmid 3x3:3_Schmid 3x3:4_Schmid 3x3:5_Schmid 3x3:6_Schmid 3x3:7_Schmid 3x3:8_Schmid 3x3:9_Schmid +-0.3743506488634663 -0.21613144789263322 -0.4584840372976439 -0.21613144789263333 -0.12478354962115536 -0.2647058823529411 0.4075413664867946 0.23529411764705865 0.4991341984846217 +0.0 -1.1032987950073291e-16 -1.1702250894369439e-16 0.0 -0.4991341984846217 -0.5294117647058824 0.0 0.47058823529411753 0.4991341984846218 +-0.3743506488634663 0.2161314478926334 0.4584840372976439 0.2161314478926334 -0.1247835496211555 -0.2647058823529412 -0.4075413664867947 0.23529411764705888 0.4991341984846218 +-0.3743506488634663 -0.21613144789263322 0.4584840372976439 -0.2161314478926334 -0.12478354962115538 0.2647058823529412 -0.4075413664867946 -0.23529411764705865 0.4991341984846217 +0.0 -1.4562117094830577e-16 1.5445457619280876e-16 0.0 -0.4991341984846217 0.5294117647058824 0.0 -0.47058823529411753 0.4991341984846218 +-0.3743506488634663 0.2161314478926334 -0.4584840372976439 0.21613144789263342 -0.12478354962115551 0.26470588235294124 0.4075413664867947 -0.23529411764705888 0.4991341984846218 +-0.06998542122237655 -0.1212183053462653 -0.04285714285714286 -0.12121830534626532 -0.20995626366712955 -0.07423074889580901 0.45714285714285724 0.7917946548886295 0.279941684889506 +-0.06998542122237653 0.12121830534626528 0.04285714285714286 0.12121830534626532 -0.20995626366712958 -0.07423074889580904 -0.45714285714285724 0.7917946548886297 0.2799416848895061 +-0.27994168488950616 6.2159540823338e-17 0.08571428571428573 -5.407625012016776e-17 1.2007339593759827e-32 1.6557402499691063e-17 -0.9142857142857143 2.0301221021717148e-16 0.27994168488950605 +-0.0699854212223766 -0.12121830534626538 0.04285714285714291 -0.12121830534626538 -0.20995626366712963 0.07423074889580909 -0.45714285714285724 -0.7917946548886295 0.2799416848895062 +-0.06998542122237648 0.12121830534626521 -0.04285714285714283 0.12121830534626538 -0.20995626366712966 0.07423074889580906 0.45714285714285724 -0.7917946548886297 0.2799416848895061 +-0.27994168488950605 6.215954082333798e-17 -0.08571428571428569 5.407625012016776e-17 -1.2007339593759827e-32 1.6557402499691063e-17 0.9142857142857143 -2.0301221021717148e-16 0.27994168488950605 +0.3104371234477526 0.17923095678901296 0.19010313741609627 0.17923095678901305 0.1034790411492508 0.1097560975609756 -0.6759222663683424 -0.39024390243902424 -0.41391616459700337 +0.0 7.68600963028337e-17 4.07612214737886e-17 0.0 0.4139161645970035 0.21951219512195125 0.0 -0.7804878048780488 -0.4139161645970034 +0.31043712344775254 -0.17923095678901305 -0.19010313741609627 -0.17923095678901302 0.10347904114925086 0.1097560975609756 0.6759222663683423 -0.3902439024390244 -0.41391616459700337 +0.3104371234477527 0.179230956789013 -0.19010313741609638 0.17923095678901313 0.10347904114925086 -0.10975609756097568 0.6759222663683424 0.3902439024390242 -0.4139161645970035 +0.0 1.3539199431344235e-16 -7.180244797305419e-17 0.0 0.4139161645970036 -0.21951219512195136 0.0 0.7804878048780487 -0.41391616459700353 +0.3104371234477525 -0.179230956789013 0.19010313741609622 -0.17923095678901313 0.10347904114925092 -0.10975609756097565 -0.6759222663683423 0.3902439024390244 -0.41391616459700337 +0.11134044285378089 0.19284730395996755 0.1363636363636364 0.19284730395996755 0.3340213285613424 0.23618874648666507 -0.36363636363636365 -0.6298366572977734 -0.44536177141512323 +0.11134044285378081 -0.1928473039599675 -0.13636363636363633 -0.19284730395996758 0.3340213285613426 0.2361887464866651 0.3636363636363637 -0.6298366572977737 -0.4453617714151233 +0.44536177141512323 -9.889017858258314e-17 -0.2727272727272727 4.301519895922435e-17 -9.551292858672588e-33 -2.634132215859942e-17 0.7272727272727272 -1.6148698540002275e-16 -0.44536177141512323 +0.1113404428537809 0.1928473039599676 -0.13636363636363644 0.1928473039599676 0.33402132856134253 -0.23618874648666516 0.36363636363636365 0.6298366572977734 -0.44536177141512323 +0.11134044285378074 -0.19284730395996735 0.13636363636363627 -0.19284730395996758 0.33402132856134253 -0.23618874648666516 -0.3636363636363636 0.6298366572977734 -0.44536177141512323 +0.44536177141512334 -9.889017858258316e-17 0.2727272727272727 -4.3015198959224354e-17 9.55129285867259e-33 -2.634132215859942e-17 -0.7272727272727273 1.6148698540002277e-16 -0.44536177141512323 diff --git a/python/tests/reference/Orientation/unitcell_cubic_0_0_0.pdf b/python/tests/reference/Orientation/unitcell_cubic_0_0_0.pdf new file mode 100644 index 000000000..f8ed27163 Binary files /dev/null and b/python/tests/reference/Orientation/unitcell_cubic_0_0_0.pdf differ diff --git a/python/tests/reference/Orientation/unitcell_cubic_0_45_0.pdf b/python/tests/reference/Orientation/unitcell_cubic_0_45_0.pdf new file mode 100644 index 000000000..1ebbb10d3 Binary files /dev/null and b/python/tests/reference/Orientation/unitcell_cubic_0_45_0.pdf differ diff --git a/python/tests/reference/Orientation/unitcell_cubic_180_45_180.pdf b/python/tests/reference/Orientation/unitcell_cubic_180_45_180.pdf new file mode 100644 index 000000000..9ecb9f44b Binary files /dev/null and b/python/tests/reference/Orientation/unitcell_cubic_180_45_180.pdf differ diff --git a/python/tests/reference/Orientation/unitcell_cubic_45_0_0.pdf b/python/tests/reference/Orientation/unitcell_cubic_45_0_0.pdf new file mode 100644 index 000000000..433a01bb2 Binary files /dev/null and b/python/tests/reference/Orientation/unitcell_cubic_45_0_0.pdf differ diff --git a/python/tests/reference/Orientation/unitcell_cubic_90_45_270.pdf b/python/tests/reference/Orientation/unitcell_cubic_90_45_270.pdf new file mode 100644 index 000000000..46d32c470 Binary files /dev/null and b/python/tests/reference/Orientation/unitcell_cubic_90_45_270.pdf differ diff --git a/python/tests/test_Lattice.py b/python/tests/test_Lattice.py deleted file mode 100644 index 93fc56445..000000000 --- a/python/tests/test_Lattice.py +++ /dev/null @@ -1,157 +0,0 @@ -import random - -import pytest -import numpy as np - -from damask import Rotation -from damask import Symmetry - -def in_FZ(system,rho): - """Non-vectorized version of 'in_FZ'.""" - rho_abs = abs(rho) - - if system == 'cubic': - return np.sqrt(2.0)-1.0 >= rho_abs[0] \ - and np.sqrt(2.0)-1.0 >= rho_abs[1] \ - and np.sqrt(2.0)-1.0 >= rho_abs[2] \ - and 1.0 >= rho_abs[0] + rho_abs[1] + rho_abs[2] - elif system == 'hexagonal': - return 1.0 >= rho_abs[0] and 1.0 >= rho_abs[1] and 1.0 >= rho_abs[2] \ - and 2.0 >= np.sqrt(3)*rho_abs[0] + rho_abs[1] \ - and 2.0 >= np.sqrt(3)*rho_abs[1] + rho_abs[0] \ - and 2.0 >= np.sqrt(3) + rho_abs[2] - elif system == 'tetragonal': - return 1.0 >= rho_abs[0] and 1.0 >= rho_abs[1] \ - and np.sqrt(2.0) >= rho_abs[0] + rho_abs[1] \ - and np.sqrt(2.0) >= rho_abs[2] + 1.0 - elif system == 'orthorhombic': - return 1.0 >= rho_abs[0] and 1.0 >= rho_abs[1] and 1.0 >= rho_abs[2] - else: - return np.all(np.isfinite(rho_abs)) - - -def in_disorientation_SST(system,rho): - """Non-vectorized version of 'in_Disorientation_SST'.""" - epsilon = 0.0 - if system == 'cubic': - return rho[0] >= rho[1]+epsilon and rho[1] >= rho[2]+epsilon and rho[2] >= epsilon - elif system == 'hexagonal': - return rho[0] >= np.sqrt(3)*(rho[1]-epsilon) and rho[1] >= epsilon and rho[2] >= epsilon - elif system == 'tetragonal': - return rho[0] >= rho[1]-epsilon and rho[1] >= epsilon and rho[2] >= epsilon - elif system == 'orthorhombic': - return rho[0] >= epsilon and rho[1] >= epsilon and rho[2] >= epsilon - else: - return True - - -def in_SST(system,vector,proper = False): - """Non-vectorized version of 'in_SST'.""" - if system == 'cubic': - basis = {'improper':np.array([ [-1. , 0. , 1. ], - [ np.sqrt(2.) , -np.sqrt(2.) , 0. ], - [ 0. , np.sqrt(3.) , 0. ] ]), - 'proper':np.array([ [ 0. , -1. , 1. ], - [-np.sqrt(2.) , np.sqrt(2.) , 0. ], - [ np.sqrt(3.) , 0. , 0. ] ]), - } - elif system == 'hexagonal': - basis = {'improper':np.array([ [ 0. , 0. , 1. ], - [ 1. , -np.sqrt(3.) , 0. ], - [ 0. , 2. , 0. ] ]), - 'proper':np.array([ [ 0. , 0. , 1. ], - [-1. , np.sqrt(3.) , 0. ], - [ np.sqrt(3.) , -1. , 0. ] ]), - } - elif system == 'tetragonal': - basis = {'improper':np.array([ [ 0. , 0. , 1. ], - [ 1. , -1. , 0. ], - [ 0. , np.sqrt(2.) , 0. ] ]), - 'proper':np.array([ [ 0. , 0. , 1. ], - [-1. , 1. , 0. ], - [ np.sqrt(2.) , 0. , 0. ] ]), - } - elif system == 'orthorhombic': - basis = {'improper':np.array([ [ 0., 0., 1.], - [ 1., 0., 0.], - [ 0., 1., 0.] ]), - 'proper':np.array([ [ 0., 0., 1.], - [-1., 0., 0.], - [ 0., 1., 0.] ]), - } - else: - return True - - v = np.array(vector,dtype=float) - if proper: - theComponents = np.around(np.dot(basis['improper'],v),12) - inSST = np.all(theComponents >= 0.0) - if not inSST: - theComponents = np.around(np.dot(basis['proper'],v),12) - inSST = np.all(theComponents >= 0.0) - else: - v[2] = abs(v[2]) - theComponents = np.around(np.dot(basis['improper'],v),12) - inSST = np.all(theComponents >= 0.0) - - return inSST - - -@pytest.fixture -def set_of_rodrigues(set_of_quaternions): - return Rotation(set_of_quaternions).as_Rodrigues(vector=True)[:200] - -class TestSymmetry: - - @pytest.mark.parametrize('system',Symmetry.crystal_systems) - def test_in_FZ_vectorize(self,set_of_rodrigues,system): - result = Symmetry(system).in_FZ(set_of_rodrigues.reshape(50,4,3)).reshape(200) - for i,r in enumerate(result): - assert r == in_FZ(system,set_of_rodrigues[i]) - - @pytest.mark.parametrize('system',Symmetry.crystal_systems) - def test_in_disorientation_SST_vectorize(self,set_of_rodrigues,system): - result = Symmetry(system).in_disorientation_SST(set_of_rodrigues.reshape(50,4,3)).reshape(200) - for i,r in enumerate(result): - assert r == in_disorientation_SST(system,set_of_rodrigues[i]) - - @pytest.mark.parametrize('proper',[True,False]) - @pytest.mark.parametrize('system',Symmetry.crystal_systems) - def test_in_SST_vectorize(self,system,proper): - vecs = np.random.rand(20,4,3) - result = Symmetry(system).in_SST(vecs,proper).reshape(20*4) - for i,r in enumerate(result): - assert r == in_SST(system,vecs.reshape(20*4,3)[i],proper) - - @pytest.mark.parametrize('invalid_symmetry',['fcc','bcc','hello']) - def test_invalid_symmetry(self,invalid_symmetry): - with pytest.raises(KeyError): - s = Symmetry(invalid_symmetry) # noqa - - def test_equal(self): - symmetry = random.choice(Symmetry.crystal_systems) - print(symmetry) - assert Symmetry(symmetry) == Symmetry(symmetry) - - def test_not_equal(self): - symmetries = random.sample(Symmetry.crystal_systems,k=2) - assert Symmetry(symmetries[0]) != Symmetry(symmetries[1]) - - @pytest.mark.parametrize('system',Symmetry.crystal_systems) - def test_in_FZ(self,system): - assert Symmetry(system).in_FZ(np.zeros(3)) - - @pytest.mark.parametrize('system',Symmetry.crystal_systems) - def test_in_disorientation_SST(self,system): - assert Symmetry(system).in_disorientation_SST(np.zeros(3)) - - @pytest.mark.parametrize('system',Symmetry.crystal_systems) - @pytest.mark.parametrize('proper',[True,False]) - def test_in_SST(self,system,proper): - assert Symmetry(system).in_SST(np.zeros(3),proper) - - @pytest.mark.parametrize('function',['in_FZ','in_disorientation_SST','in_SST']) - def test_invalid_argument(self,function): - s = Symmetry() # noqa - with pytest.raises(ValueError): - eval(f's.{function}(np.ones(4))') diff --git a/python/tests/test_Orientation.py b/python/tests/test_Orientation.py index 669f73e91..fe35bda14 100644 --- a/python/tests/test_Orientation.py +++ b/python/tests/test_Orientation.py @@ -1,131 +1,536 @@ -import os -from itertools import permutations - import pytest import numpy as np +from itertools import permutations -from damask import Table from damask import Rotation from damask import Orientation -from damask import Lattice - -n = 1000 - -def IPF_color(orientation,direction): - """TSL color of inverse pole figure for given axis (non-vectorized).""" - for o in orientation.equivalent: - pole = o.rotation@direction - inSST,color = orientation.lattice.in_SST(pole,color=True) - if inSST: break - - return color - -def inverse_pole(orientation,axis,proper=False,SST=True): - if SST: - for eq in orientation.equivalent: - pole = eq.rotation @ axis/np.linalg.norm(axis) - if orientation.lattice.in_SST(pole,proper=proper): - return pole - else: - return orientation.rotation @ axis/np.linalg.norm(axis) +from damask import Table +from damask import lattice +from damask import util @pytest.fixture def reference_dir(reference_dir_base): """Directory containing reference results.""" - return reference_dir_base/'Rotation' + return reference_dir_base/'Orientation' + +@pytest.fixture +def set_of_rodrigues(set_of_quaternions): + return Rotation(set_of_quaternions).as_Rodrigues()[:200] class TestOrientation: - @pytest.mark.parametrize('model',['Bain','KS','GT','GT_prime','NW','Pitsch']) - @pytest.mark.parametrize('lattice',['fcc','bcc']) - def test_relationship_vectorize(self,set_of_quaternions,lattice,model): - result = Orientation(set_of_quaternions[:200].reshape(50,4,4),lattice).related(model) - ref_qu = result.rotation.quaternion.reshape(-1,200,4) - for i in range(200): - single = Orientation(set_of_quaternions[i],lattice).related(model).rotation.quaternion - assert np.allclose(ref_qu[:,i,:],single) + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('shape',[None,5,(4,6)]) + def test_equal(self,lattice,shape): + R = Rotation.from_random(shape) + assert Orientation(R,lattice) == Orientation(R,lattice) - @pytest.mark.parametrize('lattice',Lattice.lattices) - def test_IPF_vectorize(self,set_of_quaternions,lattice): - direction = np.random.random(3)*2.0-1 - oris = Orientation(Rotation(set_of_quaternions),lattice)[:200] - for i,color in enumerate(oris.IPF_color(direction)): - assert np.allclose(color,IPF_color(oris[i],direction)) + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('shape',[None,5,(4,6)]) + def test_unequal(self,lattice,shape): + R = Rotation.from_random(shape) + assert not(Orientation(R,lattice) != Orientation(R,lattice)) - @pytest.mark.parametrize('SST',[False,True]) + @pytest.mark.parametrize('a,b',[ + (dict(rotation=[1,0,0,0]), + dict(rotation=[0.5,0.5,0.5,0.5])), + + (dict(rotation=[1,0,0,0],lattice='cubic'), + dict(rotation=[1,0,0,0],lattice='hexagonal')), + + (dict(rotation=[1,0,0,0],lattice='cF',a=1), + dict(rotation=[1,0,0,0],lattice='cF',a=2)), + ]) + def test_nonequal(self,a,b): + assert Orientation(**a) != Orientation(**b) + + @pytest.mark.parametrize('kwargs',[ + dict(lattice='aP', alpha=np.pi/4,beta=np.pi/3, ), + dict(lattice='mP', c=1.2,alpha=np.pi/4, gamma=np.pi/2), + dict(lattice='oP', c=1.2,alpha=np.pi/4, ), + dict(lattice='oS',a=1.0, c=2.0,alpha=np.pi/2,beta=np.pi/3, ), + dict(lattice='tP',a=1.0,b=1.2, ), + dict(lattice='tI', alpha=np.pi/3, ), + dict(lattice='hP', gamma=np.pi/2), + dict(lattice='cI',a=1.0, c=2.0,alpha=np.pi/2,beta=np.pi/2, ), + dict(lattice='cF', beta=np.pi/3, ), + ]) + def test_invalid_init(self,kwargs): + with pytest.raises(ValueError): + Orientation(**kwargs).parameters # noqa + + @pytest.mark.parametrize('kwargs',[ + dict(lattice='aP',a=1.0,b=1.1,c=1.2,alpha=np.pi/4,beta=np.pi/3,gamma=np.pi/2), + dict(lattice='mP',a=1.0,b=1.1,c=1.2, beta=np.pi/3 ), + dict(lattice='oS',a=1.0,b=1.1,c=1.2, ), + dict(lattice='tI',a=1.0, c=1.2, ), + dict(lattice='hP',a=1.0 ), + dict(lattice='cI',a=1.0, ), + ]) + def test_repr(self,kwargs): + o = Orientation.from_random(**kwargs) + assert isinstance(o.__repr__(),str) + + @pytest.mark.parametrize('kwargs',[ + dict(lattice='aP',a=1.0,b=1.1,c=1.2,alpha=np.pi/4,beta=np.pi/3,gamma=np.pi/2), + dict(lattice='mP',a=1.0,b=1.1,c=1.2, beta=np.pi/3 ), + dict(lattice='oS',a=1.0,b=1.1,c=1.2, ), + dict(lattice='tI',a=1.0, c=1.2, ), + dict(lattice='hP',a=1.0 ), + dict(lattice='cI',a=1.0, ), + ]) + def test_copy(self,kwargs): + o = Orientation.from_random(**kwargs) + p = o.copy(rotation=Rotation.from_random()) + assert o != p + + def test_from_quaternion(self): + assert np.all(Orientation.from_quaternion(q=np.array([1,0,0,0]),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_Eulers(self): + assert np.all(Orientation.from_Eulers(phi=np.zeros(3),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_axis_angle(self): + assert np.all(Orientation.from_axis_angle(axis_angle=[1,0,0,0],lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_basis(self): + assert np.all(Orientation.from_basis(basis=np.eye(3),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_matrix(self): + assert np.all(Orientation.from_matrix(R=np.eye(3),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_Rodrigues(self): + assert np.all(Orientation.from_Rodrigues(rho=np.array([0,0,1,0]),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_homochoric(self): + assert np.all(Orientation.from_homochoric(h=np.zeros(3),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_cubochoric(self): + assert np.all(Orientation.from_cubochoric(c=np.zeros(3),lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_spherical_component(self): + assert np.all(Orientation.from_spherical_component(center=Rotation(), + sigma=0.0,N=1,lattice='triclinic').as_matrix() + == np.eye(3)) + + def test_from_fiber_component(self): + r = Rotation.from_fiber_component(alpha=np.zeros(2),beta=np.zeros(2), + sigma=0.0,N=1,seed=0) + assert np.all(Orientation.from_fiber_component(alpha=np.zeros(2),beta=np.zeros(2), + sigma=0.0,N=1,seed=0,lattice='triclinic').quaternion + == r.quaternion) + + @pytest.mark.parametrize('kwargs',[ + dict(lattice='aP',a=1.0,b=1.1,c=1.2,alpha=np.pi/4.5,beta=np.pi/3.5,gamma=np.pi/2.5), + dict(lattice='mP',a=1.0,b=1.1,c=1.2, beta=np.pi/3.5), + dict(lattice='oS',a=1.0,b=1.1,c=1.2,), + dict(lattice='tI',a=1.0, c=1.2,), + dict(lattice='hP',a=1.0 ), + dict(lattice='cI',a=1.0, ), + ]) + def test_from_direction(self,kwargs): + for a,b in np.random.random((10,2,3)): + c = np.cross(b,a) + if np.all(np.isclose(c,0)): continue + o = Orientation.from_directions(uvw=a,hkl=c,**kwargs) + x = o.to_pole(uvw=a) + z = o.to_pole(hkl=c) + assert np.isclose(np.dot(x/np.linalg.norm(x),np.array([1,0,0])),1) \ + and np.isclose(np.dot(z/np.linalg.norm(z),np.array([0,0,1])),1) + + + def test_negative_angle(self): + with pytest.raises(ValueError): + Orientation(lattice='aP',a=1,b=2,c=3,alpha=45,beta=45,gamma=-45,degrees=True) # noqa + + def test_excess_angle(self): + with pytest.raises(ValueError): + Orientation(lattice='aP',a=1,b=2,c=3,alpha=45,beta=45,gamma=90.0001,degrees=True) # noqa + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('angle',[10,20,30,40]) + def test_average(self,angle,lattice): + o = Orientation.from_axis_angle(lattice=lattice,axis_angle=[[0,0,1,10],[0,0,1,angle]],degrees=True) + avg_angle = o.average().as_axis_angle(degrees=True,pair=True)[1] + assert np.isclose(avg_angle,10+(angle-10)/2.) + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + def test_reduced_equivalent(self,lattice): + i = Orientation(lattice=lattice) + o = Orientation.from_random(lattice=lattice) + eq = o.equivalent + FZ = np.argmin(abs(eq.misorientation(i.broadcast_to(len(eq))).as_axis_angle(pair=True)[1])) + assert o.reduced == eq[FZ] + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('N',[1,8,32]) + def test_disorientation(self,lattice,N): + o = Orientation.from_random(lattice=lattice,shape=N,seed=0) + p = Orientation.from_random(lattice=lattice,shape=N,seed=1) + + d,ops = o.disorientation(p,return_operators=True) + + for n in range(N): + assert np.allclose(d[n].as_quaternion(), + o[n].equivalent[ops[n][0]] + .misorientation(p[n].equivalent[ops[n][1]]) + .as_quaternion()) \ + or np.allclose((~d)[n].as_quaternion(), + o[n].equivalent[ops[n][0]] + .misorientation(p[n].equivalent[ops[n][1]]) + .as_quaternion()) + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('a,b',[ + ((2,3,2),(2,3,2)), + ((2,2),(4,4)), + ((3,1),(1,3)), + (None,None), + ]) + def test_disorientation_blending(self,lattice,a,b): + o = Orientation.from_random(lattice=lattice,shape=a,seed=0) + p = Orientation.from_random(lattice=lattice,shape=b,seed=1) + blend = util.shapeblender(o.shape,p.shape) + for loc in np.random.randint(0,blend,(10,len(blend))): + assert o[tuple(loc[:len(o.shape)])].disorientation(p[tuple(loc[-len(p.shape):])]) \ + == o.disorientation(p)[tuple(loc)] + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + def test_disorientation360(self,lattice): + o_1 = Orientation(Rotation(),lattice) + o_2 = Orientation.from_Eulers(lattice=lattice,phi=[360,0,0],degrees=True) + assert np.allclose((o_1.disorientation(o_2)).as_matrix(),np.eye(3)) + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) + def test_reduced_vectorization(self,lattice,shape): + o = Orientation.from_random(lattice=lattice,shape=shape,seed=0) + for r, theO in zip(o.reduced.flatten(),o.flatten()): + assert r == theO.reduced + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) + @pytest.mark.parametrize('vector',np.array([[1,0,0],[1,2,3],[-1,1,-1]])) @pytest.mark.parametrize('proper',[True,False]) - @pytest.mark.parametrize('lattice',Lattice.lattices) - def test_inverse_pole_vectorize(self,set_of_quaternions,lattice,SST,proper): - axis = np.random.random(3)*2.0-1 - oris = Orientation(Rotation(set_of_quaternions),lattice)[:200] - for i,pole in enumerate(oris.inverse_pole(axis,SST=SST)): - assert np.allclose(pole,inverse_pole(oris[i],axis,SST=SST)) + def test_to_SST_vectorization(self,lattice,shape,vector,proper): + o = Orientation.from_random(lattice=lattice,shape=shape,seed=0) + for r, theO in zip(o.to_SST(vector=vector,proper=proper).reshape((-1,3)),o.flatten()): + assert np.allclose(r,theO.to_SST(vector=vector,proper=proper)) + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) + @pytest.mark.parametrize('vector',np.array([[1,0,0],[1,2,3],[-1,1,-1]])) + @pytest.mark.parametrize('proper',[True,False]) + def test_IPF_color_vectorization(self,lattice,shape,vector,proper): + o = Orientation.from_random(lattice=lattice,shape=shape,seed=0) + poles = o.to_SST(vector=vector,proper=proper) + for r, theO in zip(o.IPF_color(poles,proper=proper).reshape((-1,3)),o.flatten()): + assert np.allclose(r,theO.IPF_color(theO.to_SST(vector=vector,proper=proper),proper=proper)) + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('a,b',[ + ((2,3,2),(2,3,2)), + ((2,2),(4,4)), + ((3,1),(1,3)), + (None,(3,)), + ]) + def test_to_SST_blending(self,lattice,a,b): + o = Orientation.from_random(lattice=lattice,shape=a,seed=0) + v = np.random.random(b+(3,)) + blend = util.shapeblender(o.shape,b) + for loc in np.random.randint(0,blend,(10,len(blend))): + print(f'{a}/{b} @ {loc}') + print(o[tuple(loc[:len(o.shape)])].to_SST(v[tuple(loc[-len(b):])])) + print(o.to_SST(v)[tuple(loc)]) + assert np.allclose(o[tuple(loc[:len(o.shape)])].to_SST(v[tuple(loc[-len(b):])]), + o.to_SST(v)[tuple(loc)]) @pytest.mark.parametrize('color',[{'label':'red', 'RGB':[1,0,0],'direction':[0,0,1]}, {'label':'green','RGB':[0,1,0],'direction':[0,1,1]}, {'label':'blue', 'RGB':[0,0,1],'direction':[1,1,1]}]) - @pytest.mark.parametrize('lattice',['fcc','bcc']) - def test_IPF_cubic(self,color,lattice): - cube = Orientation(Rotation(),lattice) + @pytest.mark.parametrize('proper',[True,False]) + def test_IPF_cubic(self,color,proper): + cube = Orientation(lattice='cubic') for direction in set(permutations(np.array(color['direction']))): - assert np.allclose(cube.IPF_color(np.array(direction)),np.array(color['RGB'])) + assert np.allclose(np.array(color['RGB']), + cube.IPF_color(cube.to_SST(vector=np.array(direction),proper=proper),proper=proper)) - @pytest.mark.parametrize('lattice',Lattice.lattices) - def test_IPF_equivalent(self,set_of_quaternions,lattice): - direction = np.random.random(3)*2.0-1 - for ori in Orientation(Rotation(set_of_quaternions),lattice)[:200]: - color = ori.IPF_color(direction) - for equivalent in ori.equivalent: - assert np.allclose(color,equivalent.IPF_color(direction)) + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('proper',[True,False]) + def test_IPF_equivalent(self,set_of_quaternions,lattice,proper): + direction = np.random.random(3)*2.0-1.0 + o = Orientation(rotation=set_of_quaternions,lattice=lattice).equivalent + color = o.IPF_color(o.to_SST(vector=direction,proper=proper),proper=proper) + assert np.allclose(np.broadcast_to(color[0,...],color.shape),color) - @pytest.mark.parametrize('lattice',Lattice.lattices) - def test_reduced(self,set_of_quaternions,lattice): - oris = Orientation(Rotation(set_of_quaternions),lattice) - reduced = oris.reduced - assert np.all(reduced.in_FZ) and oris.rotation.shape == reduced.rotation.shape + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + def test_in_FZ_vectorization(self,set_of_rodrigues,lattice): + result = Orientation.from_Rodrigues(rho=set_of_rodrigues.reshape((50,4,-1)),lattice=lattice).in_FZ.reshape(-1) + for r,rho in zip(result,set_of_rodrigues[:len(result)]): + assert r == Orientation.from_Rodrigues(rho=rho,lattice=lattice).in_FZ + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + def test_in_disorientation_FZ_vectorization(self,set_of_rodrigues,lattice): + result = Orientation.from_Rodrigues(rho=set_of_rodrigues.reshape((50,4,-1)), + lattice=lattice).in_disorientation_FZ.reshape(-1) + for r,rho in zip(result,set_of_rodrigues[:len(result)]): + assert r == Orientation.from_Rodrigues(rho=rho,lattice=lattice).in_disorientation_FZ + + @pytest.mark.parametrize('proper',[True,False]) + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + def test_in_SST_vectorization(self,lattice,proper): + vecs = np.random.rand(20,4,3) + result = Orientation(lattice=lattice).in_SST(vecs,proper).flatten() + for r,v in zip(result,vecs.reshape((-1,3))): + assert np.all(r == Orientation(lattice=lattice).in_SST(v,proper)) + + @pytest.mark.parametrize('invalid_lattice',['fcc','bcc','hello']) + def test_invalid_lattice_init(self,invalid_lattice): + with pytest.raises(KeyError): + Orientation(lattice=invalid_lattice) # noqa + + @pytest.mark.parametrize('invalid_family',[None,'fcc','bcc','hello']) + def test_invalid_symmetry_family(self,invalid_family): + with pytest.raises(KeyError): + o = Orientation(lattice='cubic') + o.family = invalid_family + o.symmetry_operations # noqa + + def test_missing_symmetry_equivalent(self): + with pytest.raises(ValueError): + Orientation(lattice=None).equivalent # noqa + + def test_missing_symmetry_reduced(self): + with pytest.raises(ValueError): + Orientation(lattice=None).reduced # noqa + + def test_missing_symmetry_in_FZ(self): + with pytest.raises(ValueError): + Orientation(lattice=None).in_FZ # noqa + + def test_missing_symmetry_in_disorientation_FZ(self): + with pytest.raises(ValueError): + Orientation(lattice=None).in_disorientation_FZ # noqa + + def test_missing_symmetry_disorientation(self): + with pytest.raises(ValueError): + Orientation(lattice=None).disorientation(Orientation(lattice=None)) # noqa + + def test_missing_symmetry_average(self): + with pytest.raises(ValueError): + Orientation(lattice=None).average() # noqa + + def test_missing_symmetry_to_SST(self): + with pytest.raises(ValueError): + Orientation(lattice=None).to_SST(np.zeros(3)) # noqa + + def test_missing_symmetry_immutable(self): + with pytest.raises(KeyError): + Orientation(lattice=None).immutable # noqa + + def test_missing_symmetry_basis_real(self): + with pytest.raises(KeyError): + Orientation(lattice=None).basis_real # noqa + + def test_missing_symmetry_basis_reciprocal(self): + with pytest.raises(KeyError): + Orientation(lattice=None).basis_reciprocal # noqa + + def test_double_Bravais_to_Miller(self): + with pytest.raises(KeyError): + Orientation.Bravais_to_Miller(uvtw=np.ones(4),hkil=np.ones(4)) # noqa + + def test_double_Miller_to_Bravais(self): + with pytest.raises(KeyError): + Orientation.Miller_to_Bravais(uvw=np.ones(4),hkl=np.ones(4)) # noqa + + def test_double_to_lattice(self): + with pytest.raises(KeyError): + Orientation().to_lattice(direction=np.ones(3),plane=np.ones(3)) # noqa + + def test_double_to_frame(self): + with pytest.raises(KeyError): + Orientation().to_frame(uvw=np.ones(3),hkl=np.ones(3)) # noqa + + @pytest.mark.parametrize('relation',[None,'Peter','Paul']) + def test_unknown_relation(self,relation): + with pytest.raises(KeyError): + Orientation(lattice='cF').related(relation) # noqa + + @pytest.mark.parametrize('relation,lattice,a,b,c,alpha,beta,gamma', + [ + ('Bain', 'aP',0.5,2.0,3.0,0.8,0.5,1.2), + ('KS', 'mP',1.0,2.0,3.0,np.pi/2,0.5,np.pi/2), + ('Pitsch', 'oI',0.5,1.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('Burgers','tP',0.5,0.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('GT', 'hP',1.0,None,1.6,np.pi/2,np.pi/2,2*np.pi/3), + ('Burgers','cF',1.0,1.0,None,np.pi/2,np.pi/2,np.pi/2), + ]) + def test_unknown_relation_lattice(self,relation,lattice,a,b,c,alpha,beta,gamma): + with pytest.raises(KeyError): + Orientation(lattice=lattice, + a=a,b=b,c=c, + alpha=alpha,beta=beta,gamma=gamma).related(relation) # noqa + + @pytest.mark.parametrize('lattice',Orientation.crystal_families) + @pytest.mark.parametrize('proper',[True,False]) + def test_in_SST(self,lattice,proper): + assert Orientation(lattice=lattice).in_SST(np.zeros(3),proper) + + @pytest.mark.parametrize('function',['in_SST','IPF_color']) + def test_invalid_argument(self,function): + o = Orientation(lattice='cubic') # noqa + with pytest.raises(ValueError): + eval(f'o.{function}(np.ones(4))') + + @pytest.mark.parametrize('model',lattice.relations) + def test_relationship_definition(self,model): + m,o = list(lattice.relations[model]) + assert lattice.relations[model][m].shape[:-1] == lattice.relations[model][o].shape[:-1] @pytest.mark.parametrize('model',['Bain','KS','GT','GT_prime','NW','Pitsch']) - @pytest.mark.parametrize('lattice',['fcc','bcc']) + @pytest.mark.parametrize('lattice',['cF','cI']) + def test_relationship_vectorize(self,set_of_quaternions,lattice,model): + r = Orientation(rotation=set_of_quaternions[:200].reshape((50,4,4)),lattice=lattice).related(model) + for i in range(200): + assert r.reshape((-1,200))[:,i] == Orientation(set_of_quaternions[i],lattice).related(model) + + @pytest.mark.parametrize('model',['Bain','KS','GT','GT_prime','NW','Pitsch']) + @pytest.mark.parametrize('lattice',['cF','cI']) def test_relationship_forward_backward(self,model,lattice): - ori = Orientation(Rotation.from_random(),lattice) - for i,r in enumerate(ori.related(model)): - ori2 = r.related(model)[i] - misorientation = ori.rotation.misorientation(ori2.rotation) - assert misorientation.as_axis_angle(degrees=True)[3]<1.0e-5 + o = Orientation.from_random(lattice=lattice) + for i,r in enumerate(o.related(model)): + assert o.disorientation(r.related(model)[i]).as_axis_angle(degrees=True,pair=True)[1]<1.0e-5 @pytest.mark.parametrize('model',['Bain','KS','GT','GT_prime','NW','Pitsch']) - @pytest.mark.parametrize('lattice',['fcc','bcc']) + @pytest.mark.parametrize('lattice',['cF','cI']) def test_relationship_reference(self,update,reference_dir,model,lattice): - reference = os.path.join(reference_dir,f'{lattice}_{model}.txt') - ori = Orientation(Rotation(),lattice) - eu = np.array([o.rotation.as_Eulers(degrees=True) for o in ori.related(model)]) + reference = reference_dir/f'{lattice}_{model}.txt' + o = Orientation(lattice=lattice) + eu = o.related(model).as_Eulers(degrees=True) if update: coords = np.array([(1,i+1) for i,x in enumerate(eu)]) - table = Table(eu,{'Eulers':(3,)}) - table = table.add('pos',coords) - table.save(reference) + Table(eu,{'Eulers':(3,)})\ + .add('pos',coords)\ + .save(reference) assert np.allclose(eu,Table.load(reference).get('Eulers')) - @pytest.mark.parametrize('lattice',Lattice.lattices) - def test_disorientation360(self,lattice): - R_1 = Orientation(Rotation(),lattice) - R_2 = Orientation(Rotation.from_Eulers([360,0,0],degrees=True),lattice) - assert np.allclose(R_1.disorientation(R_2).as_matrix(),np.eye(3)) + def test_basis_real(self): + for gamma in np.random.random(2**8)*np.pi: + basis = np.tril(np.random.random((3,3))+1e-6) + basis[1,:2] = basis[1,1]*np.array([np.cos(gamma),np.sin(gamma)]) + basis[2,:2] = basis[2,:2]*2-1 + lengths = np.linalg.norm(basis,axis=-1) + cosines = np.roll(np.einsum('ij,ij->i',basis,np.roll(basis,1,axis=0))/lengths/np.roll(lengths,1),1) + o = Orientation.from_random(lattice='aP', + **dict(zip(['a','b','c'],lengths)), + **dict(zip(['alpha','beta','gamma'],np.arccos(cosines))), + ) + assert np.allclose(o.to_frame(uvw=np.eye(3)),basis), 'Lattice basis disagrees with initialization' - @pytest.mark.parametrize('lattice',Lattice.lattices) - @pytest.mark.parametrize('angle',[10,20,30,40]) - def test_average(self,angle,lattice): - R_1 = Orientation(Rotation.from_axis_angle([0,0,1,10],degrees=True),lattice) - R_2 = Orientation(Rotation.from_axis_angle([0,0,1,angle],degrees=True),lattice) - avg_angle = R_1.average(R_2).rotation.as_axis_angle(degrees=True,pair=True)[1] - assert np.isclose(avg_angle,10+(angle-10)/2.) + @pytest.mark.parametrize('lattice,a,b,c,alpha,beta,gamma', + [ + ('aP',0.5,2.0,3.0,0.8,0.5,1.2), + ('mP',1.0,2.0,3.0,np.pi/2,0.5,np.pi/2), + ('oI',0.5,1.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('tP',0.5,0.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('hP',1.0,None,1.6,np.pi/2,np.pi/2,2*np.pi/3), + ('cF',1.0,1.0,None,np.pi/2,np.pi/2,np.pi/2), + ]) + def test_bases_contraction(self,lattice,a,b,c,alpha,beta,gamma): + L = Orientation(lattice=lattice, + a=a,b=b,c=c, + alpha=alpha,beta=beta,gamma=gamma) + assert np.allclose(np.eye(3),np.einsum('ik,jk',L.basis_real,L.basis_reciprocal)) - @pytest.mark.parametrize('lattice',Lattice.lattices) - def test_from_average(self,lattice): - R_1 = Orientation(Rotation.from_random(),lattice) - eqs = [r for r in R_1.equivalent] - R_2 = Orientation.from_average(eqs) - assert np.allclose(R_1.rotation.quaternion,R_2.rotation.quaternion) + @pytest.mark.parametrize('keyFrame,keyLattice',[('uvw','direction'),('hkl','plane'),]) + @pytest.mark.parametrize('vector',np.array([ + [1.,1.,1.], + [-2.,3.,0.5], + [0.,0.,1.], + [1.,1.,1.], + [2.,2.,2.], + [0.,1.,1.], + ])) + @pytest.mark.parametrize('lattice,a,b,c,alpha,beta,gamma', + [ + ('aP',0.5,2.0,3.0,0.8,0.5,1.2), + ('mP',1.0,2.0,3.0,np.pi/2,0.5,np.pi/2), + ('oI',0.5,1.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('tP',0.5,0.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('hP',1.0,1.0,1.6,np.pi/2,np.pi/2,2*np.pi/3), + ('cF',1.0,1.0,1.0,np.pi/2,np.pi/2,np.pi/2), + ]) + def test_to_frame_to_lattice(self,lattice,a,b,c,alpha,beta,gamma,vector,keyFrame,keyLattice): + L = Orientation(lattice=lattice, + a=a,b=b,c=c, + alpha=alpha,beta=beta,gamma=gamma) + assert np.allclose(vector, + L.to_frame(**{keyFrame:L.to_lattice(**{keyLattice:vector})})) + + @pytest.mark.parametrize('vector',np.array([ + [1,0,0], + [1,1,0], + [1,1,1], + [1,0,-2], + ])) + @pytest.mark.parametrize('kw_Miller,kw_Bravais',[('uvw','uvtw'),('hkl','hkil')]) + def test_Miller_Bravais_Miller(self,vector,kw_Miller,kw_Bravais): + assert np.all(vector == Orientation.Bravais_to_Miller(**{kw_Bravais:Orientation.Miller_to_Bravais(**{kw_Miller:vector})})) + + @pytest.mark.parametrize('vector',np.array([ + [1,0,-1,2], + [1,-1,0,3], + [1,1,-2,-3], + [0,0,0,1], + ])) + @pytest.mark.parametrize('kw_Miller,kw_Bravais',[('uvw','uvtw'),('hkl','hkil')]) + def test_Bravais_Miller_Bravais(self,vector,kw_Miller,kw_Bravais): + assert np.all(vector == Orientation.Miller_to_Bravais(**{kw_Miller:Orientation.Bravais_to_Miller(**{kw_Bravais:vector})})) + + @pytest.mark.parametrize('lattice,a,b,c,alpha,beta,gamma', + [ + ('aP',0.5,2.0,3.0,0.8,0.5,1.2), + ('mP',1.0,2.0,3.0,np.pi/2,0.5,np.pi/2), + ('oI',0.5,1.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('tP',0.5,0.5,3.0,np.pi/2,np.pi/2,np.pi/2), + ('hP',1.0,1.0,1.6,np.pi/2,np.pi/2,2*np.pi/3), + ('cF',1.0,1.0,1.0,np.pi/2,np.pi/2,np.pi/2), + ]) + + @pytest.mark.parametrize('kw',['uvw','hkl']) + @pytest.mark.parametrize('with_symmetry',[False,True]) + @pytest.mark.parametrize('shape',[None,1,(12,24)]) + @pytest.mark.parametrize('vector',[ + np.random.random( 3 ), + np.random.random( (4,3)), + np.random.random((4,8,3)), + ]) + def test_to_pole(self,shape,lattice,a,b,c,alpha,beta,gamma,vector,kw,with_symmetry): + o = Orientation.from_random(shape=shape, + lattice=lattice, + a=a,b=b,c=c, + alpha=alpha,beta=beta,gamma=gamma) + assert o.to_pole(**{kw:vector,'with_symmetry':with_symmetry}).shape \ + == o.shape + (o.symmetry_operations.shape if with_symmetry else ()) + vector.shape + + @pytest.mark.parametrize('lattice',['hP','cI','cF']) + def test_Schmid(self,update,reference_dir,lattice): + L = Orientation(lattice=lattice) + for mode in L.kinematics: + reference = reference_dir/f'{lattice}_{mode}.txt' + P = L.Schmid(mode) + if update: + table = Table(P.reshape(-1,9),{'Schmid':(3,3,)}) + table.save(reference) + assert np.allclose(P,Table.load(reference).get('Schmid')) diff --git a/python/tests/test_Result.py b/python/tests/test_Result.py index fb178e609..929c0ef44 100644 --- a/python/tests/test_Result.py +++ b/python/tests/test_Result.py @@ -168,15 +168,16 @@ class TestResult: @pytest.mark.parametrize('d',[[1,0,0],[0,1,0],[0,0,1]]) def test_add_IPF_color(self,default,d): - default.add_IPF_color('O',d) - loc = {'orientation': default.get_dataset_location('O'), - 'color': default.get_dataset_location('IPFcolor_[{} {} {}]'.format(*d))} - qu = default.read_dataset(loc['orientation']).view(np.double).reshape(-1,4) + default.add_IPF_color('O',np.array(d)) + loc = {'O': default.get_dataset_location('O'), + 'color': default.get_dataset_location('IPFcolor_[{} {} {}]'.format(*d))} + qu = default.read_dataset(loc['O']).view(np.double).squeeze() crystal_structure = default.get_crystal_structure() - in_memory = np.empty((qu.shape[0],3),np.uint8) - for i,q in enumerate(qu): - o = Orientation(q,crystal_structure).reduced - in_memory[i] = np.uint8(o.IPF_color(np.array(d))*255) + c = Orientation(rotation=qu, + lattice={'fcc':'cF', + 'bcc':'cI', + 'hex':'hP'}[crystal_structure]) + in_memory = np.uint8(c.IPF_color(c.to_SST(np.array(d)))*255) in_file = default.read_dataset(loc['color']) assert np.allclose(in_memory,in_file) @@ -244,13 +245,14 @@ class TestResult: in_file = default.read_dataset(loc['S'],0) assert np.allclose(in_memory,in_file) + @pytest.mark.skip(reason='requires rework of lattice.f90') @pytest.mark.parametrize('polar',[True,False]) def test_add_pole(self,default,polar): pole = np.array([1.,0.,0.]) default.add_pole('O',pole,polar) - loc = {'orientation': default.get_dataset_location('O'), - 'pole': default.get_dataset_location('p^{}_[1 0 0)'.format(u'rφ' if polar else 'xy'))} - rot = Rotation(default.read_dataset(loc['orientation']).view(np.double)) + loc = {'O': default.get_dataset_location('O'), + 'pole': default.get_dataset_location('p^{}_[1 0 0)'.format(u'rφ' if polar else 'xy'))} + rot = Rotation(default.read_dataset(loc['O']).view(np.double)) rotated_pole = rot * np.broadcast_to(pole,rot.shape+(3,)) xy = rotated_pole[:,0:2]/(1.+abs(pole[2])) in_memory = xy if not polar else \ diff --git a/python/tests/test_Rotation.py b/python/tests/test_Rotation.py index d70d650a0..9361008e5 100644 --- a/python/tests/test_Rotation.py +++ b/python/tests/test_Rotation.py @@ -771,6 +771,53 @@ class TestRotation: def test_random(self,shape): Rotation.from_random(shape) + def test_equal(self): + r = Rotation.from_random(seed=0) + assert r == r + + def test_unequal(self): + r = Rotation.from_random(seed=0) + assert not (r != r) + + def test_inversion(self): + r = Rotation.from_random(seed=0) + assert r == ~~r + + @pytest.mark.parametrize('shape',[None,1,(1,),(4,2),(1,1,1)]) + def test_shape(self,shape): + r = Rotation.from_random(shape=shape) + assert r.shape == (shape if isinstance(shape,tuple) else (shape,) if shape else ()) + + @pytest.mark.parametrize('shape',[None,1,(1,),(4,2),(3,3,2)]) + def test_append(self,shape): + r = Rotation.from_random(shape=shape) + p = Rotation.from_random(shape=shape) + s = r.append(p) + print(f'append 2x {shape} --> {s.shape}') + assert s[0,...] == r[0,...] and s[-1,...] == p[-1,...] + + @pytest.mark.parametrize('quat,standardized',[ + ([-1,0,0,0],[1,0,0,0]), + ([-0.5,-0.5,-0.5,-0.5],[0.5,0.5,0.5,0.5]), + ]) + def test_standardization(self,quat,standardized): + assert Rotation(quat)._standardize() == Rotation(standardized) + + @pytest.mark.parametrize('shape,length',[ + ((2,3,4),2), + (4,4), + ((),0) + ]) + def test_len(self,shape,length): + r = Rotation.from_random(shape=shape) + assert len(r) == length + + @pytest.mark.parametrize('shape',[(4,6),(2,3,4),(3,3,3)]) + @pytest.mark.parametrize('order',['C','F']) + def test_flatten_reshape(self,shape,order): + r = Rotation.from_random(shape=shape) + assert r == r.flatten(order).reshape(shape,order) + @pytest.mark.parametrize('function',[Rotation.from_quaternion, Rotation.from_Eulers, Rotation.from_axis_angle, @@ -848,7 +895,8 @@ class TestRotation: np.random.rand(3,3,3,3)]) def test_rotate_identity(self,data): R = Rotation() - assert np.allclose(data,R*data) + print(R,data) + assert np.allclose(data,R@data) @pytest.mark.parametrize('data',[np.random.rand(3), np.random.rand(3,3), @@ -860,6 +908,16 @@ class TestRotation: R_2 = Rotation.from_Eulers(np.array([0.,0.,phi_2])) assert np.allclose(data,R_2@(R_1@data)) + @pytest.mark.parametrize('pwr',[-10,0,1,2.5,np.pi,np.random.random()]) + def test_rotate_power(self,pwr): + R = Rotation.from_random() + axis_angle = R.as_axis_angle() + axis_angle[ 3] = (pwr*axis_angle[-1])%(2.*np.pi) + if axis_angle[3] > np.pi: + axis_angle[3] -= 2.*np.pi + axis_angle *= -1 + assert R**pwr == Rotation.from_axis_angle(axis_angle) + def test_rotate_inverse(self): R = Rotation.from_random() assert np.allclose(np.eye(3),(~R@R).as_matrix()) @@ -877,7 +935,7 @@ class TestRotation: def test_rotate_invalid_shape(self,data): R = Rotation.from_random() with pytest.raises(ValueError): - R*data + R@data @pytest.mark.parametrize('data',['does_not_work', (1,2), @@ -885,7 +943,7 @@ class TestRotation: def test_rotate_invalid_type(self,data): R = Rotation.from_random() with pytest.raises(TypeError): - R*data + R@data def test_misorientation(self): R = Rotation.from_random() @@ -898,9 +956,8 @@ class TestRotation: @pytest.mark.parametrize('angle',[10,20,30,40,50,60,70,80,90,100,120]) def test_average(self,angle): - R_1 = Rotation.from_axis_angle([0,0,1,10],degrees=True) - R_2 = Rotation.from_axis_angle([0,0,1,angle],degrees=True) - avg_angle = R_1.average(R_2).as_axis_angle(degrees=True,pair=True)[1] + R = Rotation.from_axis_angle([[0,0,1,10],[0,0,1,angle]],degrees=True) + avg_angle = R.average().as_axis_angle(degrees=True,pair=True)[1] assert np.isclose(avg_angle,10+(angle-10)/2.) diff --git a/python/tests/test_util.py b/python/tests/test_util.py index 0b5593e8e..19b6377f6 100644 --- a/python/tests/test_util.py +++ b/python/tests/test_util.py @@ -44,3 +44,52 @@ class TestUtil: selected = util.hybrid_IA(dist,N_samples) dist_sampled = np.histogram(centers[selected],bins)[0]/N_samples*np.sum(dist) assert np.sqrt(((dist - dist_sampled) ** 2).mean()) < .025 and selected.shape[0]==N_samples + + @pytest.mark.parametrize('point,normalize,answer', + [ + ([1,0,0],False,[1,0,0]), + ([1,0,0],True, [1,0,0]), + ([0,1,1],False,[0,0.5,0]), + ([0,1,1],True, [0,0.41421356,0]), + ([1,1,1],False,[0.5,0.5,0]), + ([1,1,1],True, [0.3660254, 0.3660254, 0]), + ]) + def test_project_stereographic(self,point,normalize,answer): + assert np.allclose(util.project_stereographic(np.array(point),normalize=normalize),answer) + + @pytest.mark.parametrize('fro,to,mode,answer', + [ + ((),(1,),'left',(1,)), + ((1,),(7,),'right',(1,)), + ((1,2),(1,1,2,2),'right',(1,1,2,1)), + ((1,2),(1,1,2,2),'left',(1,1,1,2)), + ((1,2,3),(1,1,2,3,4),'right',(1,1,2,3,1)), + ((10,2),(10,3,2,2,),'right',(10,1,2,1)), + ((10,2),(10,3,2,2,),'left',(10,1,1,2)), + ((2,2,3),(2,2,2,3,4),'left',(1,2,2,3,1)), + ((2,2,3),(2,2,2,3,4),'right',(2,2,1,3,1)), + ]) + def test_shapeshifter(self,fro,to,mode,answer): + assert util.shapeshifter(fro,to,mode) == answer + + @pytest.mark.parametrize('fro,to,mode', + [ + ((10,3,4),(10,3,2,2),'left'), + ((2,3),(10,3,2,2),'right'), + ]) + def test_invalid_shapeshifter(self,fro,to,mode): + with pytest.raises(ValueError): + util.shapeshifter(fro,to,mode) + + @pytest.mark.parametrize('a,b,answer', + [ + ((),(1,),(1,)), + ((1,),(),(1,)), + ((1,),(7,),(1,7)), + ((2,),(2,2),(2,2)), + ((1,2),(2,2),(1,2,2)), + ((1,2,3),(2,3,4),(1,2,3,4)), + ((1,2,3),(1,2,3),(1,2,3)), + ]) + def test_shapeblender(self,a,b,answer): + assert util.shapeblender(a,b) == answer diff --git a/src/prec.f90 b/src/prec.f90 index b2866a4f4..738775e3b 100644 --- a/src/prec.f90 +++ b/src/prec.f90 @@ -93,7 +93,8 @@ subroutine prec_init print'(a,i19)', ' Maximum value: ',huge(0) print'(/,a,i3)', ' Size of float in bit: ',storage_size(0.0_pReal) print'(a,e10.3)', ' Maximum value: ',huge(0.0_pReal) - print'(a,e10.3)', ' Minimum value: ',tiny(0.0_pReal) + print'(a,e10.3)', ' Minimum value: ',PREAL_MIN + print'(a,e10.3)', ' Epsilon value: ',PREAL_EPSILON print'(a,i3)', ' Decimal precision: ',precision(0.0_pReal) call selfTest