223031012
/
DAMASK_EICMD
mirror of https://github.com/achalhp/DAMASK_EICMD.git
1
0
Fork 0
Code Issues Projects Releases Wiki Activity

Merge branch 'account-for-floating-point-precision-in-orientation' into 'development' 

ensures that at least one orientation in the FZ is found

See merge request damask/DAMASK!337
This commit is contained in:
Philip Eisenlohr 2021-02-15 20:53:37 +00:00
parent 8933dd8562 4bcbcb34d0
commit 002cfbd085
2 changed files with 23 additions and 12 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

24
python/damask/_orientation.py
View File

@ -524,26 +524,26 @@ class Orientation(Rotation):
if self.family is None: if self.family is None:
raise ValueError('Missing crystal symmetry') raise ValueError('Missing crystal symmetry')
rho_abs = np.abs(self.as_Rodrigues_vector(compact=True)) rho_abs = np.abs(self.as_Rodrigues_vector(compact=True))*(1.-1.e-9)
with np.errstate(invalid='ignore'): with np.errstate(invalid='ignore'):
# using '*'/prod for 'and' # using '*'/prod for 'and'
if self.family == 'cubic': if self.family == 'cubic':
return (np.prod(np.sqrt(2)-1. >= rho_abs,axis=-1) * return (np.prod(np.sqrt(2)-1. >= rho_abs,axis=-1) *
(1. >= np.sum(rho_abs,axis=-1))).astype(np.bool) (1. >= np.sum(rho_abs,axis=-1))).astype(bool)
elif self.family == 'hexagonal': elif self.family == 'hexagonal':
return (np.prod(1. >= rho_abs,axis=-1) * return (np.prod(1. >= rho_abs,axis=-1) *
(2. >= np.sqrt(3)*rho_abs[...,0] + rho_abs[...,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[...,1] + rho_abs[...,0]) *
(2. >= np.sqrt(3) + rho_abs[...,2])).astype(np.bool) (2. >= np.sqrt(3) + rho_abs[...,2])).astype(bool)
elif self.family == 'tetragonal': elif self.family == 'tetragonal':
return (np.prod(1. >= rho_abs[...,:2],axis=-1) * return (np.prod(1. >= rho_abs[...,:2],axis=-1) *
(np.sqrt(2) >= rho_abs[...,0] + rho_abs[...,1]) * (np.sqrt(2) >= rho_abs[...,0] + rho_abs[...,1]) *
(np.sqrt(2) >= rho_abs[...,2] + 1.)).astype(np.bool) (np.sqrt(2) >= rho_abs[...,2] + 1.)).astype(bool)
elif self.family == 'orthorhombic': elif self.family == 'orthorhombic':
return (np.prod(1. >= rho_abs,axis=-1)).astype(np.bool) return (np.prod(1. >= rho_abs,axis=-1)).astype(bool)
elif self.family == 'monoclinic': elif self.family == 'monoclinic':
return (1. >= rho_abs[...,1]).astype(np.bool) return (1. >= rho_abs[...,1]).astype(bool)
else: else:
return np.all(np.isfinite(rho_abs),axis=-1) return np.all(np.isfinite(rho_abs),axis=-1)
@ -567,28 +567,28 @@ class Orientation(Rotation):
if self.family is None: if self.family is None:
raise ValueError('Missing crystal symmetry') raise ValueError('Missing crystal symmetry')
rho = self.as_Rodrigues_vector(compact=True) rho = self.as_Rodrigues_vector(compact=True)*(1.0-1.0e-9)
with np.errstate(invalid='ignore'): with np.errstate(invalid='ignore'):
if self.family == 'cubic': if self.family == 'cubic':
return ((rho[...,0] >= rho[...,1]) & return ((rho[...,0] >= rho[...,1]) &
(rho[...,1] >= rho[...,2]) & (rho[...,1] >= rho[...,2]) &
(rho[...,2] >= 0)).astype(np.bool) (rho[...,2] >= 0)).astype(bool)
elif self.family == 'hexagonal': elif self.family == 'hexagonal':
return ((rho[...,0] >= rho[...,1]*np.sqrt(3)) & return ((rho[...,0] >= rho[...,1]*np.sqrt(3)) &
(rho[...,1] >= 0) & (rho[...,1] >= 0) &
(rho[...,2] >= 0)).astype(np.bool) (rho[...,2] >= 0)).astype(bool)
elif self.family == 'tetragonal': elif self.family == 'tetragonal':
return ((rho[...,0] >= rho[...,1]) & return ((rho[...,0] >= rho[...,1]) &
(rho[...,1] >= 0) & (rho[...,1] >= 0) &
(rho[...,2] >= 0)).astype(np.bool) (rho[...,2] >= 0)).astype(bool)
elif self.family == 'orthorhombic': elif self.family == 'orthorhombic':
return ((rho[...,0] >= 0) & return ((rho[...,0] >= 0) &
(rho[...,1] >= 0) & (rho[...,1] >= 0) &
(rho[...,2] >= 0)).astype(np.bool) (rho[...,2] >= 0)).astype(bool)
elif self.family == 'monoclinic': elif self.family == 'monoclinic':
return ((rho[...,1] >= 0) & return ((rho[...,1] >= 0) &
(rho[...,2] >= 0)).astype(np.bool) (rho[...,2] >= 0)).astype(bool)
else: else:
return np.ones_like(rho[...,0],dtype=bool) return np.ones_like(rho[...,0],dtype=bool)

11
python/tests/test_Orientation.py
View File

@ -7,6 +7,7 @@ from damask import Orientation
from damask import Table from damask import Table
from damask import lattice from damask import lattice
from damask import util from damask import util
from damask import grid_filters
@pytest.fixture @pytest.fixture
@ -236,6 +237,16 @@ class TestOrientation:
for r, theO in zip(o.reduced.flatten(),o.flatten()): for r, theO in zip(o.reduced.flatten(),o.flatten()):
assert r == theO.reduced assert r == theO.reduced
@pytest.mark.parametrize('lattice',Orientation.crystal_families)
def test_reduced_corner_cases(self,lattice):
# test whether there is always a sym-eq rotation that falls into the FZ
N = np.random.randint(10,40)
size = np.ones(3)*np.pi**(2./3.)
grid = grid_filters.coordinates0_node([N+1,N+1,N+1],size,-size*.5)
evenly_distributed = Orientation.from_cubochoric(c=grid[:-2,:-2,:-2],lattice=lattice)
assert evenly_distributed.shape == evenly_distributed.reduced.shape
@pytest.mark.parametrize('lattice',Orientation.crystal_families) @pytest.mark.parametrize('lattice',Orientation.crystal_families)
@pytest.mark.parametrize('shape',[(1),(2,3),(4,3,2)]) @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('vector',np.array([[1,0,0],[1,2,3],[-1,1,-1]]))