From 56afc03f3cb108c63fc394e7ab2413864fac6a12 Mon Sep 17 00:00:00 2001 From: Martin Diehl Date: Wed, 27 May 2020 18:05:08 +0200 Subject: [PATCH] only vectorized version needed use single point/simple versions only for testing --- python/damask/mechanics.py | 20 ++++++-------------- python/tests/test_mechanics.py | 15 +++++++++++++++ 2 files changed, 21 insertions(+), 14 deletions(-) diff --git a/python/damask/mechanics.py b/python/damask/mechanics.py index e388b3569..484341772 100644 --- a/python/damask/mechanics.py +++ b/python/damask/mechanics.py @@ -28,12 +28,12 @@ def deviatoric_part(T): Parameters ---------- - T : numpy.ndarray of shape (:,3,3) or (3,3) + T : numpy.ndarray of shape (...,3,3) Tensor of which the deviatoric part is computed. """ - return T - _np.eye(3)*spherical_part(T) if _np.shape(T) == (3,3) else \ - T - _np.einsum('...ij,...->...ij',_np.eye(3),spherical_part(T)) + return T - _np.einsum('...ij,...->...ij',_np.eye(3),spherical_part(T)) + def eigenvalues(T_sym): """ @@ -180,22 +180,14 @@ def spherical_part(T,tensor=False): Parameters ---------- - T : numpy.ndarray of shape (:,3,3) or (3,3) + T : numpy.ndarray of shape (...,3,3) Tensor of which the hydrostatic part is computed. tensor : bool, optional Map spherical part onto identity tensor. Default is false """ - if T.shape == (3,3): - sph = _np.trace(T)/3.0 - return sph if not tensor else _np.eye(3)*sph - else: - sph = _np.trace(T,axis2=-2,axis1=-1)/3.0 - if not tensor: - return sph - else: - #return _np.einsum('ijk,i->ijk',_np.broadcast_to(_np.eye(3),(T.shape[0],3,3)),sph) - return _np.einsum('...jk,...->...jk',_np.eye(3),sph) + sph = _np.trace(T,axis2=-2,axis1=-1)/3.0 + return _np.einsum('...jk,...->...jk',_np.eye(3),sph) if tensor else sph def strain_tensor(F,t,m): diff --git a/python/tests/test_mechanics.py b/python/tests/test_mechanics.py index f7128401f..da9e83f7f 100644 --- a/python/tests/test_mechanics.py +++ b/python/tests/test_mechanics.py @@ -3,12 +3,27 @@ import numpy as np from damask import mechanics +def deviatoric_part(T): + return T - np.eye(3)*spherical_part(T) + +def spherical_part(T,tensor=False): + sph = np.trace(T)/3.0 + return sph if not tensor else np.eye(3)*sph + class TestMechanics: n = 1000 c = np.random.randint(n) + @pytest.mark.parametrize('vectorized,single',[(mechanics.deviatoric_part, deviatoric_part), + (mechanics.spherical_part, spherical_part) + ]) + def test_vectorize_1_arg_(self,vectorized,single): + test_data = np.random.rand(self.n,3,3) + for i,v in enumerate(vectorized(test_data)): + assert np.allclose(single(test_data[i]),v) + @pytest.mark.parametrize('function',[mechanics.deviatoric_part, mechanics.eigenvalues, mechanics.eigenvectors,