return rotation type (ensures proper rotation)

This commit is contained in:
Martin Diehl 2020-11-19 22:36:19 +01:00
parent a4b5c2a537
commit a87596cefc
7 changed files with 46 additions and 20 deletions

View File

@ -956,9 +956,9 @@ class Orientation(Rotation):
"""
if self.family is None or other.family is None:
raise ValueError('Missing crystal symmetry')
raise ValueError('missing crystal symmetry')
if self.family != other.family:
raise NotImplementedError('Disorientation between different crystal families not supported yet.')
raise NotImplementedError('disorientation between different crystal families')
blend = util.shapeblender(self.shape,other.shape)
s = self.equivalent

View File

@ -944,7 +944,7 @@ class Result:
@staticmethod
def _add_rotational(F):
return {
'data': mechanics.rotational(F['data']),
'data': mechanics.rotational(F['data']).as_matrix(),
'label': f"R({F['label']})",
'meta': {
'Unit': F['meta']['Unit'],

View File

@ -61,18 +61,21 @@ class Rotation:
elif np.array(rotation).shape[-1] == 4:
self.quaternion = np.array(rotation)
else:
raise ValueError('"rotation" is neither a Rotation nor a quaternion')
raise TypeError('"rotation" is neither a Rotation nor a quaternion')
def __repr__(self):
"""Represent rotation as unit quaternion, rotation matrix, and Bunge-Euler angles."""
return 'As 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_Euler_angles(degrees=True)),
])
if self == Rotation():
return 'Rotation()'
else:
return f'Quaternions {self.shape}:\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_Euler_angles(degrees=True)),
])
# ToDo: Check difference __copy__ vs __deepcopy__

View File

@ -1,6 +1,7 @@
"""Finite-strain continuum mechanics."""
from . import tensor as _tensor
from . import _rotation
import numpy as _np
@ -107,11 +108,11 @@ def rotational(T):
Returns
-------
R : numpy.ndarray of shape (...,3,3)
Rotational part.
R : damask.Rotation of shape (...)
Rotational part of the vector.
"""
return _polar_decomposition(T,'R')[0]
return _rotation.Rotation.from_matrix(_polar_decomposition(T,'R')[0])
def strain(F,t,m):
@ -260,7 +261,7 @@ def _polar_decomposition(T,requested):
output.append(_np.einsum('...ji,...jk',R,T))
if len(output) == 0:
raise ValueError('Output needs to be out of V,R,U')
raise ValueError('output needs to be out of V, R, U')
return tuple(output)

View File

@ -265,7 +265,7 @@ class TestResult:
default.add_rotational('F')
loc = {'F': default.get_dataset_location('F'),
'R(F)': default.get_dataset_location('R(F)')}
in_memory = mechanics.rotational(default.read_dataset(loc['F'],0))
in_memory = mechanics.rotational(default.read_dataset(loc['F'],0)).as_matrix()
in_file = default.read_dataset(loc['R(F)'],0)
assert np.allclose(in_memory,in_file)

View File

@ -682,6 +682,9 @@ class TestRotation:
for v,o in zip(vec,ori):
assert np.allclose(func(v,normalize=True).as_quaternion(),o.as_quaternion())
def test_invalid_init(self):
with pytest.raises(TypeError):
Rotation(np.ones(3))
@pytest.mark.parametrize('degrees',[True,False])
def test_Eulers(self,set_of_rotations,degrees):
@ -831,6 +834,13 @@ class TestRotation:
with pytest.raises(ValueError):
function(invalid_shape)
def test_invalid_shape_parallel(self):
invalid_a = np.random.random(np.random.randint(8,32,(3)))
invalid_b = np.random.random(np.random.randint(8,32,(3)))
with pytest.raises(ValueError):
Rotation.from_parallel(invalid_a,invalid_b)
@pytest.mark.parametrize('fr,to',[(Rotation.from_quaternion,'as_quaternion'),
(Rotation.from_axis_angle,'as_axis_angle'),
(Rotation.from_Rodrigues_vector, 'as_Rodrigues_vector'),

View File

@ -3,6 +3,7 @@ import numpy as np
from damask import tensor
from damask import mechanics
from damask import Rotation
def stress_Cauchy(P,F):
sigma = 1.0/np.linalg.det(F) * np.dot(P,F.T)
@ -96,7 +97,6 @@ class TestMechanics:
@pytest.mark.parametrize('vectorized,single',[(mechanics.maximum_shear, maximum_shear),
(mechanics.equivalent_stress_Mises, equivalent_stress_Mises),
(mechanics.equivalent_strain_Mises, equivalent_strain_Mises),
(mechanics.rotational, rotational),
(mechanics.stretch_left, stretch_left),
(mechanics.stretch_right, stretch_right),
])
@ -106,6 +106,14 @@ class TestMechanics:
for i,v in enumerate(np.reshape(vectorized(epsilon_vec),vectorized(epsilon).shape)):
assert np.allclose(single(epsilon[i]),v)
def test_vectorize_rotation(self):
epsilon = Rotation.from_random(self.n).as_matrix()
epsilon_vec = np.reshape(epsilon,(self.n//10,10,3,3))
for i,v in enumerate(np.reshape(mechanics.rotational(epsilon_vec).as_matrix(),
mechanics.rotational(epsilon).as_matrix().shape)):
assert np.allclose(rotational(epsilon[i]),v)
@pytest.mark.parametrize('vectorized,single',[(mechanics.stress_Cauchy, stress_Cauchy),
(mechanics.stress_second_Piola_Kirchhoff, stress_second_Piola_Kirchhoff)
])
@ -138,7 +146,7 @@ class TestMechanics:
def test_polar_decomposition(self):
"""F = RU = VR."""
F = np.broadcast_to(np.eye(3),[self.n,3,3])*np.random.rand(self.n,3,3)
R = mechanics.rotational(F)
R = mechanics.rotational(F).as_matrix()
V = mechanics.stretch_left(F)
U = mechanics.stretch_right(F)
assert np.allclose(np.matmul(R,U),
@ -162,7 +170,7 @@ class TestMechanics:
@pytest.mark.parametrize('t',['V','U'])
def test_strain_rotation(self,m,t):
"""Ensure that pure rotation results in no strain."""
F = mechanics.rotational(np.random.rand(self.n,3,3))
F = Rotation.from_random(self.n).as_matrix()
assert np.allclose(mechanics.strain(F,t,m),
0.0)
@ -173,7 +181,7 @@ class TestMechanics:
Should be +1, but random F might contain a reflection.
"""
x = np.random.rand(self.n,3,3)
assert np.allclose(np.abs(np.linalg.det(mechanics.rotational(x))),
assert np.allclose(np.abs(np.linalg.det(mechanics._polar_decomposition(x,'R')[0])),
1.0)
def test_deviatoric_Mises(self):
@ -203,3 +211,7 @@ class TestMechanics:
"""Ensure that sherical stress has max shear of 0.0."""
A = tensor.spherical(tensor.symmetric(np.random.rand(self.n,3,3)),True)
assert np.allclose(mechanics.maximum_shear(A),0.0)
def test_invalid_decomposition(self):
with pytest.raises(ValueError):
mechanics._polar_decomposition(np.random.rand(10,3,3),'A')