Merge branch 'development' into 36-faster-file-handling-for-material-config-use-stream-access-instead-of-line-wise-reading

This commit is contained in:
Martin Diehl 2018-08-03 08:57:44 +02:00
commit 38183fdbf6
36 changed files with 174 additions and 214 deletions

View File

@ -12,21 +12,38 @@ echo + Send to damask@mpie.de for support
echo + view with \'cat $OUTFILE\' echo + view with \'cat $OUTFILE\'
echo =========================================== echo ===========================================
function firstLevel {
echo -e '\n\n=============================================================================================='
echo $1
echo ==============================================================================================
}
function secondLevel {
echo ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
echo $1
echo ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
}
function thirdLevel {
echo -e '\n----------------------------------------------------------------------------------------------'
echo $1
echo ----------------------------------------------------------------------------------------------
}
function getDetails { function getDetails {
if which $1 &> /dev/null; then if which $1 &> /dev/null; then
echo ---------------------------------------------------------------------------------------------- secondLevel $1:
echo $1:
echo ----------------------------------------------------------------------------------------------
echo + location: echo + location:
which $1 which $1
echo + $1 $2: echo + $1 $2:
$1 $2 $1 $2
echo -e '\n'
else else
echo $1 not found echo $1 not found
fi fi
echo
} }
# redirect STDOUT and STDERR to logfile # redirect STDOUT and STDERR to logfile
# https://stackoverflow.com/questions/11229385/redirect-all-output-in-a-bash-script-when-using-set-x^ # https://stackoverflow.com/questions/11229385/redirect-all-output-in-a-bash-script-when-using-set-x^
exec > $OUTFILE 2>&1 exec > $OUTFILE 2>&1
@ -38,28 +55,18 @@ DAMASK_ROOT="$( cd "$( dirname "${BASH_SOURCE[0]}" )" && pwd )"
echo XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX echo XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
echo System report for \'$(hostname)\' created on $(date '+%Y-%m-%d %H:%M:%S') by \'$(whoami)\' echo System report for \'$(hostname)\' created on $(date '+%Y-%m-%d %H:%M:%S') by \'$(whoami)\'
echo XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX echo XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX
echo
echo ============================================================================================== firstLevel "DAMASK settings"
echo DAMASK settings secondLevel "DAMASK_ROOT"
echo ==============================================================================================
echo ----------------------------------------------------------------------------------------------
echo DAMASK_ROOT:
echo ----------------------------------------------------------------------------------------------
echo $DAMASK_ROOT echo $DAMASK_ROOT
echo echo
echo ---------------------------------------------------------------------------------------------- secondLevel "Version"
echo Version:
echo ----------------------------------------------------------------------------------------------
cat VERSION cat VERSION
echo echo
echo ---------------------------------------------------------------------------------------------- secondLevel "Settings in CONFIG"
echo Settings in CONFIG:
echo ----------------------------------------------------------------------------------------------
cat CONFIG cat CONFIG
echo
echo ============================================================================================== firstLevel "System"
echo System
echo ==============================================================================================
uname -a uname -a
echo echo
echo PATH: $PATH echo PATH: $PATH
@ -69,74 +76,52 @@ echo SHELL: $SHELL
echo PETSC_ARCH: $PETSC_ARCH echo PETSC_ARCH: $PETSC_ARCH
echo PETSC_DIR: $PETSC_DIR echo PETSC_DIR: $PETSC_DIR
ls $PETSC_DIR/lib ls $PETSC_DIR/lib
echo
echo ==============================================================================================
echo Python
echo ==============================================================================================
firstLevel "Python"
DEFAULT_PYTHON=python2.7 DEFAULT_PYTHON=python2.7
for executable in python python2 python3 python2.7; do for executable in python python2 python3 python2.7; do
getDetails $executable '--version' getDetails $executable '--version'
done done
echo ---------------------------------------------------------------------------------------------- secondLevel "Details on $DEFAULT_PYTHON:"
echo Details on $DEFAULT_PYTHON:
echo ----------------------------------------------------------------------------------------------
echo $(ls -la $(which $DEFAULT_PYTHON)) echo $(ls -la $(which $DEFAULT_PYTHON))
for module in numpy scipy;do for module in numpy scipy;do
echo -e '\n----------------------------------------------------------------------------------------------' thirdLevel $module
echo $module
echo ----------------------------------------------------------------------------------------------
$DEFAULT_PYTHON -c "import $module; \ $DEFAULT_PYTHON -c "import $module; \
print('Version: {}'.format($module.__version__)); \ print('Version: {}'.format($module.__version__)); \
print('Location: {}'.format($module.__file__))" print('Location: {}'.format($module.__file__))"
done done
echo ---------------------------------------------------------------------------------------------- thirdLevel vtk
echo vtk
echo ----------------------------------------------------------------------------------------------
$DEFAULT_PYTHON -c "import vtk; \ $DEFAULT_PYTHON -c "import vtk; \
print('Version: {}'.format(vtk.vtkVersion.GetVTKVersion())); \ print('Version: {}'.format(vtk.vtkVersion.GetVTKVersion())); \
print('Location: {}'.format(vtk.__file__))" print('Location: {}'.format(vtk.__file__))"
echo ---------------------------------------------------------------------------------------------- thirdLevel h5py
echo h5py
echo ----------------------------------------------------------------------------------------------
$DEFAULT_PYTHON -c "import h5py; \ $DEFAULT_PYTHON -c "import h5py; \
print('Version: {}'.format(h5py.version.version)); \ print('Version: {}'.format(h5py.version.version)); \
print('Location: {}'.format(h5py.__file__))" print('Location: {}'.format(h5py.__file__))"
echo
echo ============================================================================================== firstLevel "GNU Compiler Collection"
echo GCC
echo ==============================================================================================
for executable in gcc g++ gfortran ;do for executable in gcc g++ gfortran ;do
getDetails $executable '--version' getDetails $executable '--version'
done done
echo
echo ============================================================================================== firstLevel "Intel Compiler Suite"
echo Intel Compiler Suite
echo ==============================================================================================
for executable in icc icpc ifort ;do for executable in icc icpc ifort ;do
getDetails $executable '--version' getDetails $executable '--version'
done done
echo
echo ============================================================================================== firstLevel "MPI Wrappers"
echo MPI Wrappers
echo ==============================================================================================
for executable in mpicc mpiCC mpic++ mpicpc mpicxx mpifort mpif90 mpif77; do for executable in mpicc mpiCC mpic++ mpicpc mpicxx mpifort mpif90 mpif77; do
getDetails $executable '-show' getDetails $executable '-show'
done done
echo
echo ============================================================================================== firstLevel "MPI Launchers"
echo MPI Launchers
echo ==============================================================================================
for executable in mpirun mpiexec; do for executable in mpirun mpiexec; do
getDetails $executable '--version' getDetails $executable '--version'
done done
echo
echo ============================================================================================== firstLevel "Abaqus"
echo Abaqus
echo ==============================================================================================
cd installation/mods_Abaqus # to have the right environment file cd installation/mods_Abaqus # to have the right environment file
for executable in abaqus abq2016 abq2017; do for executable in abaqus abq2017 abq2018; do
getDetails $executable 'information=all' getDetails $executable 'information=all'
done done
cd ../.. cd ../..

@ -1 +1 @@
Subproject commit d1d46580823d2896059b9514ddc975f9fe5f6b1f Subproject commit 50eb21714e2f501b111bb62096ebb6a5bfc6708a

View File

@ -1 +1 @@
v2.0.2-228-ge08b7325 v2.0.2-261-gbc3f6ae9

View File

@ -0,0 +1,8 @@
#! /usr/bin/env bash
if [ $1x != 3to2x ]; then
echo 'python2.7 to python'
find . -name '*.py' | xargs sed -i 's/usr\/bin\/env python2.7/usr\/bin\/env python/g'
else
echo 'python to python2.7'
find . -name '*.py' | xargs sed -i 's/usr\/bin\/env python/usr\/bin\/env python2.7/g'
fi

View File

@ -36,8 +36,8 @@ class bcolors:
def srepr(arg,glue = '\n'): def srepr(arg,glue = '\n'):
"""Joins arguments as individual lines""" """Joins arguments as individual lines"""
if (not hasattr(arg, "strip") and if (not hasattr(arg, "strip") and
hasattr(arg, "__getitem__") or (hasattr(arg, "__getitem__") or
hasattr(arg, "__iter__")): hasattr(arg, "__iter__"))):
return glue.join(str(x) for x in arg) return glue.join(str(x) for x in arg)
return arg if isinstance(arg,str) else repr(arg) return arg if isinstance(arg,str) else repr(arg)
@ -59,9 +59,9 @@ def report_geom(info,
what = ['grid','size','origin','homogenization','microstructures']): what = ['grid','size','origin','homogenization','microstructures']):
"""Reports (selected) geometry information""" """Reports (selected) geometry information"""
output = { output = {
'grid' : 'grid a b c: {}'.format(' x '.join(map(str,info['grid' ]))), 'grid' : 'grid a b c: {}'.format(' x '.join(list(map(str,info['grid' ])))),
'size' : 'size x y z: {}'.format(' x '.join(map(str,info['size' ]))), 'size' : 'size x y z: {}'.format(' x '.join(list(map(str,info['size' ])))),
'origin' : 'origin x y z: {}'.format(' : '.join(map(str,info['origin']))), 'origin' : 'origin x y z: {}'.format(' : '.join(list(map(str,info['origin'])))),
'homogenization' : 'homogenization: {}'.format(info['homogenization']), 'homogenization' : 'homogenization: {}'.format(info['homogenization']),
'microstructures' : 'microstructures: {}'.format(info['microstructures']), 'microstructures' : 'microstructures: {}'.format(info['microstructures']),
} }
@ -93,8 +93,10 @@ def execute(cmd,
stdout = subprocess.PIPE, stdout = subprocess.PIPE,
stderr = subprocess.PIPE, stderr = subprocess.PIPE,
stdin = subprocess.PIPE) stdin = subprocess.PIPE)
out,error = [i.replace(b"\x08",b"") for i in (process.communicate() if streamIn is None out,error = [i for i in (process.communicate() if streamIn is None
else process.communicate(streamIn.read()))] else process.communicate(streamIn.read().encode('utf-8')))]
out = out.decode('utf-8').replace('\x08','')
error = error.decode('utf-8').replace('\x08','')
os.chdir(initialPath) os.chdir(initialPath)
if process.returncode != 0: raise RuntimeError('{} failed with returncode {}'.format(cmd,process.returncode)) if process.returncode != 0: raise RuntimeError('{} failed with returncode {}'.format(cmd,process.returncode))
return out,error return out,error
@ -103,9 +105,9 @@ def coordGridAndSize(coordinates):
"""Determines grid count and overall physical size along each dimension of an ordered array of coordinates""" """Determines grid count and overall physical size along each dimension of an ordered array of coordinates"""
dim = coordinates.shape[1] dim = coordinates.shape[1]
coords = [np.unique(coordinates[:,i]) for i in range(dim)] coords = [np.unique(coordinates[:,i]) for i in range(dim)]
mincorner = np.array(map(min,coords)) mincorner = np.array(list(map(min,coords)))
maxcorner = np.array(map(max,coords)) maxcorner = np.array(list(map(max,coords)))
grid = np.array(map(len,coords),'i') grid = np.array(list(map(len,coords)),'i')
size = grid/np.maximum(np.ones(dim,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1) size = grid/np.maximum(np.ones(dim,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1)
size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 equal to smallest among other ones size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 equal to smallest among other ones
return grid,size return grid,size

View File

@ -75,8 +75,8 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
F = np.array(map(float,table.data[column[options.defgrad]:column[options.defgrad]+9]),'d').reshape(3,3) F = np.array(list(map(float,table.data[column[options.defgrad]:column[options.defgrad]+9])),'d').reshape(3,3)
P = np.array(map(float,table.data[column[options.stress ]:column[options.stress ]+9]),'d').reshape(3,3) P = np.array(list(map(float,table.data[column[options.stress ]:column[options.stress ]+9])),'d').reshape(3,3)
table.data_append(list(1.0/np.linalg.det(F)*np.dot(P,F.T).reshape(9))) # [Cauchy] = (1/det(F)) * [P].[F_transpose] table.data_append(list(1.0/np.linalg.det(F)*np.dot(P,F.T).reshape(9))) # [Cauchy] = (1/det(F)) * [P].[F_transpose]
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line

View File

@ -282,19 +282,12 @@ for name in filenames:
table.data_readArray([options.defgrad,options.pos]) table.data_readArray([options.defgrad,options.pos])
table.data_rewind() table.data_rewind()
if len(table.data.shape) < 2: table.data.shape += (1,) # expand to 2D shape
if table.data[:,9:].shape[1] < 3: if table.data[:,9:].shape[1] < 3:
table.data = np.hstack((table.data, table.data = np.hstack((table.data,
np.zeros((table.data.shape[0], np.zeros((table.data.shape[0],
3-table.data[:,9:].shape[1]),dtype='f'))) # fill coords up to 3D with zeros 3-table.data[:,9:].shape[1]),dtype='f'))) # fill coords up to 3D with zeros
coords = [np.unique(table.data[:,9+i]) for i in range(3)] grid,size = damask.util.coordGridAndSize(table.data[:,9:12])
mincorner = np.array(map(min,coords))
maxcorner = np.array(map(max,coords))
grid = np.array(map(len,coords),'i')
size = grid/np.maximum(np.ones(3,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1)
size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 set to smallest among other spacings
N = grid.prod() N = grid.prod()
if N != len(table.data): errors.append('data count {} does not match grid {}x{}x{}.'.format(N,*grid)) if N != len(table.data): errors.append('data count {} does not match grid {}x{}x{}.'.format(N,*grid))

View File

@ -138,7 +138,6 @@ for name in filenames:
# --------------- figure out size and grid --------------------------------------------------------- # --------------- figure out size and grid ---------------------------------------------------------
table.data_readArray() table.data_readArray()
grid,size = damask.util.coordGridAndSize(table.data[:,table.label_indexrange(options.pos)]) grid,size = damask.util.coordGridAndSize(table.data[:,table.label_indexrange(options.pos)])
# ------------------------------------------ process value field ----------------------------------- # ------------------------------------------ process value field -----------------------------------

View File

@ -58,7 +58,7 @@ for name in filenames:
errors = [] errors = []
remarks = [] remarks = []
for type, data in items.iteritems(): for type, data in items.items():
for what in data['labels']: for what in data['labels']:
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type))
@ -81,12 +81,11 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
for type, data in items.iteritems(): for type, data in items.items():
for column in data['column']: for column in data['column']:
table.data_append(determinant(map(float,table.data[column: table.data_append(determinant(list(map(float,table.data[column: column+data['dim']]))))
column+data['dim']])))
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line
# ------------------------------------------ output finalization ----------------------------------- # ------------------------------------------ output finalization -----------------------------------
table.close() # close input ASCII table (works for stdin) table.close() # close input ASCII table (works for stdin)

View File

@ -66,7 +66,7 @@ for name in filenames:
remarks = [] remarks = []
column = {} column = {}
for type, data in items.iteritems(): for type, data in items.items():
for what in data['labels']: for what in data['labels']:
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}.'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}.'.format(what,type))
@ -83,7 +83,7 @@ for name in filenames:
# ------------------------------------------ assemble header -------------------------------------- # ------------------------------------------ assemble header --------------------------------------
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:])) table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
for type, data in items.iteritems(): for type, data in items.items():
for label in data['active']: for label in data['active']:
table.labels_append(['{}_dev({})'.format(i+1,label) for i in range(data['dim'])] + \ table.labels_append(['{}_dev({})'.format(i+1,label) for i in range(data['dim'])] + \
(['sph({})'.format(label)] if options.spherical else [])) # extend ASCII header with new labels (['sph({})'.format(label)] if options.spherical else [])) # extend ASCII header with new labels
@ -93,10 +93,10 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
for type, data in items.iteritems(): for type, data in items.items():
for column in data['column']: for column in data['column']:
table.data_append(deviator(map(float,table.data[column: table.data_append(deviator(list(map(float,table.data[column:
column+data['dim']]),options.spherical)) column+data['dim']])),options.spherical))
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line
# ------------------------------------------ output finalization ----------------------------------- # ------------------------------------------ output finalization -----------------------------------

View File

@ -168,13 +168,7 @@ for name in filenames:
np.zeros((table.data.shape[0], np.zeros((table.data.shape[0],
3-table.data[:,9:].shape[1]),dtype='f'))) # fill coords up to 3D with zeros 3-table.data[:,9:].shape[1]),dtype='f'))) # fill coords up to 3D with zeros
coords = [np.unique(table.data[:,9+i]) for i in range(3)] grid,size = damask.util.coordGridAndSize(table.data[:,9:12])
mincorner = np.array(map(min,coords))
maxcorner = np.array(map(max,coords))
grid = np.array(map(len,coords),'i')
size = grid/np.maximum(np.ones(3,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1)
size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 set to smallest among other spacings
N = grid.prod() N = grid.prod()
if N != len(table.data): errors.append('data count {} does not match grid {}x{}x{}.'.format(N,*grid)) if N != len(table.data): errors.append('data count {} does not match grid {}x{}x{}.'.format(N,*grid))

View File

@ -88,9 +88,9 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
for column in columns: for column in columns:
table.data_append(E_hkl(map(float,table.data[column:column+3]),options.hkl)) table.data_append(E_hkl(list(map(float,table.data[column:column+3])),options.hkl))
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line
# ------------------------------------------ output finalization ----------------------------------- # ------------------------------------------ output finalization -----------------------------------
table.close() # close ASCII tables table.close() # close ASCII tables

View File

@ -102,7 +102,7 @@ parser.add_option('-t',
help = 'feature type {{{}}} '.format(', '.join(map(lambda x:'/'.join(x['names']),features))) ) help = 'feature type {{{}}} '.format(', '.join(map(lambda x:'/'.join(x['names']),features))) )
parser.add_option('-n', parser.add_option('-n',
'--neighborhood', '--neighborhood',
dest = 'neighborhood', choices = neighborhoods.keys(), metavar = 'string', dest = 'neighborhood', choices = list(neighborhoods.keys()), metavar = 'string',
help = 'neighborhood type [neumann] {{{}}}'.format(', '.join(neighborhoods.keys()))) help = 'neighborhood type [neumann] {{{}}}'.format(', '.join(neighborhoods.keys())))
parser.add_option('-s', parser.add_option('-s',
'--scale', '--scale',
@ -151,10 +151,8 @@ for name in filenames:
remarks = [] remarks = []
column = {} column = {}
coordDim = table.label_dimension(options.pos) if not 3 >= table.label_dimension(options.pos) >= 1:
if not 3 >= coordDim >= 1:
errors.append('coordinates "{}" need to have one, two, or three dimensions.'.format(options.pos)) errors.append('coordinates "{}" need to have one, two, or three dimensions.'.format(options.pos))
else: coordCol = table.label_index(options.pos)
if table.label_dimension(options.id) != 1: errors.append('grain identifier {} not found.'.format(options.id)) if table.label_dimension(options.id) != 1: errors.append('grain identifier {} not found.'.format(options.id))
else: idCol = table.label_index(options.id) else: idCol = table.label_index(options.id)
@ -178,11 +176,7 @@ for name in filenames:
table.data_readArray() table.data_readArray()
coords = [np.unique(table.data[:,coordCol+i]) for i in range(coordDim)] grid,size = damask.util.coordGridAndSize(table.data[:,table.label_indexrange(options.pos)])
mincorner = np.array(map(min,coords))
maxcorner = np.array(map(max,coords))
grid = np.array(map(len,coords)+[1]*(3-len(coords)),'i')
N = grid.prod() N = grid.prod()
if N != len(table.data): errors.append('data count {} does not match grid {}.'.format(N,'x'.join(map(str,grid)))) if N != len(table.data): errors.append('data count {} does not match grid {}.'.format(N,'x'.join(map(str,grid))))

View File

@ -83,7 +83,7 @@ for name in filenames:
if table.label_dimension(options.pos) != 3: errors.append('coordinates {} are not a vector.'.format(options.pos)) if table.label_dimension(options.pos) != 3: errors.append('coordinates {} are not a vector.'.format(options.pos))
else: colCoord = table.label_index(options.pos) else: colCoord = table.label_index(options.pos)
for type, data in items.iteritems(): for type, data in items.items():
for what in (data['labels'] if data['labels'] is not None else []): for what in (data['labels'] if data['labels'] is not None else []):
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}.'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}.'.format(what,type))
@ -100,7 +100,7 @@ for name in filenames:
# ------------------------------------------ assemble header -------------------------------------- # ------------------------------------------ assemble header --------------------------------------
table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:])) table.info_append(scriptID + '\t' + ' '.join(sys.argv[1:]))
for type, data in items.iteritems(): for type, data in items.items():
for label in data['active']: for label in data['active']:
table.labels_append(['Gauss{}({})'.format(options.sigma,label)]) # extend ASCII header with new labels table.labels_append(['Gauss{}({})'.format(options.sigma,label)]) # extend ASCII header with new labels
table.head_write() table.head_write()
@ -114,7 +114,7 @@ for name in filenames:
# ------------------------------------------ process value field ----------------------------------- # ------------------------------------------ process value field -----------------------------------
stack = [table.data] stack = [table.data]
for type, data in items.iteritems(): for type, data in items.items():
for i,label in enumerate(data['active']): for i,label in enumerate(data['active']):
stack.append(ndimage.filters.gaussian_filter(table.data[:,data['column'][i]], stack.append(ndimage.filters.gaussian_filter(table.data[:,data['column'][i]],
options.sigma,options.order, options.sigma,options.order,

View File

@ -116,18 +116,18 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
if inputtype == 'eulers': if inputtype == 'eulers':
o = damask.Orientation(Eulers = np.array(map(float,table.data[column:column+3]))*toRadians, o = damask.Orientation(Eulers = np.array(list(map(float,table.data[column:column+3])))*toRadians,
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'matrix': elif inputtype == 'matrix':
o = damask.Orientation(matrix = np.array(map(float,table.data[column:column+9])).reshape(3,3).transpose(), o = damask.Orientation(matrix = np.array(list(map(float,table.data[column:column+9]))).reshape(3,3).transpose(),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'frame': elif inputtype == 'frame':
o = damask.Orientation(matrix = np.array(map(float,table.data[column[0]:column[0]+3] + \ o = damask.Orientation(matrix = np.array(list(map(float,table.data[column[0]:column[0]+3] + \
table.data[column[1]:column[1]+3] + \ table.data[column[1]:column[1]+3] + \
table.data[column[2]:column[2]+3])).reshape(3,3), table.data[column[2]:column[2]+3]))).reshape(3,3),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'quaternion': elif inputtype == 'quaternion':
o = damask.Orientation(quaternion = np.array(map(float,table.data[column:column+4])), o = damask.Orientation(quaternion = np.array(list(map(float,table.data[column:column+4]))),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
table.data_append(o.IPFcolor(pole)) table.data_append(o.IPFcolor(pole))

View File

@ -70,7 +70,7 @@ for name in filenames:
errors = [] errors = []
remarks = [] remarks = []
for type, data in items.iteritems(): for type, data in items.items():
for what in data['labels']: for what in data['labels']:
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type))
@ -94,7 +94,7 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
for type, data in items.iteritems(): for type, data in items.items():
for column in data['column']: for column in data['column']:
table.data_append(Mises(type, table.data_append(Mises(type,
np.array(table.data[column:column+data['dim']],'d').reshape(data['shape']))) np.array(table.data[column:column+data['dim']],'d').reshape(data['shape'])))

View File

@ -80,7 +80,7 @@ parser.set_defaults(output = [],
(options, filenames) = parser.parse_args() (options, filenames) = parser.parse_args()
options.output = map(lambda x: x.lower(), options.output) options.output = list(map(lambda x: x.lower(), options.output))
if options.output == [] or (not set(options.output).issubset(set(outputChoices))): if options.output == [] or (not set(options.output).issubset(set(outputChoices))):
parser.error('output must be chosen from {}.'.format(', '.join(outputChoices))) parser.error('output must be chosen from {}.'.format(', '.join(outputChoices)))
@ -147,21 +147,21 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
if inputtype == 'eulers': if inputtype == 'eulers':
o = damask.Orientation(Eulers = np.array(map(float,table.data[column:column+3]))*toRadians, o = damask.Orientation(Eulers = np.array(list(map(float,table.data[column:column+3])))*toRadians,
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'rodrigues': elif inputtype == 'rodrigues':
o = damask.Orientation(Rodrigues= np.array(map(float,table.data[column:column+3])), o = damask.Orientation(Rodrigues= np.array(list(map(float,table.data[column:column+3]))),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'matrix': elif inputtype == 'matrix':
o = damask.Orientation(matrix = np.array(map(float,table.data[column:column+9])).reshape(3,3).transpose(), o = damask.Orientation(matrix = np.array(list(map(float,table.data[column:column+9]))).reshape(3,3).transpose(),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'frame': elif inputtype == 'frame':
o = damask.Orientation(matrix = np.array(map(float,table.data[column[0]:column[0]+3] + \ o = damask.Orientation(matrix = np.array(list(map(float,table.data[column[0]:column[0]+3] + \
table.data[column[1]:column[1]+3] + \ table.data[column[1]:column[1]+3] + \
table.data[column[2]:column[2]+3])).reshape(3,3), table.data[column[2]:column[2]+3]))).reshape(3,3),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
elif inputtype == 'quaternion': elif inputtype == 'quaternion':
o = damask.Orientation(quaternion = np.array(map(float,table.data[column:column+4])), o = damask.Orientation(quaternion = np.array(list(map(float,table.data[column:column+4]))),
symmetry = options.symmetry).reduced() symmetry = options.symmetry).reduced()
o.quaternion = r*o.quaternion*R # apply additional lab and crystal frame rotations o.quaternion = r*o.quaternion*R # apply additional lab and crystal frame rotations

View File

@ -75,8 +75,8 @@ for name in filenames:
# ------------------------------------------ process data ------------------------------------------ # ------------------------------------------ process data ------------------------------------------
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
F = np.array(map(float,table.data[column[options.defgrad]:column[options.defgrad]+9]),'d').reshape(3,3) F = np.array(list(map(float,table.data[column[options.defgrad]:column[options.defgrad]+9])),'d').reshape(3,3)
P = np.array(map(float,table.data[column[options.stress ]:column[options.stress ]+9]),'d').reshape(3,3) P = np.array(list(map(float,table.data[column[options.stress ]:column[options.stress ]+9])),'d').reshape(3,3)
table.data_append(list(np.dot(np.linalg.inv(F),P).reshape(9))) # [S] =[P].[F-1] table.data_append(list(np.dot(np.linalg.inv(F),P).reshape(9))) # [S] =[P].[F-1]
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line

View File

@ -120,15 +120,15 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
if inputtype == 'eulers': if inputtype == 'eulers':
o = damask.Orientation(Eulers = np.array(map(float,table.data[column:column+3]))*toRadians) o = damask.Orientation(Eulers = np.array(list(map(float,table.data[column:column+3])))*toRadians)
elif inputtype == 'matrix': elif inputtype == 'matrix':
o = damask.Orientation(matrix = np.array(map(float,table.data[column:column+9])).reshape(3,3).transpose()) o = damask.Orientation(matrix = np.array(list(map(float,table.data[column:column+9]))).reshape(3,3).transpose())
elif inputtype == 'frame': elif inputtype == 'frame':
o = damask.Orientation(matrix = np.array(map(float,table.data[column[0]:column[0]+3] + \ o = damask.Orientation(matrix = np.array(list(map(float,table.data[column[0]:column[0]+3] + \
table.data[column[1]:column[1]+3] + \ table.data[column[1]:column[1]+3] + \
table.data[column[2]:column[2]+3])).reshape(3,3)) table.data[column[2]:column[2]+3]))).reshape(3,3))
elif inputtype == 'quaternion': elif inputtype == 'quaternion':
o = damask.Orientation(quaternion = np.array(map(float,table.data[column:column+4]))) o = damask.Orientation(quaternion = np.array(list(map(float,table.data[column:column+4]))))
rotatedPole = o.quaternion*pole # rotate pole according to crystal orientation rotatedPole = o.quaternion*pole # rotate pole according to crystal orientation
(x,y) = rotatedPole[0:2]/(1.+abs(pole[2])) # stereographic projection (x,y) = rotatedPole[0:2]/(1.+abs(pole[2])) # stereographic projection

View File

@ -252,15 +252,15 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
if inputtype == 'eulers': if inputtype == 'eulers':
o = damask.Orientation(Eulers = np.array(map(float,table.data[column:column+3]))*toRadians,) o = damask.Orientation(Eulers = np.array(list(map(float,table.data[column:column+3])))*toRadians,)
elif inputtype == 'matrix': elif inputtype == 'matrix':
o = damask.Orientation(matrix = np.array(map(float,table.data[column:column+9])).reshape(3,3).transpose(),) o = damask.Orientation(matrix = np.array(list(map(float,table.data[column:column+9]))).reshape(3,3).transpose(),)
elif inputtype == 'frame': elif inputtype == 'frame':
o = damask.Orientation(matrix = np.array(map(float,table.data[column[0]:column[0]+3] + \ o = damask.Orientation(matrix = np.array(list(map(float,table.data[column[0]:column[0]+3] + \
table.data[column[1]:column[1]+3] + \ table.data[column[1]:column[1]+3] + \
table.data[column[2]:column[2]+3])).reshape(3,3),) table.data[column[2]:column[2]+3]))).reshape(3,3),)
elif inputtype == 'quaternion': elif inputtype == 'quaternion':
o = damask.Orientation(quaternion = np.array(map(float,table.data[column:column+4])),) o = damask.Orientation(quaternion = np.array(list(map(float,table.data[column:column+4]))),)
rotForce = o.quaternion.conjugated() * force rotForce = o.quaternion.conjugated() * force
rotNormal = o.quaternion.conjugated() * normal rotNormal = o.quaternion.conjugated() * normal

View File

@ -58,7 +58,7 @@ for name in filenames:
errors = [] errors = []
remarks = [] remarks = []
for type, data in items.iteritems(): for type, data in items.items():
for what in data['labels']: for what in data['labels']:
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type))
@ -84,9 +84,9 @@ for name in filenames:
outputAlive = True outputAlive = True
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
for type, data in items.iteritems(): for type, data in items.items():
for column in data['column']: for column in data['column']:
(u,v) = np.linalg.eigh(np.array(map(float,table.data[column:column+data['dim']])).reshape(data['shape'])) (u,v) = np.linalg.eigh(np.array(list(map(float,table.data[column:column+data['dim']]))).reshape(data['shape']))
if options.rh and np.dot(np.cross(v[:,0], v[:,1]), v[:,2]) < 0.0 : v[:, 2] *= -1.0 # ensure right-handed eigenvector basis if options.rh and np.dot(np.cross(v[:,0], v[:,1]), v[:,2]) < 0.0 : v[:, 2] *= -1.0 # ensure right-handed eigenvector basis
table.data_append(list(u)) # vector of max,mid,min eigval table.data_append(list(u)) # vector of max,mid,min eigval
table.data_append(list(v.transpose().reshape(data['dim']))) # 3x3=9 combo vector of max,mid,min eigvec coordinates table.data_append(list(v.transpose().reshape(data['dim']))) # 3x3=9 combo vector of max,mid,min eigvec coordinates

View File

@ -101,7 +101,7 @@ for name in filenames:
errors = [] errors = []
remarks = [] remarks = []
for type, data in items.iteritems(): for type, data in items.items():
for what in data['labels']: for what in data['labels']:
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}...'.format(what,type))
@ -132,7 +132,7 @@ for name in filenames:
while outputAlive and table.data_read(): # read next data line of ASCII table while outputAlive and table.data_read(): # read next data line of ASCII table
for column in items['tensor']['column']: # loop over all requested defgrads for column in items['tensor']['column']: # loop over all requested defgrads
F = np.array(map(float,table.data[column:column+items['tensor']['dim']]),'d').reshape(items['tensor']['shape']) F = np.array(list(map(float,table.data[column:column+items['tensor']['dim']])),'d').reshape(items['tensor']['shape'])
(U,S,Vh) = np.linalg.svd(F) # singular value decomposition (U,S,Vh) = np.linalg.svd(F) # singular value decomposition
R = np.dot(U,Vh) # rotation of polar decomposition R = np.dot(U,Vh) # rotation of polar decomposition
stretch['U'] = np.dot(np.linalg.inv(R),F) # F = RU stretch['U'] = np.dot(np.linalg.inv(R),F) # F = RU

View File

@ -76,7 +76,6 @@ for name in filenames:
remarks = [] remarks = []
if table.label_dimension(options.pos) != 3: errors.append('coordinates {} are not a vector.'.format(options.pos)) if table.label_dimension(options.pos) != 3: errors.append('coordinates {} are not a vector.'.format(options.pos))
else: colCoord = table.label_index(options.pos)
if remarks != []: damask.util.croak(remarks) if remarks != []: damask.util.croak(remarks)
if errors != []: if errors != []:
@ -94,14 +93,7 @@ for name in filenames:
table.data_readArray() table.data_readArray()
if (any(options.grid) == 0 or any(options.size) == 0.0): if (any(options.grid) == 0 or any(options.size) == 0.0):
coords = [np.unique(table.data[:,colCoord+i]) for i in range(3)] grid,size = damask.util.coordGridAndSize(table.data[:,table.label_indexrange(options.pos)])
mincorner = np.array(map(min,coords))
maxcorner = np.array(map(max,coords))
grid = np.array(map(len,coords),'i')
size = grid/np.maximum(np.ones(3,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1)
size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 set to smallest among other spacings
delta = size/np.maximum(np.ones(3,'d'), grid)
origin = mincorner - 0.5*delta # shift from cell center to corner
else: else:
grid = np.array(options.grid,'i') grid = np.array(options.grid,'i')
@ -129,16 +121,15 @@ for name in filenames:
#--- generate grid -------------------------------------------------------------------------------- #--- generate grid --------------------------------------------------------------------------------
if colCoord: x = (0.5 + shift[0] + np.arange(packedGrid[0],dtype=float))/packedGrid[0]*size[0] + origin[0]
x = (0.5 + shift[0] + np.arange(packedGrid[0],dtype=float))/packedGrid[0]*size[0] + origin[0] y = (0.5 + shift[1] + np.arange(packedGrid[1],dtype=float))/packedGrid[1]*size[1] + origin[1]
y = (0.5 + shift[1] + np.arange(packedGrid[1],dtype=float))/packedGrid[1]*size[1] + origin[1] z = (0.5 + shift[2] + np.arange(packedGrid[2],dtype=float))/packedGrid[2]*size[2] + origin[2]
z = (0.5 + shift[2] + np.arange(packedGrid[2],dtype=float))/packedGrid[2]*size[2] + origin[2]
xx = np.tile( x, packedGrid[1]* packedGrid[2]) xx = np.tile( x, packedGrid[1]* packedGrid[2])
yy = np.tile(np.repeat(y,packedGrid[0] ),packedGrid[2]) yy = np.tile(np.repeat(y,packedGrid[0] ),packedGrid[2])
zz = np.repeat(z,packedGrid[0]*packedGrid[1]) zz = np.repeat(z,packedGrid[0]*packedGrid[1])
table.data[:,colCoord:colCoord+3] = np.squeeze(np.dstack((xx,yy,zz))) table.data[:,table.label_indexragen(options.pos)] = np.squeeze(np.dstack((xx,yy,zz)))
# ------------------------------------------ output result ----------------------------------------- # ------------------------------------------ output result -----------------------------------------

View File

@ -64,7 +64,6 @@ for name in filenames:
remarks = [] remarks = []
if table.label_dimension(options.pos) != 3: errors.append('coordinates "{}" are not a vector.'.format(options.pos)) if table.label_dimension(options.pos) != 3: errors.append('coordinates "{}" are not a vector.'.format(options.pos))
else: colCoord = table.label_index(options.pos)
colElem = table.label_index('elem') colElem = table.label_index('elem')
@ -79,12 +78,7 @@ for name in filenames:
table.data_readArray(options.pos) table.data_readArray(options.pos)
table.data_rewind() table.data_rewind()
coords = [np.unique(table.data[:,i]) for i in range(3)] grid,size = damask.util.coordGridAndSize(table.data)
mincorner = np.array(map(min,coords))
maxcorner = np.array(map(max,coords))
grid = np.array(map(len,coords),'i')
size = grid/np.maximum(np.ones(3,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1)
size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 set to smallest among other spacings
packing = np.array(options.packing,'i') packing = np.array(options.packing,'i')
outSize = grid*packing outSize = grid*packing
@ -113,7 +107,7 @@ for name in filenames:
for c in range(outSize[2]): for c in range(outSize[2]):
for b in range(outSize[1]): for b in range(outSize[1]):
for a in range(outSize[0]): for a in range(outSize[0]):
data[a,b,c,colCoord:colCoord+3] = [a+0.5,b+0.5,c+0.5]*elementSize data[a,b,c,table.label_indexrange(options.pos)] = [a+0.5,b+0.5,c+0.5]*elementSize
if colElem != -1: data[a,b,c,colElem] = elem if colElem != -1: data[a,b,c,colElem] = elem
table.data = data[a,b,c,:].tolist() table.data = data[a,b,c,:].tolist()
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line

View File

@ -73,7 +73,7 @@ for name in filenames:
remarks = [] remarks = []
column = {} column = {}
for type, data in items.iteritems(): for type, data in items.items():
for what in data['labels']: for what in data['labels']:
dim = table.label_dimension(what) dim = table.label_dimension(what)
if dim != data['dim']: remarks.append('column {} is not a {}.'.format(what,type)) if dim != data['dim']: remarks.append('column {} is not a {}.'.format(what,type))
@ -100,13 +100,13 @@ for name in filenames:
for column in items[datatype]['column']: # loop over all requested labels for column in items[datatype]['column']: # loop over all requested labels
table.data[column:column+items[datatype]['dim']] = \ table.data[column:column+items[datatype]['dim']] = \
q * np.array(map(float,table.data[column:column+items[datatype]['dim']])) q * np.array(list(map(float,table.data[column:column+items[datatype]['dim']])))
datatype = 'tensor' datatype = 'tensor'
for column in items[datatype]['column']: # loop over all requested labels for column in items[datatype]['column']: # loop over all requested labels
table.data[column:column+items[datatype]['dim']] = \ table.data[column:column+items[datatype]['dim']] = \
np.dot(R,np.dot(np.array(map(float,table.data[column:column+items[datatype]['dim']])).\ np.dot(R,np.dot(np.array(list(map(float,table.data[column:column+items[datatype]['dim']]))).\
reshape(items[datatype]['shape']),R.transpose())).reshape(items[datatype]['dim']) reshape(items[datatype]['shape']),R.transpose())).reshape(items[datatype]['dim'])
outputAlive = table.data_write() # output processed line outputAlive = table.data_write() # output processed line

View File

@ -421,8 +421,6 @@ for filename in filenames:
meshActor.GetProperty().SetOpacity(0.2) meshActor.GetProperty().SetOpacity(0.2)
meshActor.GetProperty().SetColor(1.0,1.0,0) meshActor.GetProperty().SetColor(1.0,1.0,0)
meshActor.GetProperty().BackfaceCullingOn() meshActor.GetProperty().BackfaceCullingOn()
# meshActor.GetProperty().SetEdgeColor(1,1,0.5)
# meshActor.GetProperty().EdgeVisibilityOn()
boxpoints = vtk.vtkPoints() boxpoints = vtk.vtkPoints()
for n in range(8): for n in range(8):

View File

@ -82,7 +82,7 @@ for name in filenames:
[coords[i][j-1] + coords[i][j] for j in range(1,len(coords[i]))] + \ [coords[i][j-1] + coords[i][j] for j in range(1,len(coords[i]))] + \
[3.0 * coords[i][-1] - coords[i][-1 - int(len(coords[i]) > 1)]]) for i in range(3)] [3.0 * coords[i][-1] - coords[i][-1 - int(len(coords[i]) > 1)]]) for i in range(3)]
grid = np.array(map(len,coords),'i') grid = np.array(list(map(len,coords)),'i')
N = grid.prod() if options.mode == 'point' else (grid-1).prod() N = grid.prod() if options.mode == 'point' else (grid-1).prod()
if N != len(table.data): if N != len(table.data):

View File

@ -55,9 +55,9 @@ parser.set_defaults(canal = 25e-6,
(options,filename) = parser.parse_args() (options,filename) = parser.parse_args()
if np.any(options.grid < 2): if np.any(np.array(options.grid) < 2):
parser('invalid grid a b c.') parser('invalid grid a b c.')
if np.any(options.size <= 0.0): if np.any(np.array(options.size) <= 0.0):
parser('invalid size x y z.') parser('invalid size x y z.')
# --- open input files ---------------------------------------------------------------------------- # --- open input files ----------------------------------------------------------------------------
@ -114,12 +114,8 @@ for y in range(info['grid'][1]):
info['microstructures'] += 1 info['microstructures'] += 1
#--- report --------------------------------------------------------------------------------------- #--- report ---------------------------------------------------------------------------------------
damask.util.croak(['grid a b c: %s'%(' x '.join(map(str,info['grid']))), damask.util.report_geom(info,['grid','size','origin','homogenization','microstructures'])
'size x y z: %s'%(' x '.join(map(str,info['size']))),
'origin x y z: %s'%(' : '.join(map(str,info['origin']))),
'homogenization: %i'%info['homogenization'],
'microstructures: %i'%info['microstructures']])
# -------------------------------------- switch according to task ----------------------------------
formatwidth = 1+int(math.floor(math.log10(info['microstructures']-1))) formatwidth = 1+int(math.floor(math.log10(info['microstructures']-1)))
header = [scriptID + ' ' + ' '.join(sys.argv[1:])] header = [scriptID + ' ' + ' '.join(sys.argv[1:])]
header.append('<microstructure>') header.append('<microstructure>')

View File

@ -1,7 +1,7 @@
#!/usr/bin/env python2.7 #!/usr/bin/env python2.7
# -*- coding: UTF-8 no BOM -*- # -*- coding: UTF-8 no BOM -*-
import os,sys,math,types,time import os,sys,math,time
import scipy.spatial, numpy as np import scipy.spatial, numpy as np
from optparse import OptionParser from optparse import OptionParser
import damask import damask
@ -152,7 +152,7 @@ for name in filenames:
continue continue
table.data_readArray([options.pos] \ table.data_readArray([options.pos] \
+ ([label] if isinstance(label, types.StringTypes) else label) \ + (label if isinstance(label, list) else [label]) \
+ ([options.phase] if options.phase else [])) + ([options.phase] if options.phase else []))
if coordDim == 2: if coordDim == 2:
@ -165,9 +165,9 @@ for name in filenames:
# --------------- figure out size and grid --------------------------------------------------------- # --------------- figure out size and grid ---------------------------------------------------------
coords = [np.unique(table.data[:,i]) for i in range(3)] coords = [np.unique(table.data[:,i]) for i in range(3)]
mincorner = np.array(map(min,coords)) mincorner = np.array(list(map(min,coords)))
maxcorner = np.array(map(max,coords)) maxcorner = np.array(list(map(max,coords)))
grid = np.array(map(len,coords),'i') grid = np.array(list(map(len,coords)),'i')
size = grid/np.maximum(np.ones(3,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1) size = grid/np.maximum(np.ones(3,'d'), grid-1.0) * (maxcorner-mincorner) # size from edge to edge = dim * n/(n-1)
size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 set to smallest among other spacings size = np.where(grid > 1, size, min(size[grid > 1]/grid[grid > 1])) # spacing for grid==1 set to smallest among other spacings
delta = size/np.maximum(np.ones(3,'d'), grid) delta = size/np.maximum(np.ones(3,'d'), grid)

View File

@ -15,8 +15,7 @@ scriptID = ' '.join([scriptName,damask.version])
def meshgrid2(*arrs): def meshgrid2(*arrs):
"""Code inspired by http://stackoverflow.com/questions/1827489/numpy-meshgrid-in-3d""" """Code inspired by http://stackoverflow.com/questions/1827489/numpy-meshgrid-in-3d"""
arrs = tuple(reversed(arrs)) arrs = tuple(reversed(arrs))
arrs = tuple(arrs) lens = np.array(list(map(len, arrs)))
lens = np.array(map(len, arrs))
dim = len(arrs) dim = len(arrs)
ans = [] ans = []
for i, arr in enumerate(arrs): for i, arr in enumerate(arrs):

View File

@ -92,7 +92,7 @@ for name in filenames:
} }
substituted = np.copy(microstructure) substituted = np.copy(microstructure)
for k, v in sub.iteritems(): substituted[microstructure==k] = v # substitute microstructure indices for k, v in sub.items(): substituted[microstructure==k] = v # substitute microstructure indices
substituted += options.microstructure # shift microstructure indices substituted += options.microstructure # shift microstructure indices

View File

@ -270,7 +270,7 @@ for name in filenames:
ODF['limit'] = np.radians(limits[1,:]) # right hand limits in radians ODF['limit'] = np.radians(limits[1,:]) # right hand limits in radians
ODF['center'] = 0.0 if all(limits[0,:]<1e-8) else 0.5 # vertex or cell centered ODF['center'] = 0.0 if all(limits[0,:]<1e-8) else 0.5 # vertex or cell centered
ODF['interval'] = np.array(map(len,[np.unique(table.data[:,i]) for i in range(3)]),'i') # steps are number of distict values ODF['interval'] = np.array(list(map(len,[np.unique(table.data[:,i]) for i in range(3)])),'i') # steps are number of distict values
ODF['nBins'] = ODF['interval'].prod() ODF['nBins'] = ODF['interval'].prod()
ODF['delta'] = np.radians(np.array(limits[1,0:3]-limits[0,0:3])/(ODF['interval']-1)) # step size ODF['delta'] = np.radians(np.array(limits[1,0:3]-limits[0,0:3])/(ODF['interval']-1)) # step size

View File

@ -77,7 +77,14 @@ def mesh(r,d):
"%f %f %f"%(-d[0],d[1],d[2]), "%f %f %f"%(-d[0],d[1],d[2]),
"%f %f %f"%(-d[0],d[1],0.0), "%f %f %f"%(-d[0],d[1],0.0),
"*add_elements", "*add_elements",
range(1,9), "1",
"2",
"3",
"4",
"5",
"6",
"7",
"8",
"*sub_divisions", "*sub_divisions",
"%i %i %i"%(r[2],r[1],r[0]), "%i %i %i"%(r[2],r[1],r[0]),
"*subdivide_elements", "*subdivide_elements",
@ -201,7 +208,7 @@ if options.port:
except: except:
parser.error('no valid Mentat release found.') parser.error('no valid Mentat release found.')
# --- loop over input files ------------------------------------------------------------------------- # --- loop over input files ------------------------------------------------------------------------
if filenames == []: filenames = [None] if filenames == []: filenames = [None]
@ -236,7 +243,7 @@ for name in filenames:
# --- read data ------------------------------------------------------------------------------------ # --- read data ------------------------------------------------------------------------------------
microstructure = table.microstructure_read(info['grid']).reshape(info['grid'].prod(),order='F') # read microstructure microstructure = table.microstructure_read(info['grid']).reshape(info['grid'].prod(),order='F') # read microstructure
cmds = [\ cmds = [\
init(), init(),

View File

@ -344,7 +344,7 @@ def rcbParser(content,M,size,tolerance,idcolumn,segmentcolumn):
else: else:
myNeighbors[grainNeighbors[leg][side]] = 1 myNeighbors[grainNeighbors[leg][side]] = 1
if myNeighbors: # do I have any neighbors (i.e., non-bounding box segment) if myNeighbors: # do I have any neighbors (i.e., non-bounding box segment)
candidateGrains = sorted(myNeighbors.iteritems(), key=lambda p: (p[1],p[0]), reverse=True) # sort grain counting candidateGrains = sorted(myNeighbors.items(), key=lambda p: (p[1],p[0]), reverse=True) # sort grain counting
# most frequent one not yet seen? # most frequent one not yet seen?
rcData['grainMapping'].append(candidateGrains[0 if candidateGrains[0][0] not in rcData['grainMapping'] else 1][0]) # must be me then rcData['grainMapping'].append(candidateGrains[0 if candidateGrains[0][0] not in rcData['grainMapping'] else 1][0]) # must be me then
# special case of bi-crystal situation... # special case of bi-crystal situation...

View File

@ -667,6 +667,7 @@ function getStrings(this,key,defaultVal,requiredShape,raw)
endif endif
else notAllocated else notAllocated
if (whole) then if (whole) then
str = item%string%val(item%string%pos(4):)
getStrings = [getStrings,str] getStrings = [getStrings,str]
else else
do i=2_pInt,item%string%pos(1) do i=2_pInt,item%string%pos(1)

View File

@ -160,7 +160,7 @@ module math
math_rotate_forward33, & math_rotate_forward33, &
math_rotate_backward33, & math_rotate_backward33, &
math_rotate_forward3333, & math_rotate_forward3333, &
math_limit math_clip
private :: & private :: &
math_check, & math_check, &
halton halton
@ -1363,16 +1363,16 @@ pure function math_RtoEuler(R)
sqhk =sqrt(R(1,3)*R(1,3)+R(2,3)*R(2,3)) sqhk =sqrt(R(1,3)*R(1,3)+R(2,3)*R(2,3))
! calculate PHI ! calculate PHI
math_RtoEuler(2) = acos(math_limit(R(3,3)/sqhkl,-1.0_pReal, 1.0_pReal)) math_RtoEuler(2) = acos(math_clip(R(3,3)/sqhkl,-1.0_pReal, 1.0_pReal))
if((math_RtoEuler(2) < 1.0e-8_pReal) .or. (pi-math_RtoEuler(2) < 1.0e-8_pReal)) then if((math_RtoEuler(2) < 1.0e-8_pReal) .or. (pi-math_RtoEuler(2) < 1.0e-8_pReal)) then
math_RtoEuler(3) = 0.0_pReal math_RtoEuler(3) = 0.0_pReal
math_RtoEuler(1) = acos(math_limit(R(1,1)/squvw, -1.0_pReal, 1.0_pReal)) math_RtoEuler(1) = acos(math_clip(R(1,1)/squvw, -1.0_pReal, 1.0_pReal))
if(R(2,1) > 0.0_pReal) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1) if(R(2,1) > 0.0_pReal) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1)
else else
math_RtoEuler(3) = acos(math_limit(R(2,3)/sqhk, -1.0_pReal, 1.0_pReal)) math_RtoEuler(3) = acos(math_clip(R(2,3)/sqhk, -1.0_pReal, 1.0_pReal))
if(R(1,3) < 0.0) math_RtoEuler(3) = 2.0_pReal*pi-math_RtoEuler(3) if(R(1,3) < 0.0) math_RtoEuler(3) = 2.0_pReal*pi-math_RtoEuler(3)
math_RtoEuler(1) = acos(math_limit(-R(3,2)/sin(math_RtoEuler(2)), -1.0_pReal, 1.0_pReal)) math_RtoEuler(1) = acos(math_clip(-R(3,2)/sin(math_RtoEuler(2)), -1.0_pReal, 1.0_pReal))
if(R(3,1) < 0.0) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1) if(R(3,1) < 0.0) math_RtoEuler(1) = 2.0_pReal*pi-math_RtoEuler(1)
end if end if
@ -1654,7 +1654,7 @@ pure function math_qToEuler(qPassive)
math_qToEuler(2) = acos(1.0_pReal-2.0_pReal*(q(2)**2+q(3)**2)) math_qToEuler(2) = acos(1.0_pReal-2.0_pReal*(q(2)**2+q(3)**2))
if (abs(math_qToEuler(2)) < 1.0e-6_pReal) then if (abs(math_qToEuler(2)) < 1.0e-6_pReal) then
math_qToEuler(1) = sign(2.0_pReal*acos(math_limit(q(1),-1.0_pReal, 1.0_pReal)),q(4)) math_qToEuler(1) = sign(2.0_pReal*acos(math_clip(q(1),-1.0_pReal, 1.0_pReal)),q(4))
math_qToEuler(3) = 0.0_pReal math_qToEuler(3) = 0.0_pReal
else else
math_qToEuler(1) = atan2(+q(1)*q(3)+q(2)*q(4), q(1)*q(2)-q(3)*q(4)) math_qToEuler(1) = atan2(+q(1)*q(3)+q(2)*q(4), q(1)*q(2)-q(3)*q(4))
@ -1681,7 +1681,7 @@ pure function math_qToAxisAngle(Q)
real(pReal) :: halfAngle, sinHalfAngle real(pReal) :: halfAngle, sinHalfAngle
real(pReal), dimension(4) :: math_qToAxisAngle real(pReal), dimension(4) :: math_qToAxisAngle
halfAngle = acos(math_limit(Q(1),-1.0_pReal,1.0_pReal)) halfAngle = acos(math_clip(Q(1),-1.0_pReal,1.0_pReal))
sinHalfAngle = sin(halfAngle) sinHalfAngle = sin(halfAngle)
smallRotation: if (sinHalfAngle <= 1.0e-4_pReal) then smallRotation: if (sinHalfAngle <= 1.0e-4_pReal) then
@ -1741,7 +1741,7 @@ real(pReal) pure function math_EulerMisorientation(EulerA,EulerB)
cosTheta = (math_trace33(math_mul33x33(math_EulerToR(EulerB), & cosTheta = (math_trace33(math_mul33x33(math_EulerToR(EulerB), &
transpose(math_EulerToR(EulerA)))) - 1.0_pReal) * 0.5_pReal transpose(math_EulerToR(EulerA)))) - 1.0_pReal) * 0.5_pReal
math_EulerMisorientation = acos(math_limit(cosTheta,-1.0_pReal,1.0_pReal)) math_EulerMisorientation = acos(math_clip(cosTheta,-1.0_pReal,1.0_pReal))
end function math_EulerMisorientation end function math_EulerMisorientation
@ -2052,7 +2052,7 @@ function math_eigenvectorBasisSym33(m)
EB(3,3,3)=1.0_pReal EB(3,3,3)=1.0_pReal
else threeSimilarEigenvalues else threeSimilarEigenvalues
rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal
phi=acos(math_limit(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal)) phi=acos(math_clip(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal))
values = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* & values = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* &
[cos(phi/3.0_pReal), & [cos(phi/3.0_pReal), &
cos((phi+2.0_pReal*PI)/3.0_pReal), & cos((phi+2.0_pReal*PI)/3.0_pReal), &
@ -2117,7 +2117,7 @@ function math_eigenvectorBasisSym33_log(m)
EB(3,3,3)=1.0_pReal EB(3,3,3)=1.0_pReal
else threeSimilarEigenvalues else threeSimilarEigenvalues
rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal
phi=acos(math_limit(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal)) phi=acos(math_clip(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal))
values = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* & values = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* &
[cos(phi/3.0_pReal), & [cos(phi/3.0_pReal), &
cos((phi+2.0_pReal*PI)/3.0_pReal), & cos((phi+2.0_pReal*PI)/3.0_pReal), &
@ -2229,7 +2229,7 @@ function math_eigenvaluesSym33(m)
math_eigenvaluesSym33 = math_eigenvaluesSym(m) math_eigenvaluesSym33 = math_eigenvaluesSym(m)
else else
rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal rho=sqrt(-3.0_pReal*P**3.0_pReal)/9.0_pReal
phi=acos(math_limit(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal)) phi=acos(math_clip(-Q/rho*0.5_pReal,-1.0_pReal,1.0_pReal))
math_eigenvaluesSym33 = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* & math_eigenvaluesSym33 = 2.0_pReal*rho**(1.0_pReal/3.0_pReal)* &
[cos(phi/3.0_pReal), & [cos(phi/3.0_pReal), &
cos((phi+2.0_pReal*PI)/3.0_pReal), & cos((phi+2.0_pReal*PI)/3.0_pReal), &
@ -2614,7 +2614,7 @@ end function math_rotate_forward3333
!> @brief limits a scalar value to a certain range (either one or two sided) !> @brief limits a scalar value to a certain range (either one or two sided)
! Will return NaN if left > right ! Will return NaN if left > right
!-------------------------------------------------------------------------------------------------- !--------------------------------------------------------------------------------------------------
real(pReal) pure function math_limit(a, left, right) real(pReal) pure function math_clip(a, left, right)
use, intrinsic :: & use, intrinsic :: &
IEEE_arithmetic IEEE_arithmetic
@ -2623,14 +2623,14 @@ real(pReal) pure function math_limit(a, left, right)
real(pReal), intent(in), optional :: left, right real(pReal), intent(in), optional :: left, right
math_limit = min ( & math_clip = min ( &
max (merge(left, -huge(a), present(left)), a), & max (merge(left, -huge(a), present(left)), a), &
merge(right, huge(a), present(right)) & merge(right, huge(a), present(right)) &
) )
if (present(left) .and. present(right)) & if (present(left) .and. present(right)) &
math_limit = merge (IEEE_value(1.0_pReal,IEEE_quiet_NaN),math_limit, left>right) math_clip = merge (IEEE_value(1.0_pReal,IEEE_quiet_NaN),math_clip, left>right)
end function math_limit end function math_clip
end module math end module math