better have tensor functions available

allows in-memory evaluation of results
This commit is contained in:
Martin Diehl 2019-10-18 20:50:03 +02:00
parent 1a34a6f7b5
commit 3336cfc3da
3 changed files with 32 additions and 7 deletions

View File

@ -6,6 +6,7 @@ with open(os.path.join(os.path.dirname(__file__),'VERSION')) as f:
name = 'damask'
# classes
from .environment import Environment # noqa
from .asciitable import ASCIItable # noqa
@ -14,8 +15,11 @@ from .colormaps import Colormap, Color # noqa
from .orientation import Symmetry, Lattice, Rotation, Orientation # noqa
from .dadf5 import DADF5 # noqa
#from .block import Block # only one class
from .geom import Geom # noqa
from .solver import Solver # noqa
from .test import Test # noqa
from .util import extendableOption # noqa
# functions in modules
from . import mechanics

View File

@ -7,6 +7,7 @@ import numpy as np
from . import util
from . import version
from . import mechanics
# ------------------------------------------------------------------
class DADF5():
@ -379,10 +380,8 @@ class DADF5():
Resulting tensor is symmetrized as the Cauchy stress should be symmetric.
"""
def Cauchy(F,P):
sigma = np.einsum('i,ijk,ilk->ijl',1.0/np.linalg.det(F['data']),P['data'],F['data'])
sigma = (sigma + np.transpose(sigma,(0,2,1)))*0.5 # enforce symmetry
return {
'data' : sigma,
'data' : mechanics.Cauchy(F['data']),P['data']),
'label' : 'sigma',
'meta' : {
'Unit' : P['meta']['Unit'],
@ -529,13 +528,12 @@ class DADF5():
def add_deviator(self,x):
"""Adds the deviator of a tensor."""
def deviator(x):
d = x['data']
if not np.all(np.array(d.shape[1:]) == np.array([3,3])):
if not np.all(np.array(x['data'].shape[1:]) == np.array([3,3])):
raise ValueError
return {
'data' : d - np.einsum('ijk,i->ijk',np.broadcast_to(np.eye(3),[d.shape[0],3,3]),np.trace(d,axis1=1,axis2=2)/3.0),
'data' : mechanics.deviator(x['data']),
'label' : 'dev({})'.format(x['label']),
'meta' : {
'Unit' : x['meta']['Unit'],

View File

@ -0,0 +1,23 @@
import numpy as np
def Cauchy(F,P):
if np.shape(F) == np.shape(P) == (3,3):
sigma = 1.0/np.linalg.det(F) * np.dot(F,P)
return (sigma+sigma.T)*0.5
else:
sigma = np.einsum('i,ijk,ilk->ijl',1.0/np.linalg.det(F),P,F)
return (sigma + np.transpose(sigma,(0,2,1)))*0.5
def deviator(x):
if np.shape(x) == (3,3):
return x - np.eye(3)*np.trace(x)/3.0
else:
return d - np.einsum('ijk,i->ijk',np.broadcast_to(np.eye(3),[d.shape[0],3,3]),np.trace(d,axis1=1,axis2=2)/3.0)
def spherical(x):
if np.shape(x) == (3,3):
return np.trace(x)/3.0
else:
return np.trace(x,axis1=1,axis2=2)/3.0