clearer names
This commit is contained in:
Martin Diehl
parent
4970f22e13
commit
15c3cab549
|
|
@ -21,21 +21,21 @@ def Cauchy(P,F):
|
|||
return symmetric(sigma)
|
||||
|
||||
|
||||
def deviatoric_part(x):
|
||||
def deviatoric_part(T):
|
||||
"""
|
||||
Return deviatoric part of a tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the deviatoric part is computed.
|
||||
|
||||
"""
|
||||
return x - np.eye(3)*spherical_part(x) if np.shape(x) == (3,3) else \
|
||||
x - np.einsum('ijk,i->ijk',np.broadcast_to(np.eye(3),[x.shape[0],3,3]),spherical_part(x))
|
||||
return T - np.eye(3)*spherical_part(T) if np.shape(T) == (3,3) else \
|
||||
T - np.einsum('ijk,i->ijk',np.broadcast_to(np.eye(3),[T.shape[0],3,3]),spherical_part(T))
|
||||
|
||||
|
||||
def eigenvalues(x):
|
||||
def eigenvalues(T_sym):
|
||||
"""
|
||||
Return the eigenvalues, i.e. principal components, of a symmetric tensor.
|
||||
|
||||
|
|
@ -44,14 +44,14 @@ def eigenvalues(x):
|
|||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T_sym : numpy.array of shape (:,3,3) or (3,3)
|
||||
Symmetric tensor of which the eigenvalues are computed.
|
||||
|
||||
"""
|
||||
return np.linalg.eigvalsh(symmetric(x))
|
||||
return np.linalg.eigvalsh(symmetric(T_sym))
|
||||
|
||||
|
||||
def eigenvectors(x,RHS=False):
|
||||
def eigenvectors(T_sym,RHS=False):
|
||||
"""
|
||||
Return eigenvectors of a symmetric tensor.
|
||||
|
||||
|
|
@ -59,47 +59,47 @@ def eigenvectors(x,RHS=False):
|
|||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T_sym : numpy.array of shape (:,3,3) or (3,3)
|
||||
Symmetric tensor of which the eigenvectors are computed.
|
||||
RHS: bool, optional
|
||||
Enforce right-handed coordinate system. Default is False.
|
||||
|
||||
"""
|
||||
(u,v) = np.linalg.eigh(symmetric(x))
|
||||
(u,v) = np.linalg.eigh(symmetric(T_sym))
|
||||
|
||||
if RHS:
|
||||
if np.shape(x) == (3,3):
|
||||
if np.shape(T_sym) == (3,3):
|
||||
if np.linalg.det(v) < 0.0: v[:,2] *= -1.0
|
||||
else:
|
||||
v[np.linalg.det(v) < 0.0,:,2] *= -1.0
|
||||
return v
|
||||
|
||||
|
||||
def left_stretch(x):
|
||||
def left_stretch(T):
|
||||
"""
|
||||
Return the left stretch of a tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the left stretch is computed.
|
||||
|
||||
"""
|
||||
return __polar_decomposition(x,'V')[0]
|
||||
return __polar_decomposition(T,'V')[0]
|
||||
|
||||
|
||||
def maximum_shear(x):
|
||||
def maximum_shear(T_sym):
|
||||
"""
|
||||
Return the maximum shear component of a symmetric tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T_sym : numpy.array of shape (:,3,3) or (3,3)
|
||||
Symmetric tensor of which the maximum shear is computed.
|
||||
|
||||
"""
|
||||
w = eigenvalues(x)
|
||||
return (w[0] - w[2])*0.5 if np.shape(x) == (3,3) else \
|
||||
w = eigenvalues(T_sym)
|
||||
return (w[0] - w[2])*0.5 if np.shape(T_sym) == (3,3) else \
|
||||
(w[:,0] - w[:,2])*0.5
|
||||
|
||||
|
||||
|
|
@ -147,53 +147,54 @@ def PK2(P,F):
|
|||
S = np.einsum('ijk,ikl->ijl',np.linalg.inv(F),P)
|
||||
return symmetric(S)
|
||||
|
||||
def right_stretch(x):
|
||||
|
||||
def right_stretch(T):
|
||||
"""
|
||||
Return the right stretch of a tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the right stretch is computed.
|
||||
|
||||
"""
|
||||
return __polar_decomposition(x,'U')[0]
|
||||
return __polar_decomposition(T,'U')[0]
|
||||
|
||||
|
||||
def rotational_part(x):
|
||||
def rotational_part(T):
|
||||
"""
|
||||
Return the rotational part of a tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the rotational part is computed.
|
||||
|
||||
"""
|
||||
return __polar_decomposition(x,'R')[0]
|
||||
return __polar_decomposition(T,'R')[0]
|
||||
|
||||
|
||||
def spherical_part(x,tensor=False):
|
||||
def spherical_part(T,tensor=False):
|
||||
"""
|
||||
Return spherical (hydrostatic) part of a tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the hydrostatic part is computed.
|
||||
tensor : bool, optional
|
||||
Map spherical part onto identity tensor. Default is false
|
||||
|
||||
"""
|
||||
if x.shape == (3,3):
|
||||
sph = np.trace(x)/3.0
|
||||
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(x,axis1=1,axis2=2)/3.0
|
||||
sph = np.trace(T,axis1=1,axis2=2)/3.0
|
||||
if not tensor:
|
||||
return sph
|
||||
else:
|
||||
return np.einsum('ijk,i->ijk',np.broadcast_to(np.eye(3),(x.shape[0],3,3)),sph)
|
||||
return np.einsum('ijk,i->ijk',np.broadcast_to(np.eye(3),(T.shape[0],3,3)),sph)
|
||||
|
||||
|
||||
def strain_tensor(F,t,m):
|
||||
|
|
@ -234,73 +235,73 @@ def strain_tensor(F,t,m):
|
|||
eps
|
||||
|
||||
|
||||
def symmetric(x):
|
||||
def symmetric(T):
|
||||
"""
|
||||
Return the symmetrized tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the symmetrized values are computed.
|
||||
|
||||
"""
|
||||
return (x+transpose(x))*0.5
|
||||
return (T+transpose(T))*0.5
|
||||
|
||||
|
||||
def transpose(x):
|
||||
def transpose(T):
|
||||
"""
|
||||
Return the transpose of a tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the transpose is computed.
|
||||
|
||||
"""
|
||||
return x.T if np.shape(x) == (3,3) else \
|
||||
np.transpose(x,(0,2,1))
|
||||
return T.T if np.shape(T) == (3,3) else \
|
||||
np.transpose(T,(0,2,1))
|
||||
|
||||
|
||||
def __polar_decomposition(x,requested):
|
||||
def __polar_decomposition(T,requested):
|
||||
"""
|
||||
Singular value decomposition.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T : numpy.array of shape (:,3,3) or (3,3)
|
||||
Tensor of which the singular values are computed.
|
||||
requested : iterable of str
|
||||
Requested outputs: ‘R’ for the rotation tensor,
|
||||
‘V’ for left stretch tensor and ‘U’ for right stretch tensor.
|
||||
|
||||
"""
|
||||
u, s, vh = np.linalg.svd(x)
|
||||
R = np.dot(u,vh) if np.shape(x) == (3,3) else \
|
||||
u, s, vh = np.linalg.svd(T)
|
||||
R = np.dot(u,vh) if np.shape(T) == (3,3) else \
|
||||
np.einsum('ijk,ikl->ijl',u,vh)
|
||||
|
||||
output = []
|
||||
if 'R' in requested:
|
||||
output.append(R)
|
||||
if 'V' in requested:
|
||||
output.append(np.dot(x,R.T) if np.shape(x) == (3,3) else np.einsum('ijk,ilk->ijl',x,R))
|
||||
output.append(np.dot(T,R.T) if np.shape(T) == (3,3) else np.einsum('ijk,ilk->ijl',T,R))
|
||||
if 'U' in requested:
|
||||
output.append(np.dot(R.T,x) if np.shape(x) == (3,3) else np.einsum('ikj,ikl->ijl',R,x))
|
||||
output.append(np.dot(R.T,T) if np.shape(T) == (3,3) else np.einsum('ikj,ikl->ijl',R,T))
|
||||
|
||||
return tuple(output)
|
||||
|
||||
|
||||
def __Mises(x,s):
|
||||
def __Mises(T_sym,s):
|
||||
"""
|
||||
Base equation for Mises equivalent of a stres or strain tensor.
|
||||
|
||||
Parameters
|
||||
----------
|
||||
x : numpy.array of shape (:,3,3) or (3,3)
|
||||
T_sym : numpy.array of shape (:,3,3) or (3,3)
|
||||
Symmetric tensor of which the von Mises equivalent is computed.
|
||||
s : float
|
||||
Scaling factor (2/3 for strain, 3/2 for stress).
|
||||
|
||||
|
||||
"""
|
||||
d = deviatoric_part(x)
|
||||
return np.sqrt(s*(np.sum(d**2.0))) if np.shape(x) == (3,3) else \
|
||||
d = deviatoric_part(T_sym)
|
||||
return np.sqrt(s*(np.sum(d**2.0))) if np.shape(T_sym) == (3,3) else \
|
||||
np.sqrt(s*np.einsum('ijk->i',d**2.0))
|
||||
|
|
|
|||
Loading…
Reference in New Issue