Merge branch 'misc-improvements' of magit1.mpie.de:/damask/DAMASK into misc-improvements

This commit is contained in:
Martin Diehl 2020-05-25 16:40:23 +02:00
commit 9d7158b51a
20 changed files with 674 additions and 326 deletions

@ -1 +1 @@
Subproject commit 72d526e5750366a9efe4d1fd9d92e0d1ecd2cd38 Subproject commit 8dde2a68538b7cffbe9d370e2b60be90a31627ab

View File

@ -1 +1 @@
v2.0.3-2504-gcee9daff v2.0.3-2514-g873b9fa8

View File

@ -2,6 +2,7 @@
import os import os
import sys import sys
from io import StringIO
from optparse import OptionParser from optparse import OptionParser
import numpy as np import numpy as np
@ -15,91 +16,82 @@ scriptID = ' '.join([scriptName,damask.version])
slipSystems = { slipSystems = {
'fcc': 'fcc':
np.array([ np.array([
# Slip direction Plane normal [+0,+1,-1 , +1,+1,+1],
[ 0, 1,-1, 1, 1, 1, ], [-1,+0,+1 , +1,+1,+1],
[-1, 0, 1, 1, 1, 1, ], [+1,-1,+0 , +1,+1,+1],
[ 1,-1, 0, 1, 1, 1, ], [+0,-1,-1 , -1,-1,+1],
[ 0,-1,-1, -1,-1, 1, ], [+1,+0,+1 , -1,-1,+1],
[ 1, 0, 1, -1,-1, 1, ], [-1,+1,+0 , -1,-1,+1],
[-1, 1, 0, -1,-1, 1, ], [+0,-1,+1 , +1,-1,-1],
[ 0,-1, 1, 1,-1,-1, ], [-1,+0,-1 , +1,-1,-1],
[-1, 0,-1, 1,-1,-1, ], [+1,+1,+0 , +1,-1,-1],
[ 1, 1, 0, 1,-1,-1, ], [+0,+1,+1 , -1,+1,-1],
[ 0, 1, 1, -1, 1,-1, ], [+1,+0,-1 , -1,+1,-1],
[ 1, 0,-1, -1, 1,-1, ], [-1,-1,+0 , -1,+1,-1],
[-1,-1, 0, -1, 1,-1, ], ],'d'),
],'f'),
'bcc': 'bcc':
np.array([ np.array([
# Slip system <111>{110} [+1,-1,+1 , +0,+1,+1],
[ 1,-1, 1, 0, 1, 1, ], [-1,-1,+1 , +0,+1,+1],
[-1,-1, 1, 0, 1, 1, ], [+1,+1,+1 , +0,-1,+1],
[ 1, 1, 1, 0,-1, 1, ], [-1,+1,+1 , +0,-1,+1],
[-1, 1, 1, 0,-1, 1, ], [-1,+1,+1 , +1,+0,+1],
[-1, 1, 1, 1, 0, 1, ], [-1,-1,+1 , +1,+0,+1],
[-1,-1, 1, 1, 0, 1, ], [+1,+1,+1 , -1,+0,+1],
[ 1, 1, 1, -1, 0, 1, ], [+1,-1,+1 , -1,+0,+1],
[ 1,-1, 1, -1, 0, 1, ], [-1,+1,+1 , +1,+1,+0],
[-1, 1, 1, 1, 1, 0, ], [-1,+1,-1 , +1,+1,+0],
[-1, 1,-1, 1, 1, 0, ], [+1,+1,+1 , -1,+1,+0],
[ 1, 1, 1, -1, 1, 0, ], [+1,+1,-1 , -1,+1,+0],
[ 1, 1,-1, -1, 1, 0, ], [-1,+1,+1 , +2,+1,+1],
# Slip system <111>{112} [+1,+1,+1 , -2,+1,+1],
[-1, 1, 1, 2, 1, 1, ], [+1,+1,-1 , +2,-1,+1],
[ 1, 1, 1, -2, 1, 1, ], [+1,-1,+1 , +2,+1,-1],
[ 1, 1,-1, 2,-1, 1, ], [+1,-1,+1 , +1,+2,+1],
[ 1,-1, 1, 2, 1,-1, ], [+1,+1,-1 , -1,+2,+1],
[ 1,-1, 1, 1, 2, 1, ], [+1,+1,+1 , +1,-2,+1],
[ 1, 1,-1, -1, 2, 1, ], [-1,+1,+1 , +1,+2,-1],
[ 1, 1, 1, 1,-2, 1, ], [+1,+1,-1 , +1,+1,+2],
[-1, 1, 1, 1, 2,-1, ], [+1,-1,+1 , -1,+1,+2],
[ 1, 1,-1, 1, 1, 2, ], [-1,+1,+1 , +1,-1,+2],
[ 1,-1, 1, -1, 1, 2, ], [+1,+1,+1 , +1,+1,-2],
[-1, 1, 1, 1,-1, 2, ], ],'d'),
[ 1, 1, 1, 1, 1,-2, ],
],'f'),
'hex': 'hex':
np.array([ np.array([
# Basal systems <11.0>{00.1} (independent of c/a-ratio, Bravais notation (4 coordinate base)) [+2,-1,-1,+0 , +0,+0,+0,+1],
[ 2, -1, -1, 0, 0, 0, 0, 1, ], [-1,+2,-1,+0 , +0,+0,+0,+1],
[-1, 2, -1, 0, 0, 0, 0, 1, ], [-1,-1,+2,+0 , +0,+0,+0,+1],
[-1, -1, 2, 0, 0, 0, 0, 1, ], [+2,-1,-1,+0 , +0,+1,-1,+0],
# 1st type prismatic systems <11.0>{10.0} (independent of c/a-ratio) [-1,+2,-1,+0 , -1,+0,+1,+0],
[ 2, -1, -1, 0, 0, 1, -1, 0, ], [-1,-1,+2,+0 , +1,-1,+0,+0],
[-1, 2, -1, 0, -1, 0, 1, 0, ], [-1,+1,+0,+0 , +1,+1,-2,+0],
[-1, -1, 2, 0, 1, -1, 0, 0, ], [+0,-1,+1,+0 , -2,+1,+1,+0],
# 2nd type prismatic systems <10.0>{11.0} -- a slip; plane normals independent of c/a-ratio [+1,+0,-1,+0 , +1,-2,+1,+0],
[ 0, 1, -1, 0, 2, -1, -1, 0, ], [-1,+2,-1,+0 , +1,+0,-1,+1],
[-1, 0, 1, 0, -1, 2, -1, 0, ], [-2,+1,+1,+0 , +0,+1,-1,+1],
[ 1, -1, 0, 0, -1, -1, 2, 0, ], [-1,-1,+2,+0 , -1,+1,+0,+1],
# 1st type 1st order pyramidal systems <11.0>{-11.1} -- plane normals depend on the c/a-ratio [+1,-2,+1,+0 , -1,+0,+1,+1],
[ 2, -1, -1, 0, 0, 1, -1, 1, ], [+2,-1,-1,+0 , +0,-1,+1,+1],
[-1, 2, -1, 0, -1, 0, 1, 1, ], [+1,+1,-2,+0 , +1,-1,+0,+1],
[-1, -1, 2, 0, 1, -1, 0, 1, ], [-2,+1,+1,+3 , +1,+0,-1,+1],
[ 1, 1, -2, 0, -1, 1, 0, 1, ], [-1,-1,+2,+3 , +1,+0,-1,+1],
[-2, 1, 1, 0, 0, -1, 1, 1, ], [-1,-1,+2,+3 , +0,+1,-1,+1],
[ 1, -2, 1, 0, 1, 0, -1, 1, ], [+1,-2,+1,+3 , +0,+1,-1,+1],
# pyramidal system: c+a slip <11.3>{-10.1} -- plane normals depend on the c/a-ratio [+1,-2,+1,+3 , -1,+1,+0,+1],
[ 2, -1, -1, 3, -1, 1, 0, 1, ], [+2,-1,-1,+3 , -1,+1,+0,+1],
[ 1, -2, 1, 3, -1, 1, 0, 1, ], [+2,-1,-1,+3 , -1,+0,+1,+1],
[-1, -1, 2, 3, 1, 0, -1, 1, ], [+1,+1,-2,+3 , -1,+0,+1,+1],
[-2, 1, 1, 3, 1, 0, -1, 1, ], [+1,+1,-2,+3 , +0,-1,+1,+1],
[-1, 2, -1, 3, 0, -1, 1, 1, ], [-1,+2,-1,+3 , +0,-1,+1,+1],
[ 1, 1, -2, 3, 0, -1, 1, 1, ], [-1,+2,-1,+3 , +1,-1,+0,+1],
[-2, 1, 1, 3, 1, -1, 0, 1, ], [-2,+1,+1,+3 , +1,-1,+0,+1],
[-1, 2, -1, 3, 1, -1, 0, 1, ], [-1,-1,+2,+3 , +1,+1,-2,+2],
[ 1, 1, -2, 3, -1, 0, 1, 1, ], [+1,-2,+1,+3 , -1,+2,-1,+2],
[ 2, -1, -1, 3, -1, 0, 1, 1, ], [+2,-1,-1,+3 , -2,+1,+1,+2],
[ 1, -2, 1, 3, 0, 1, -1, 1, ], [+1,+1,-2,+3 , -1,-1,+2,+2],
[-1, -1, 2, 3, 0, 1, -1, 1, ], [-1,+2,-1,+3 , +1,-2,+1,+2],
# pyramidal system: c+a slip <11.3>{-1-1.2} -- as for hexagonal ice (Castelnau et al. 1996, similar to twin system found below) [-2,+1,+1,+3 , +2,-1,-1,+2],
[ 2, -1, -1, 3, -2, 1, 1, 2, ], # sorted according to similar twin system ],'d'),
[-1, 2, -1, 3, 1, -2, 1, 2, ], # <11.3>{-1-1.2} shear = 2((c/a)^2-2)/(3 c/a)
[-1, -1, 2, 3, 1, 1, -2, 2, ],
[-2, 1, 1, 3, 2, -1, -1, 2, ],
[ 1, -2, 1, 3, -1, 2, -1, 2, ],
[ 1, 1, -2, 3, -1, -1, 2, 2, ],
],'f'),
} }
# -------------------------------------------------------------------- # --------------------------------------------------------------------
@ -111,11 +103,11 @@ Add columns listing Schmid factors (and optional trace vector of selected system
""", version = scriptID) """, version = scriptID)
latticeChoices = ('fcc','bcc','hex') lattice_choices = list(slipSystems.keys())
parser.add_option('-l', parser.add_option('-l',
'--lattice', '--lattice',
dest = 'lattice', type = 'choice', choices = latticeChoices, metavar='string', dest = 'lattice', type = 'choice', choices = lattice_choices, metavar='string',
help = 'type of lattice structure [%default] {}'.format(latticeChoices)) help = 'type of lattice structure [%default] {}'.format(lattice_choices))
parser.add_option('--covera', parser.add_option('--covera',
dest = 'CoverA', type = 'float', metavar = 'float', dest = 'CoverA', type = 'float', metavar = 'float',
help = 'C over A ratio for hexagonal systems [%default]') help = 'C over A ratio for hexagonal systems [%default]')
@ -138,31 +130,29 @@ parser.add_option('-o',
parser.set_defaults(force = (0.0,0.0,1.0), parser.set_defaults(force = (0.0,0.0,1.0),
quaternion='orientation', quaternion='orientation',
normal = None, normal = None,
lattice = latticeChoices[0], lattice = lattice_choices[0],
CoverA = np.sqrt(8./3.), CoverA = np.sqrt(8./3.),
) )
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
if filenames == []: filenames = [None]
force = np.array(options.force) force = np.array(options.force)/np.linalg.norm(options.force)
force /= np.linalg.norm(force)
if options.normal is not None: if options.normal is not None:
normal = np.array(options.normal) normal = np.array(options.normal)/np.linalg.norm(options.ormal)
normal /= np.linalg.norm(normal)
if abs(np.dot(force,normal)) > 1e-3: if abs(np.dot(force,normal)) > 1e-3:
parser.error('stress plane normal not orthogonal to force direction') parser.error('stress plane normal not orthogonal to force direction')
else: else:
normal = force normal = force
slip_direction = np.zeros((len(slipSystems[options.lattice]),3),'f')
slip_normal = np.zeros_like(slip_direction)
if options.lattice in ['bcc','fcc']:
if options.lattice in latticeChoices[:2]:
slip_direction = slipSystems[options.lattice][:,:3] slip_direction = slipSystems[options.lattice][:,:3]
slip_normal = slipSystems[options.lattice][:,3:] slip_normal = slipSystems[options.lattice][:,3:]
elif options.lattice == latticeChoices[2]: elif options.lattice == 'hex':
slip_direction = np.zeros((len(slipSystems['hex']),3),'d')
slip_normal = np.zeros_like(slip_direction)
# convert 4 Miller index notation of hex to orthogonal 3 Miller index notation # convert 4 Miller index notation of hex to orthogonal 3 Miller index notation
for i in range(len(slip_direction)): for i in range(len(slip_direction)):
slip_direction[i] = np.array([slipSystems['hex'][i,0]*1.5, slip_direction[i] = np.array([slipSystems['hex'][i,0]*1.5,
@ -174,52 +164,29 @@ elif options.lattice == latticeChoices[2]:
slipSystems['hex'][i,7]/options.CoverA, slipSystems['hex'][i,7]/options.CoverA,
]) ])
slip_direction /= np.tile(np.linalg.norm(slip_direction,axis=1),(3,1)).T slip_direction /= np.linalg.norm(slip_direction,axis=1,keepdims=True)
slip_normal /= np.tile(np.linalg.norm(slip_normal ,axis=1),(3,1)).T slip_normal /= np.linalg.norm(slip_normal, axis=1,keepdims=True)
# --- loop over input files ------------------------------------------------------------------------ labels = ['S[{direction[0]:.1g}_{direction[1]:.1g}_{direction[2]:.1g}]'
'({normal[0]:.1g}_{normal[1]:.1g}_{normal[2]:.1g})'\
if filenames == []: filenames = [None] .format(normal = theNormal, direction = theDirection,
) for theNormal,theDirection in zip(slip_normal,slip_direction)]
for name in filenames: for name in filenames:
try:
table = damask.ASCIItable(name = name)
except IOError:
continue
damask.util.report(scriptName,name) damask.util.report(scriptName,name)
# ------------------------------------------ read header ------------------------------------------ table = damask.Table.from_ASCII(StringIO(''.join(sys.stdin.read())) if name is None else name)
table.head_read() o = damask.Rotation.from_quaternion(table.get(options.quaternion))
# ------------------------------------------ sanity checks ---------------------------------------- force = np.broadcast_to(force, o.shape+(3,))
if not table.label_dimension(options.quaternion) == 4: normal = np.broadcast_to(normal,o.shape+(3,))
damask.util.croak('input {} does not have dimension 4.'.format(options.quaternion)) slip_direction = np.broadcast_to(slip_direction,o.shape+slip_direction.shape)
table.close(dismiss = True) # close ASCIItable and remove empty file slip_normal = np.broadcast_to(slip_normal, o.shape+slip_normal.shape)
continue S = np.abs(np.einsum('ijk,ik->ij',slip_direction,(o@force))*
np.einsum('ijk,ik->ij',slip_normal, (o@normal)))
column = table.label_index(options.quaternion) for i,label in enumerate(labels):
table.add(label,S[:,i],scriptID+' '+' '.join(sys.argv[1:]))
# ------------------------------------------ assemble header --------------------------------------- table.to_ASCII(sys.stdout if name is None else name)
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
table.labels_append(['S[{direction[0]:.1g}_{direction[1]:.1g}_{direction[2]:.1g}]'
'({normal[0]:.1g}_{normal[1]:.1g}_{normal[2]:.1g})'\
.format(normal = theNormal,
direction = theDirection,
) for theNormal,theDirection in zip(slip_normal,slip_direction)])
table.head_write()
# ------------------------------------------ process data ------------------------------------------
outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table
o = damask.Rotation(np.array(list(map(float,table.data[column:column+4]))))
table.data_append( np.abs( np.sum(slip_direction * (o * force) ,axis=1) \
* np.sum(slip_normal * (o * normal),axis=1)))
outputAlive = table.data_write() # output processed line
# ------------------------------------------ output finalization -----------------------------------
table.close() # close ASCII tables

View File

@ -40,35 +40,22 @@ parser.add_option('-f','--fill',
parser.set_defaults(offset = (0,0,0)) parser.set_defaults(offset = (0,0,0))
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
if filenames == []: filenames = [None] if filenames == []: filenames = [None]
for name in filenames: for name in filenames:
damask.util.report(scriptName,name) damask.util.report(scriptName,name)
geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name) geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name)
origin = geom.get_origin()
size = geom.get_size()
old = new = geom.get_grid()
offset = np.asarray(options.offset) offset = np.asarray(options.offset)
if options.grid is not None: if options.grid is not None:
new = np.maximum(1, grid = np.maximum(1,
np.array([int(o*float(n.lower().replace('x',''))) if n.lower().endswith('x') \ np.array([int(o*float(n.lower().replace('x',''))) if n.lower().endswith('x') \
else int(n) for o,n in zip(old,options.grid)],dtype=int)) else int(n) for o,n in zip(geom.grid,options.grid)],dtype=int))
else:
grid = np.array(options.grid,dtype=int)
canvas = np.full(new,options.fill if options.fill is not None damask.util.croak(geom.canvas(grid,np.asarray(options.offset),options.fill))
else np.nanmax(geom.microstructure)+1,geom.microstructure.dtype)
l = np.clip( offset, 0,np.minimum(old +offset,new)) # noqa
r = np.clip( offset+old,0,np.minimum(old*2+offset,new))
L = np.clip(-offset, 0,np.minimum(new -offset,old))
R = np.clip(-offset+new,0,np.minimum(new*2-offset,old))
canvas[l[0]:r[0],l[1]:r[1],l[2]:r[2]] = geom.microstructure[L[0]:R[0],L[1]:R[1],L[2]:R[2]]
damask.util.croak(geom.update(canvas,origin=origin+offset*size/old,rescale=True))
geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:])) geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:]))
geom.to_file(sys.stdout if name is None else name,pack=False) geom.to_file(sys.stdout if name is None else name,pack=False)

View File

@ -22,8 +22,6 @@ Renumber sorted microstructure indices to 1,...,N.
""", version=scriptID) """, version=scriptID)
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
if filenames == []: filenames = [None] if filenames == []: filenames = [None]
for name in filenames: for name in filenames:

View File

@ -5,7 +5,6 @@ import sys
from io import StringIO from io import StringIO
from optparse import OptionParser from optparse import OptionParser
from scipy import ndimage
import numpy as np import numpy as np
import damask import damask
@ -52,6 +51,7 @@ parser.add_option('-f', '--fill',
parser.set_defaults(degrees = False) parser.set_defaults(degrees = False)
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
if filenames == []: filenames = [None]
if [options.rotation,options.eulers,options.matrix,options.quaternion].count(None) < 3: if [options.rotation,options.eulers,options.matrix,options.quaternion].count(None) < 3:
parser.error('more than one rotation specified.') parser.error('more than one rotation specified.')
@ -67,32 +67,11 @@ if options.matrix is not None:
if options.eulers is not None: if options.eulers is not None:
rot = damask.Rotation.from_Eulers(np.array(options.eulers),degrees=options.degrees) rot = damask.Rotation.from_Eulers(np.array(options.eulers),degrees=options.degrees)
eulers = rot.as_Eulers(degrees=True)
if filenames == []: filenames = [None]
for name in filenames: for name in filenames:
damask.util.report(scriptName,name) damask.util.report(scriptName,name)
geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name) geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name)
size = geom.get_size() damask.util.croak(geom.rotate(rot,options.fill))
grid = geom.get_grid()
origin = geom.get_origin()
microstructure = geom.get_microstructure()
fill = np.nanmax(microstructure)+1 if options.fill is None else options.fill
dtype = float if np.isnan(fill) or int(fill) != fill or microstructure.dtype==np.float else int
# These rotations are always applied in the reference coordinate system, i.e. (z,x,z) not (z,x',z'')
# this seems to be ok, see https://www.cs.utexas.edu/~theshark/courses/cs354/lectures/cs354-14.pdf
microstructure = ndimage.rotate(microstructure,eulers[2],(0,1),order=0,
prefilter=False,output=dtype,cval=fill) # rotation around z
microstructure = ndimage.rotate(microstructure,eulers[1],(1,2),order=0,
prefilter=False,output=dtype,cval=fill) # rotation around x
microstructure = ndimage.rotate(microstructure,eulers[0],(0,1),order=0,
prefilter=False,output=dtype,cval=fill) # rotation around z
damask.util.croak(geom.update(microstructure,origin=origin-(np.asarray(microstructure.shape)-grid)/2*size/grid,rescale=True))
geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:])) geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:]))
geom.to_file(sys.stdout if name is None else name,pack=False) geom.to_file(sys.stdout if name is None else name,pack=False)

View File

@ -40,22 +40,16 @@ parser.set_defaults(origin = (0.0,0.0,0.0),
) )
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
if filenames == []: filenames = [None]
sub = list(map(int,options.substitute)) sub = list(map(int,options.substitute))
if filenames == []: filenames = [None]
for name in filenames: for name in filenames:
damask.util.report(scriptName,name) damask.util.report(scriptName,name)
geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name) geom = damask.Geom.from_file(StringIO(''.join(sys.stdin.read())) if name is None else name)
geom.renumber(sub[0::2],sub[1::2],origin=geom.origin+options.origin)
substituted = geom.get_microstructure() geom.microstructure+= options.microstructure
for old,new in zip(sub[0::2],sub[1::2]): substituted[geom.microstructure==old] = new # substitute microstructure indices damask.util.croak(geom)
substituted += options.microstructure # constant shift
damask.util.croak(geom.update(substituted,origin=geom.get_origin()+options.origin))
geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:])) geom.add_comments(scriptID + ' ' + ' '.join(sys.argv[1:]))
geom.to_file(sys.stdout if name is None else name,pack=False) geom.to_file(sys.stdout if name is None else name,pack=False)

View File

@ -1,4 +1,4 @@
"""Main aggregator.""" """Tools for pre and post processing of DAMASK simulations."""
import os as _os import os as _os
import re as _re import re as _re

View File

@ -581,3 +581,58 @@ class Geom:
#self.add_comments('geom.py:renumber v{}'.format(version) #self.add_comments('geom.py:renumber v{}'.format(version)
return self.update(renumbered) return self.update(renumbered)
def rotate(self,R,fill=None):
"""Rotate microstructure (pad if required)."""
if fill is None: fill = np.nanmax(self.microstructure) + 1
dtype = float if np.isnan(fill) or int(fill) != fill or self.microstructure.dtype==np.float else int
Eulers = R.as_Eulers(degrees=True)
microstructure_in = self.get_microstructure()
# These rotations are always applied in the reference coordinate system, i.e. (z,x,z) not (z,x',z'')
# see https://www.cs.utexas.edu/~theshark/courses/cs354/lectures/cs354-14.pdf
for angle,axes in zip(Eulers[::-1], [(0,1),(1,2),(0,1)]):
microstructure_out = ndimage.rotate(microstructure_in,angle,axes,order=0,
prefilter=False,output=dtype,cval=fill)
if np.prod(microstructure_in.shape) == np.prod(microstructure_out.shape):
# avoid scipy interpolation errors for rotations close to multiples of 90°
microstructure_in = np.rot90(microstructure_in,k=np.rint(angle/90.).astype(int),axes=axes)
else:
microstructure_in = microstructure_out
origin = self.origin-(np.asarray(microstructure_in.shape)-self.grid)*.5 * self.size/self.grid
#self.add_comments('geom.py:rotate v{}'.format(version)
return self.update(microstructure_in,origin=origin,rescale=True)
def canvas(self,grid=None,offset=None,fill=None):
"""Crop or enlarge/pad microstructure."""
if fill is None: fill = np.nanmax(self.microstructure) + 1
if offset is None: offset = 0
dtype = float if int(fill) != fill or self.microstructure.dtype==np.float else int
canvas = np.full(self.grid if grid is None else grid,
fill if fill is not None else np.nanmax(self.microstructure)+1,dtype)
l = np.clip( offset, 0,np.minimum(self.grid +offset,grid)) # noqa
r = np.clip( offset+self.grid,0,np.minimum(self.grid*2+offset,grid))
L = np.clip(-offset, 0,np.minimum(grid -offset,self.grid))
R = np.clip(-offset+grid, 0,np.minimum(grid*2 -offset,self.grid))
canvas[l[0]:r[0],l[1]:r[1],l[2]:r[2]] = self.microstructure[L[0]:R[0],L[1]:R[1],L[2]:R[2]]
#self.add_comments('geom.py:canvas v{}'.format(version)
return self.update(canvas,origin=self.origin+offset*self.size/self.grid,rescale=True)
def substitute(self,from_microstructure,to_microstructure):
"""Substitude microstructure indices."""
substituted = self.get_microstructure()
for from_ms,to_ms in zip(from_microstructure,to_microstructure):
substituted[self.microstructure==from_ms] = to_ms
#self.add_comments('geom.py:substitute v{}'.format(version)
return self.update(substituted)

View File

@ -229,19 +229,20 @@ class Symmetry:
Return inverse pole figure color if requested. Return inverse pole figure color if requested.
Bases are computed from Bases are computed from
basis = {'cubic' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red >>> basis = {'cubic' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red
[1.,0.,1.]/np.sqrt(2.), # direction of green ... [1.,0.,1.]/np.sqrt(2.), # direction of green
[1.,1.,1.]/np.sqrt(3.)]).T), # direction of blue ... [1.,1.,1.]/np.sqrt(3.)]).T), # direction of blue
'hexagonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red ... 'hexagonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red
[1.,0.,0.], # direction of green ... [1.,0.,0.], # direction of green
[np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T), # direction of blue ... [np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T), # direction of blue
'tetragonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red ... 'tetragonal' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red
[1.,0.,0.], # direction of green ... [1.,0.,0.], # direction of green
[1.,1.,0.]/np.sqrt(2.)]).T), # direction of blue ... [1.,1.,0.]/np.sqrt(2.)]).T), # direction of blue
'orthorhombic' : np.linalg.inv(np.array([[0.,0.,1.], # direction of red ... 'orthorhombic': np.linalg.inv(np.array([[0.,0.,1.], # direction of red
[1.,0.,0.], # direction of green ... [1.,0.,0.], # direction of green
[0.,1.,0.]]).T), # direction of blue ... [0.,1.,0.]]).T), # direction of blue
} ... }
""" """
if self.lattice == 'cubic': if self.lattice == 'cubic':
basis = {'improper':np.array([ [-1. , 0. , 1. ], basis = {'improper':np.array([ [-1. , 0. , 1. ],

View File

@ -8,6 +8,7 @@ class Orientation:
Crystallographic orientation. Crystallographic orientation.
A crystallographic orientation contains a rotation and a lattice. A crystallographic orientation contains a rotation and a lattice.
""" """
__slots__ = ['rotation','lattice'] __slots__ = ['rotation','lattice']
@ -49,8 +50,10 @@ class Orientation:
Disorientation between myself and given other orientation. Disorientation between myself and given other orientation.
Rotation axis falls into SST if SST == True. Rotation axis falls into SST if SST == True.
(Currently requires same symmetry for both orientations.
Look into A. Heinz and P. Neumann 1991 for cases with differing sym.) Currently requires same symmetry for both orientations.
Look into A. Heinz and P. Neumann 1991 for cases with differing sym.
""" """
if self.lattice.symmetry != other.lattice.symmetry: if self.lattice.symmetry != other.lattice.symmetry:
raise NotImplementedError('disorientation between different symmetry classes not supported yet.') raise NotImplementedError('disorientation between different symmetry classes not supported yet.')

View File

@ -322,9 +322,10 @@ class Result:
Return groups that contain all requested datasets. Return groups that contain all requested datasets.
Only groups within Only groups within
- inc?????/constituent/*_*/* - inc*/constituent/*/*
- inc?????/materialpoint/*_*/* - inc*/materialpoint/*/*
- inc?????/geometry/* - inc*/geometry/*
are considered as they contain user-relevant data. are considered as they contain user-relevant data.
Single strings will be treated as list with one entry. Single strings will be treated as list with one entry.
@ -599,9 +600,6 @@ class Result:
@staticmethod @staticmethod
def _add_deviator(T): def _add_deviator(T):
if not T['data'].shape[1:] == (3,3):
raise ValueError
return { return {
'data': mechanics.deviatoric_part(T['data']), 'data': mechanics.deviatoric_part(T['data']),
'label': 's_{}'.format(T['label']), 'label': 's_{}'.format(T['label']),
@ -832,14 +830,12 @@ class Result:
pole = np.array(p) pole = np.array(p)
unit_pole = pole/np.linalg.norm(pole) unit_pole = pole/np.linalg.norm(pole)
m = util.scale_to_coprime(pole) m = util.scale_to_coprime(pole)
coords = np.empty((len(q['data']),2)) rot = Rotation(q['data'].view(np.double).reshape(-1,4))
for i,qu in enumerate(q['data']):
o = Rotation(np.array([qu['w'],qu['x'],qu['y'],qu['z']]))
rotatedPole = o*unit_pole # rotate pole according to crystal orientation
(x,y) = rotatedPole[0:2]/(1.+abs(unit_pole[2])) # stereographic projection
coords[i] = [np.sqrt(x*x+y*y),np.arctan2(y,x)] if polar else [x,y]
rotatedPole = rot @ np.broadcast_to(unit_pole,rot.shape+(3,)) # rotate pole according to crystal orientation
xy = rotatedPole[:,0:2]/(1.+abs(unit_pole[2])) # stereographic projection
coords = xy if not polar else \
np.block([np.sqrt(xy[:,0:1]*xy[:,0:1]+xy[:,1:2]*xy[:,1:2]),np.arctan2(xy[:,1:2],xy[:,0:1])])
return { return {
'data': coords, 'data': coords,
'label': 'p^{}_[{} {} {})'.format(u'' if polar else 'xy',*m), 'label': 'p^{}_[{} {} {})'.format(u'' if polar else 'xy',*m),
@ -869,8 +865,6 @@ class Result:
@staticmethod @staticmethod
def _add_rotational_part(F): def _add_rotational_part(F):
if not F['data'].shape[1:] == (3,3):
raise ValueError
return { return {
'data': mechanics.rotational_part(F['data']), 'data': mechanics.rotational_part(F['data']),
'label': 'R({})'.format(F['label']), 'label': 'R({})'.format(F['label']),
@ -895,9 +889,6 @@ class Result:
@staticmethod @staticmethod
def _add_spherical(T): def _add_spherical(T):
if not T['data'].shape[1:] == (3,3):
raise ValueError
return { return {
'data': mechanics.spherical_part(T['data']), 'data': mechanics.spherical_part(T['data']),
'label': 'p_{}'.format(T['label']), 'label': 'p_{}'.format(T['label']),
@ -922,9 +913,6 @@ class Result:
@staticmethod @staticmethod
def _add_strain_tensor(F,t,m): def _add_strain_tensor(F,t,m):
if not F['data'].shape[1:] == (3,3):
raise ValueError
return { return {
'data': mechanics.strain_tensor(F['data'],t,m), 'data': mechanics.strain_tensor(F['data'],t,m),
'label': 'epsilon_{}^{}({})'.format(t,m,F['label']), 'label': 'epsilon_{}^{}({})'.format(t,m,F['label']),
@ -956,9 +944,6 @@ class Result:
@staticmethod @staticmethod
def _add_stretch_tensor(F,t): def _add_stretch_tensor(F,t):
if not F['data'].shape[1:] == (3,3):
raise ValueError
return { return {
'data': mechanics.left_stretch(F['data']) if t == 'V' else mechanics.right_stretch(F['data']), 'data': mechanics.left_stretch(F['data']) if t == 'V' else mechanics.right_stretch(F['data']),
'label': '{}({})'.format(t,F['label']), 'label': '{}({})'.format(t,F['label']),
@ -1041,7 +1026,7 @@ class Result:
pool.join() pool.join()
def write_XMDF(self): def write_XDMF(self):
""" """
Write XDMF file to directly visualize data in DADF5 file. Write XDMF file to directly visualize data in DADF5 file.
@ -1049,7 +1034,7 @@ class Result:
Selection is not taken into account. Selection is not taken into account.
""" """
if len(self.constituents) != 1 or not self.structured: if len(self.constituents) != 1 or not self.structured:
raise NotImplementedError raise NotImplementedError('XDMF only available for grid results with 1 constituent.')
xdmf=ET.Element('Xdmf') xdmf=ET.Element('Xdmf')
xdmf.attrib={'Version': '2.0', xdmf.attrib={'Version': '2.0',
@ -1217,14 +1202,6 @@ class Result:
################################################################################################### ###################################################################################################
# BEGIN DEPRECATED # BEGIN DEPRECATED
def _time_to_inc(self,start,end):
selected = []
for i,time in enumerate(self.times):
if start <= time <= end:
selected.append(self.increments[i])
return selected
def set_by_time(self,start,end): def set_by_time(self,start,end):
""" """
Set active increments based on start and end time. Set active increments based on start and end time.
@ -1237,4 +1214,4 @@ class Result:
end time (included) end time (included)
""" """
self._manage_selection('set','increments',self._time_to_inc(start,end)) self._manage_selection('set','times',self.times_in_range(start,end))

View File

@ -13,30 +13,31 @@ class Rotation:
u""" u"""
Orientation stored with functionality for conversion to different representations. Orientation stored with functionality for conversion to different representations.
The following conventions apply:
- coordinate frames are right-handed.
- a rotation angle ω is taken to be positive for a counterclockwise rotation
when viewing from the end point of the rotation axis towards the origin.
- rotations will be interpreted in the passive sense.
- Euler angle triplets are implemented using the Bunge convention,
with the angular ranges as [0, 2π],[0, π],[0, 2π].
- the rotation angle ω is limited to the interval [0, π].
- the real part of a quaternion is positive, Re(q) > 0
- P = -1 (as default).
Examples
--------
Rotate vector "a" (defined in coordinate system "A") to
coordinates "b" expressed in system "B":
- b = Q @ a
- b = np.dot(Q.asMatrix(),a)
References References
---------- ----------
D. Rowenhorst et al., Modelling and Simulation in Materials Science and Engineering 23:083501, 2015 D. Rowenhorst et al., Modelling and Simulation in Materials Science and Engineering 23:083501, 2015
https://doi.org/10.1088/0965-0393/23/8/083501 https://doi.org/10.1088/0965-0393/23/8/083501
Conventions
-----------
Convention 1: Coordinate frames are right-handed.
Convention 2: A rotation angle ω is taken to be positive for a counterclockwise rotation
when viewing from the end point of the rotation axis towards the origin.
Convention 3: Rotations will be interpreted in the passive sense.
Convention 4: Euler angle triplets are implemented using the Bunge convention,
with the angular ranges as [0, 2π],[0, π],[0, 2π].
Convention 5: The rotation angle ω is limited to the interval [0, π].
Convention 6: the real part of a quaternion is positive, Re(q) > 0
Convention 7: P = -1 (as default).
Usage
-----
Vector "a" (defined in coordinate system "A") is passively rotated
resulting in new coordinates "b" when expressed in system "B".
b = Q @ a
b = np.dot(Q.as_matrix(),a)
""" """
__slots__ = ['quaternion'] __slots__ = ['quaternion']
@ -160,10 +161,10 @@ class Rotation:
if self.shape == (): if self.shape == ():
q = np.broadcast_to(self.quaternion,shape+(4,)) q = np.broadcast_to(self.quaternion,shape+(4,))
else: else:
q = np.block([np.broadcast_to(self.quaternion[...,0:1],shape+(1,)), q = np.block([np.broadcast_to(self.quaternion[...,0:1],shape).reshape(shape+(1,)),
np.broadcast_to(self.quaternion[...,1:2],shape+(1,)), np.broadcast_to(self.quaternion[...,1:2],shape).reshape(shape+(1,)),
np.broadcast_to(self.quaternion[...,2:3],shape+(1,)), np.broadcast_to(self.quaternion[...,2:3],shape).reshape(shape+(1,)),
np.broadcast_to(self.quaternion[...,3:4],shape+(1,))]) np.broadcast_to(self.quaternion[...,3:4],shape).reshape(shape+(1,))])
return self.__class__(q) return self.__class__(q)
@ -537,7 +538,7 @@ class Rotation:
) )
# reduce Euler angles to definition range # reduce Euler angles to definition range
eu[np.abs(eu)<1.e-6] = 0.0 eu[np.abs(eu)<1.e-6] = 0.0
eu = np.where(eu<0, (eu+2.0*np.pi)%np.array([2.0*np.pi,np.pi,2.0*np.pi]),eu) eu = np.where(eu<0, (eu+2.0*np.pi)%np.array([2.0*np.pi,np.pi,2.0*np.pi]),eu) # needed?
return eu return eu
@staticmethod @staticmethod

View File

@ -120,9 +120,9 @@ class VTK:
Parameters Parameters
---------- ----------
fname : str fname : str
Filename for reading. Valid extensions are *.vtr, *.vtu, *.vtp, and *.vtk. Filename for reading. Valid extensions are .vtr, .vtu, .vtp, and .vtk.
dataset_type : str, optional dataset_type : str, optional
Name of the vtk.vtkDataSet subclass when opening an *.vtk file. Valid types are vtkRectilinearGrid, Name of the vtk.vtkDataSet subclass when opening an .vtk file. Valid types are vtkRectilinearGrid,
vtkUnstructuredGrid, and vtkPolyData. vtkUnstructuredGrid, and vtkPolyData.
""" """

View File

@ -0,0 +1,85 @@
4 header
grid a 8 b 10 c 8
size x 8e-06 y 1e-05 z 8e-06
origin x 0.0 y -2.5e-06 z -2e-06
homogenization 1
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 30 14 18 42 42 42
42 42 29 13 17 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 31 15 19 42 42 42
42 6 10 2 2 42 42 42
42 1 2 2 2 2 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 32 16 20 42 42 42
42 7 11 2 2 42 42 42
42 7 11 2 2 42 42 42
42 2 2 2 2 2 42 42
42 42 2 2 31 35 42 42
42 42 22 26 10 14 1 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 4 2 2 2 2 42 42
42 42 2 2 32 36 42 42
42 42 24 28 12 16 1 42
42 42 42 7 7 1 1 42
42 42 42 7 7 1 1 42
42 42 42 1 1 1 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 25 29 13 17 1 42
42 42 42 8 8 1 1 42
42 42 42 1 1 1 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 1 1 1 42 42
42 42 42 1 1 1 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42

View File

@ -0,0 +1,104 @@
4 header
grid a 11 b 11 c 9
size x 1.1e-05 y 1.1000000000000001e-05 z 9e-06
origin x -1.5e-06 y -3.0000000000000005e-06 z -2.4999999999999998e-06
homogenization 1
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 1 2 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 1 42 42 42 42 42 42
42 42 42 42 1 5 42 42 42 42 42
42 42 42 1 7 4 42 42 42 42 42
42 42 42 42 42 27 42 42 42 42 42
42 42 42 42 42 42 2 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 1 1 42 42 42 42 42 42
42 42 42 1 1 9 29 42 42 42 42
42 42 1 1 11 8 28 2 42 42 42
42 42 42 1 10 31 2 42 42 42 42
42 42 42 42 30 2 2 2 42 42 42
42 42 42 42 42 42 2 1 42 42 42
42 42 42 42 42 42 42 1 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 1 42 42 42 42 42 42 42
42 42 42 1 1 42 42 42 42 42 42
42 42 1 16 36 12 32 42 42 42 42
42 42 42 15 35 2 2 2 42 42 42
42 42 42 42 2 2 2 11 3 42 42
42 42 42 42 42 42 10 6 42 42 42
42 42 42 42 42 42 42 6 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 1 42 42 42 42 42 42 42
42 42 42 1 17 42 42 42 42 42 42
42 42 42 16 36 37 2 42 42 42 42
42 42 42 42 39 2 2 12 42 42 42
42 42 42 38 2 2 2 11 8 42 42
42 42 42 42 2 2 14 30 42 42 42
42 42 42 42 42 42 13 30 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 40 42 42 42 42 42 42
42 42 42 42 42 2 42 42 42 42 42
42 42 42 42 42 2 2 15 42 42 42
42 42 42 42 42 2 18 42 42 42 42
42 42 42 42 42 42 17 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 2 20 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42
42 42 42 42 42 42 42 42 42 42 42

View File

@ -5,12 +5,13 @@ import pytest
import numpy as np import numpy as np
from damask import Geom from damask import Geom
from damask import Rotation
def geom_equal(a,b): def geom_equal(a,b):
return np.all(a.get_microstructure() == b.get_microstructure()) and \ return np.all(a.get_microstructure() == b.get_microstructure()) and \
np.all(a.get_size() == b.get_size()) and \ np.all(a.get_grid() == b.get_grid()) and \
np.all(a.get_grid() == b.get_grid()) np.allclose(a.get_size(), b.get_size())
@pytest.fixture @pytest.fixture
def default(): def default():
@ -36,6 +37,7 @@ class TestGeom:
default.get_size(), default.get_size(),
default.get_origin() default.get_origin()
) )
print(modified)
assert geom_equal(modified,default) assert geom_equal(modified,default)
@ -57,6 +59,22 @@ class TestGeom:
new = Geom.from_file(tmpdir.join('default.geom')) new = Geom.from_file(tmpdir.join('default.geom'))
assert geom_equal(new,default) assert geom_equal(new,default)
def test_invalid_combination(self,default):
with pytest.raises(ValueError):
default.update(default.microstructure[1:,1:,1:],size=np.ones(3), rescale=True)
def test_invalid_size(self,default):
with pytest.raises(ValueError):
default.update(default.microstructure[1:,1:,1:],size=np.ones(2))
def test_invalid_microstructure(self,default):
with pytest.raises(ValueError):
default.update(default.microstructure[1])
def test_invalid_homogenization(self,default):
with pytest.raises(TypeError):
default.set_homogenization(homogenization=0)
@pytest.mark.parametrize('directions,reflect',[ @pytest.mark.parametrize('directions,reflect',[
(['x'], False), (['x'], False),
(['x','y','z'],True), (['x','y','z'],True),
@ -98,6 +116,52 @@ class TestGeom:
if update: modified.to_file(reference) if update: modified.to_file(reference)
assert geom_equal(modified,Geom.from_file(reference)) assert geom_equal(modified,Geom.from_file(reference))
def test_renumber(self,default):
modified = copy.deepcopy(default)
microstructure = modified.get_microstructure()
for m in np.unique(microstructure):
microstructure[microstructure==m] = microstructure.max() + np.random.randint(1,30)
modified.update(microstructure)
assert not geom_equal(modified,default)
modified.renumber()
assert geom_equal(modified,default)
def test_substitute(self,default):
modified = copy.deepcopy(default)
microstructure = modified.get_microstructure()
offset = np.random.randint(1,500)
microstructure+=offset
modified.update(microstructure)
assert not geom_equal(modified,default)
modified.substitute(np.arange(default.microstructure.max())+1+offset,
np.arange(default.microstructure.max())+1)
assert geom_equal(modified,default)
@pytest.mark.parametrize('axis_angle',[np.array([1,0,0,86.7]), np.array([0,1,0,90.4]), np.array([0,0,1,90]),
np.array([1,0,0,175]),np.array([0,-1,0,178]),np.array([0,0,1,180])])
def test_rotate360(self,default,axis_angle):
modified = copy.deepcopy(default)
for i in range(np.rint(360/axis_angle[3]).astype(int)):
modified.rotate(Rotation.from_axis_angle(axis_angle,degrees=True))
assert geom_equal(modified,default)
@pytest.mark.parametrize('Eulers',[[32.0,68.0,21.0],
[0.0,32.0,240.0]])
def test_rotate(self,default,update,reference_dir,Eulers):
modified = copy.deepcopy(default)
modified.rotate(Rotation.from_Eulers(Eulers,degrees=True))
tag = 'Eulers={}-{}-{}'.format(*Eulers)
reference = os.path.join(reference_dir,'rotate_{}.geom'.format(tag))
if update: modified.to_file(reference)
assert geom_equal(modified,Geom.from_file(reference))
def test_canvas(self,default):
grid_add = np.random.randint(0,30,(3))
modified = copy.deepcopy(default)
modified.canvas(modified.grid + grid_add)
e = default.grid
assert np.all(modified.microstructure[:e[0],:e[1],:e[2]] == default.microstructure)
@pytest.mark.parametrize('periodic',[True,False]) @pytest.mark.parametrize('periodic',[True,False])
def test_tessellation_approaches(self,periodic): def test_tessellation_approaches(self,periodic):
grid = np.random.randint(10,20,3) grid = np.random.randint(10,20,3)

View File

@ -0,0 +1,41 @@
import random
import pytest
import numpy as np
from damask import Symmetry
class TestSymmetry:
@pytest.mark.parametrize('invalid_symmetry',['fcc','bcc','hello'])
def test_invalid_symmetry(self,invalid_symmetry):
with pytest.raises(KeyError):
s = Symmetry(invalid_symmetry) # noqa
def test_equal(self):
symmetry = random.choice(Symmetry.lattices)
print(symmetry)
assert Symmetry(symmetry) == Symmetry(symmetry)
def test_not_equal(self):
symmetries = random.sample(Symmetry.lattices,k=2)
assert Symmetry(symmetries[0]) != Symmetry(symmetries[1])
@pytest.mark.parametrize('lattice',Symmetry.lattices)
def test_inFZ(self,lattice):
assert Symmetry(lattice).inFZ(np.zeros(3))
@pytest.mark.parametrize('lattice',Symmetry.lattices)
def test_inDisorientationSST(self,lattice):
assert Symmetry(lattice).inDisorientationSST(np.zeros(3))
@pytest.mark.parametrize('lattice',Symmetry.lattices)
@pytest.mark.parametrize('proper',[True,False])
def test_inSST(self,lattice,proper):
assert Symmetry(lattice).inSST(np.zeros(3),proper)
@pytest.mark.parametrize('function',['inFZ','inDisorientationSST'])
def test_invalid_argument(self,function):
s = Symmetry() # noqa
with pytest.raises(ValueError):
eval('s.{}(np.ones(4))'.format(function))

View File

@ -4,6 +4,7 @@ import os
import pytest import pytest
import numpy as np import numpy as np
import damask
from damask import Result from damask import Result
from damask import mechanics from damask import mechanics
@ -13,7 +14,7 @@ def default(tmp_path,reference_dir):
fname = '12grains6x7x8_tensionY.hdf5' fname = '12grains6x7x8_tensionY.hdf5'
shutil.copy(os.path.join(reference_dir,fname),tmp_path) shutil.copy(os.path.join(reference_dir,fname),tmp_path)
f = Result(os.path.join(tmp_path,fname)) f = Result(os.path.join(tmp_path,fname))
f.set_by_time(20.0,20.0) f.pick('times',20.0)
return f return f
@pytest.fixture @pytest.fixture
@ -28,12 +29,60 @@ class TestResult:
print(default) print(default)
def test_time_increments(self,default): def test_pick_all(self,default):
shape = default.read_dataset(default.get_dataset_location('F'),0).shape default.pick('increments',True)
default.set_by_time(0.0,20.0) a = default.get_dataset_location('F')
for i in default.iterate('increments'): default.pick('increments','*')
assert shape == default.read_dataset(default.get_dataset_location('F'),0).shape b = default.get_dataset_location('F')
default.pick('increments',default.incs_in_range(0,np.iinfo(int).max))
c = default.get_dataset_location('F')
default.pick('times',True)
d = default.get_dataset_location('F')
default.pick('times','*')
e = default.get_dataset_location('F')
default.pick('times',default.times_in_range(0.0,np.inf))
f = default.get_dataset_location('F')
assert a == b == c == d == e ==f
@pytest.mark.parametrize('what',['increments','times','constituents']) # ToDo: discuss materialpoints
def test_pick_none(self,default,what):
default.pick(what,False)
a = default.get_dataset_location('F')
default.pick(what,[])
b = default.get_dataset_location('F')
assert a == b == []
@pytest.mark.parametrize('what',['increments','times','constituents']) # ToDo: discuss materialpoints
def test_pick_more(self,default,what):
default.pick(what,False)
default.pick_more(what,'*')
a = default.get_dataset_location('F')
default.pick(what,True)
b = default.get_dataset_location('F')
assert a == b
@pytest.mark.parametrize('what',['increments','times','constituents']) # ToDo: discuss materialpoints
def test_pick_less(self,default,what):
default.pick(what,True)
default.pick_less(what,'*')
a = default.get_dataset_location('F')
default.pick(what,False)
b = default.get_dataset_location('F')
assert a == b == []
def test_pick_invalid(self,default):
with pytest.raises(AttributeError):
default.pick('invalid',True)
def test_add_invalid(self,default):
with pytest.raises(Exception):
default.add_calculation('#invalid#*2')
def test_add_absolute(self,default): def test_add_absolute(self,default):
default.add_absolute('Fe') default.add_absolute('Fe')
@ -95,6 +144,20 @@ class TestResult:
in_file = default.read_dataset(loc['v(sigma)'],0) in_file = default.read_dataset(loc['v(sigma)'],0)
assert np.allclose(in_memory,in_file) assert np.allclose(in_memory,in_file)
@pytest.mark.parametrize('d',[[1,0,0],[0,1,0],[0,0,1]])
def test_add_IPFcolor(self,default,d):
default.add_IPFcolor('orientation',d)
loc = {'orientation': default.get_dataset_location('orientation'),
'color': default.get_dataset_location('IPFcolor_[{} {} {}]'.format(*d))}
qu = default.read_dataset(loc['orientation']).view(np.double).reshape(-1,4)
crystal_structure = default.get_crystal_structure()
in_memory = np.empty((qu.shape[0],3),np.uint8)
for i,q in enumerate(qu):
o = damask.Orientation(q,crystal_structure).reduced()
in_memory[i] = np.uint8(o.IPFcolor(np.array(d))*255)
in_file = default.read_dataset(loc['color'])
assert np.allclose(in_memory,in_file)
def test_add_maximum_shear(self,default): def test_add_maximum_shear(self,default):
default.add_Cauchy('P','F') default.add_Cauchy('P','F')
default.add_maximum_shear('sigma') default.add_maximum_shear('sigma')
@ -143,6 +206,20 @@ class TestResult:
in_file = default.read_dataset(loc['S'],0) in_file = default.read_dataset(loc['S'],0)
assert np.allclose(in_memory,in_file) assert np.allclose(in_memory,in_file)
@pytest.mark.parametrize('polar',[True,False])
def test_add_pole(self,default,polar):
pole = np.array([1.,0.,0.])
default.add_pole('orientation',pole,polar)
loc = {'orientation': default.get_dataset_location('orientation'),
'pole': default.get_dataset_location('p^{}_[1 0 0)'.format(u'' if polar else 'xy'))}
rot = damask.Rotation(default.read_dataset(loc['orientation']).view(np.double))
rotated_pole = rot * np.broadcast_to(pole,rot.shape+(3,))
xy = rotated_pole[:,0:2]/(1.+abs(pole[2]))
in_memory = xy if not polar else \
np.block([np.sqrt(xy[:,0:1]*xy[:,0:1]+xy[:,1:2]*xy[:,1:2]),np.arctan2(xy[:,1:2],xy[:,0:1])])
in_file = default.read_dataset(loc['pole'])
assert np.allclose(in_memory,in_file)
def test_add_rotational_part(self,default): def test_add_rotational_part(self,default):
default.add_rotational_part('F') default.add_rotational_part('F')
loc = {'F': default.get_dataset_location('F'), loc = {'F': default.get_dataset_location('F'),
@ -185,3 +262,7 @@ class TestResult:
in_memory = mechanics.left_stretch(default.read_dataset(loc['F'],0)) in_memory = mechanics.left_stretch(default.read_dataset(loc['F'],0))
in_file = default.read_dataset(loc['V(F)'],0) in_file = default.read_dataset(loc['V(F)'],0)
assert np.allclose(in_memory,in_file) assert np.allclose(in_memory,in_file)
@pytest.mark.parametrize('output',['F',[],['F','P']])
def test_vtk(self,default,output):
default.to_vtk(output)

View File

@ -866,12 +866,23 @@ class TestRotation:
with pytest.raises(ValueError): with pytest.raises(ValueError):
function(invalid_shape) function(invalid_shape)
@pytest.mark.parametrize('shape',[None,(3,),(4,2)])
def test_broadcast(self,shape):
rot = Rotation.from_random(shape)
new_shape = tuple(np.random.randint(8,32,(3))) if shape is None else \
rot.shape + (np.random.randint(8,32),)
rot_broadcast = rot.broadcast_to(tuple(new_shape))
for i in range(rot_broadcast.shape[-1]):
assert np.allclose(rot_broadcast.quaternion[...,i,:], rot.quaternion)
@pytest.mark.parametrize('function,invalid',[(Rotation.from_quaternion, np.array([-1,0,0,0])), @pytest.mark.parametrize('function,invalid',[(Rotation.from_quaternion, np.array([-1,0,0,0])),
(Rotation.from_quaternion, np.array([1,1,1,0])), (Rotation.from_quaternion, np.array([1,1,1,0])),
(Rotation.from_Eulers, np.array([1,4,0])), (Rotation.from_Eulers, np.array([1,4,0])),
(Rotation.from_axis_angle, np.array([1,0,0,4])), (Rotation.from_axis_angle, np.array([1,0,0,4])),
(Rotation.from_axis_angle, np.array([1,1,0,1])), (Rotation.from_axis_angle, np.array([1,1,0,1])),
(Rotation.from_matrix, np.random.rand(3,3)), (Rotation.from_matrix, np.random.rand(3,3)),
(Rotation.from_matrix, np.array([[1,1,0],[1,2,0],[0,0,1]])),
(Rotation.from_Rodrigues, np.array([1,0,0,-1])), (Rotation.from_Rodrigues, np.array([1,0,0,-1])),
(Rotation.from_Rodrigues, np.array([1,1,0,1])), (Rotation.from_Rodrigues, np.array([1,1,0,1])),
(Rotation.from_homochoric, np.array([2,2,2])), (Rotation.from_homochoric, np.array([2,2,2])),