Pre-Processing

The following classes and modules support the setup of DAMASK simulations.

Config

class damask.Config[source]

YAML-based configuration.

Attributes
is_complete

Check for completeness.

is_valid

Check for valid file layout.

Methods

clear()

copy()

Create deep copy.

delete(keys)

Remove configuration keys.

fromkeys(iterable[, value])

Create a new dictionary with keys from iterable and values set to value.

get(key[, default])

Return the value for key if key is in the dictionary, else default.

items()

keys()

load(fname)

Load from yaml file.

pop(k[,d])

If key is not found, d is returned if given, otherwise KeyError is raised

popitem(/)

Remove and return a (key, value) pair as a 2-tuple.

save(fname, **kwargs)

Save to yaml file.

setdefault(key[, default])

Insert key with a value of default if key is not in the dictionary.

update([E, ]**F)

If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]

values()

copy()

Create deep copy.

delete(keys)[source]

Remove configuration keys.

Parameters
keysiterable or scalar

Label of the key(s) to remove.

Returns
updateddamask.Config

Updated configuration.

classmethod load(fname)[source]

Load from yaml file.

Parameters
fnamefile, str, or pathlib.Path

Filename or file for writing.

Returns
loadeddamask.Config

Configuration from file.

save(fname, **kwargs)[source]

Save to yaml file.

Parameters
fnamefile, str, or pathlib.Path

Filename or file for writing.

**kwargsdict

Keyword arguments parsed to yaml.dump.

abstract property is_complete

Check for completeness.

abstract property is_valid

Check for valid file layout.

clear() None.  Remove all items from D.
fromkeys(iterable, value=None, /)

Create a new dictionary with keys from iterable and values set to value.

get(key, default=None, /)

Return the value for key if key is in the dictionary, else default.

items() a set-like object providing a view on D's items
keys() a set-like object providing a view on D's keys
pop(k[, d]) v, remove specified key and return the corresponding value.

If key is not found, d is returned if given, otherwise KeyError is raised

popitem(/)

Remove and return a (key, value) pair as a 2-tuple.

Pairs are returned in LIFO (last-in, first-out) order. Raises KeyError if the dict is empty.

setdefault(key, default=None, /)

Insert key with a value of default if key is not in the dictionary.

Return the value for key if key is in the dictionary, else default.

update([E, ]**F) None.  Update D from dict/iterable E and F.

If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]

values() an object providing a view on D's values

ConfigMaterial

class damask.ConfigMaterial(d=None)[source]

Material configuration.

Manipulate material configurations for storage in YAML format. A complete material configuration file has the entries ‘material’, ‘phase’, and ‘homogenization’. For use in DAMASK, it needs to be stored as ‘material.yaml’.

Attributes
is_complete

Check for completeness.

is_valid

Check for valid content.

Methods

clear()

copy()

Create deep copy.

delete(keys)

Remove configuration keys.

from_table(table, **kwargs)

Generate from an ASCII table.

fromkeys(iterable[, value])

Create a new dictionary with keys from iterable and values set to value.

get(key[, default])

Return the value for key if key is in the dictionary, else default.

items()

keys()

load([fname])

Load from yaml file.

load_DREAM3D(fname[, grain_data, cell_data, ...])

Load DREAM.3D (HDF5) file.

material_add(**kwargs)

Add material entries.

material_rename_homogenization(mapping[, ID])

Change homogenization name in material.

material_rename_phase(mapping[, ID, constituent])

Change phase name in material.

pop(k[,d])

If key is not found, d is returned if given, otherwise KeyError is raised

popitem(/)

Remove and return a (key, value) pair as a 2-tuple.

save([fname])

Save to yaml file.

setdefault(key[, default])

Insert key with a value of default if key is not in the dictionary.

update([E, ]**F)

If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]

values()

save(fname='material.yaml', **kwargs)[source]

Save to yaml file.

Parameters
fnamefile, str, or pathlib.Path, optional

Filename or file for writing. Defaults to ‘material.yaml’.

**kwargs

Keyword arguments parsed to yaml.dump.

classmethod load(fname='material.yaml')[source]

Load from yaml file.

Parameters
fnamefile, str, or pathlib.Path, optional

Filename or file to read from. Defaults to ‘material.yaml’.

Returns
loadeddamask.ConfigMaterial

Material configuration from file.

static load_DREAM3D(fname, grain_data=None, cell_data=None, cell_ensemble_data='CellEnsembleData', phases='Phases', Euler_angles='EulerAngles', phase_names='PhaseName', base_group=None)[source]

Load DREAM.3D (HDF5) file.

Data in DREAM.3D files can be stored per cell (‘CellData’) and/or per grain (‘Grain Data’). Per default, cell-wise data is assumed.

damask.Grid.load_DREAM3D allows to get the corresponding geometry for the grid solver.

Parameters
fnamestr

Filename of the DREAM.3D (HDF5) file.

grain_datastr

Name of the group (folder) containing grain-wise data. Defaults to None, in which case cell-wise data is used.

cell_datastr

Name of the group (folder) containing cell-wise data. Defaults to None in wich case it is automatically detected.

cell_ensemble_datastr

Name of the group (folder) containing data of cell ensembles. This group is used to inquire the name of the phases. Phases will get numeric IDs if this group is not found. Defaults to ‘CellEnsembleData’.

phasesstr

Name of the dataset containing the phase ID (cell-wise or grain-wise). Defaults to ‘Phases’.

Euler_anglesstr

Name of the dataset containing the crystallographic orientation as Euler angles in radians (cell-wise or grain-wise). Defaults to ‘EulerAngles’.

phase_namesstr

Name of the dataset containing the phase names. Phases will get numeric IDs if this dataset is not found. Defaults to ‘PhaseName’.

base_groupstr

Path to the group (folder) that contains geometry (_SIMPL_GEOMETRY), and grain- or cell-wise data. Defaults to None, in which case it is set as the path that contains _SIMPL_GEOMETRY/SPACING.

Returns
loadeddamask.ConfigMaterial

Material configuration from file.

Notes

Homogenization and phase entries are emtpy and need to be defined separately.

static from_table(table, **kwargs)[source]

Generate from an ASCII table.

Parameters
tabledamask.Table

Table that contains material information.

**kwargs

Keyword arguments where the key is the name and the value specifies the label of the data column in the table.

Returns
newdamask.ConfigMaterial

Material configuration from values in table.

Examples

>>> import damask
>>> import damask.ConfigMaterial as cm
>>> t = damask.Table.load('small.txt')
>>> t
    pos  pos  pos   qu   qu    qu    qu   phase    homog
0    0    0    0  0.19  0.8   0.24 -0.51  Aluminum SX
1    1    0    0  0.8   0.19  0.24 -0.51  Steel    SX
1    1    1    0  0.8   0.19  0.24 -0.51  Steel    SX
>>> cm.from_table(t,O='qu',phase='phase',homogenization='homog')
material:
  - constituents:
      - O: [0.19, 0.8, 0.24, -0.51]
        v: 1.0
        phase: Aluminum
    homogenization: SX
  - constituents:
      - O: [0.8, 0.19, 0.24, -0.51]
        v: 1.0
        phase: Steel
    homogenization: SX
homogenization: {}
phase: {}
property is_complete

Check for completeness.

Only the general file layout is considered. This check does not consider whether parameters for a particular phase/homogenization model are missing.

Returns
completebool

Whether the material.yaml definition is complete.

property is_valid

Check for valid content.

Only the generic file content is considered. This check does not consider whether parameters for a particular phase/homogenization mode are out of bounds.

Returns
validbool

Whether the material.yaml definition is valid.

material_rename_phase(mapping, ID=None, constituent=None)[source]

Change phase name in material.

Parameters
mapping: dictionary

Mapping from old name to new name

ID: list of ints, optional

Limit renaming to selected material IDs.

constituent: list of ints, optional

Limit renaming to selected constituents.

Returns
updateddamask.ConfigMaterial

Updated material configuration.

material_rename_homogenization(mapping, ID=None)[source]

Change homogenization name in material.

Parameters
mapping: dictionary

Mapping from old name to new name

ID: list of ints, optional

Limit renaming to selected homogenization IDs.

Returns
updateddamask.ConfigMaterial

Updated material configuration.

material_add(**kwargs)[source]

Add material entries.

Parameters
**kwargs

Key-value pairs.

Returns
updateddamask.ConfigMaterial

Updated material configuration.

Examples

Create a dual-phase steel microstructure for micromechanical simulations:

>>> import numpy as np
>>> import damask
>>> m = damask.ConfigMaterial()
>>> m = m.material_add(phase = ['Ferrite','Martensite'],
...                    O = damask.Rotation.from_random(2),
...                    homogenization = 'SX')
>>> m
material:
  - constituents:
      - O: [0.577764, -0.146299, -0.617669, 0.513010]
        v: 1.0
        phase: Ferrite
    homogenization: SX
  - constituents:
      - O: [0.184176, 0.340305, 0.737247, 0.553840]
        v: 1.0
        phase: Martensite
    homogenization: SX
homogenization: {}
phase: {}

Create a duplex stainless steel microstructure for forming simulations:

>>> import numpy as np
>>> import damask
>>> m = damask.ConfigMaterial()
>>> m = m.material_add(phase = np.array(['Austenite','Ferrite']).reshape(1,2),
...                    O = damask.Rotation.from_random((2,2)),
...                    v = np.array([0.2,0.8]).reshape(1,2),
...                    homogenization = 'Taylor')
>>> m
material:
  - constituents:
      - phase: Austenite
        O: [0.659802978293224, 0.6953785848195171, 0.22426295326327111, -0.17554139512785227]
        v: 0.2
      - phase: Ferrite
        O: [0.49356745891301596, 0.2841806579193434, -0.7487679215072818, -0.339085707289975]
        v: 0.8
    homogenization: Taylor
  - constituents:
      - phase: Austenite
        O: [0.26542221365204055, 0.7268854930702071, 0.4474726435701472, -0.44828201137283735]
        v: 0.2
      - phase: Ferrite
        O: [0.6545817158479885, -0.08004812803625233, -0.6226561293931374, 0.4212059104577611]
        v: 0.8
    homogenization: Taylor
homogenization: {}
phase: {}
clear() None.  Remove all items from D.
copy()

Create deep copy.

delete(keys)

Remove configuration keys.

Parameters
keysiterable or scalar

Label of the key(s) to remove.

Returns
updateddamask.Config

Updated configuration.

fromkeys(iterable, value=None, /)

Create a new dictionary with keys from iterable and values set to value.

get(key, default=None, /)

Return the value for key if key is in the dictionary, else default.

items() a set-like object providing a view on D's items
keys() a set-like object providing a view on D's keys
pop(k[, d]) v, remove specified key and return the corresponding value.

If key is not found, d is returned if given, otherwise KeyError is raised

popitem(/)

Remove and return a (key, value) pair as a 2-tuple.

Pairs are returned in LIFO (last-in, first-out) order. Raises KeyError if the dict is empty.

setdefault(key, default=None, /)

Insert key with a value of default if key is not in the dictionary.

Return the value for key if key is in the dictionary, else default.

update([E, ]**F) None.  Update D from dict/iterable E and F.

If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]

values() an object providing a view on D's values

Grid

class damask.Grid(material, size, origin=[0.0, 0.0, 0.0], comments=[])[source]

Geometry definition for grid solvers.

Create and manipulate geometry definitions for storage as VTK image data files (‘.vti’ extension). A grid contains the material ID (referring to the entry in ‘material.yaml’) and the physical size.

Attributes
N_materials

Number of (unique) material indices within grid.

cells

Number of cells in x,y,z direction.

comments

Comments, e.g.

material

Material indices.

origin

Coordinates of grid origin in meter.

size

Physical size of grid in meter.

Methods

add_primitive(dimension, center, exponent[, ...])

Insert a primitive geometric object at a given position.

canvas([cells, offset, fill])

Crop or enlarge/pad grid.

clean([stencil, selection, periodic])

Smooth grid by selecting most frequent material index within given stencil at each location.

copy()

Create deep copy.

flip(directions)

Flip grid along given directions.

from_Laguerre_tessellation(cells, size, ...)

Create grid from Laguerre tessellation.

from_Voronoi_tessellation(cells, size, seeds)

Create grid from Voronoi tessellation.

from_minimal_surface(cells, size, surface[, ...])

Create grid from definition of triply periodic minimal surface.

from_table(table, coordinates, labels)

Create grid from ASCII table.

get_grain_boundaries([periodic, directions])

Create VTK unstructured grid containing grain boundaries.

load(fname)

Load from VTK image data file.

load_ASCII(fname)

Load from geom file.

load_DREAM3D(fname[, feature_IDs, ...])

Load DREAM.3D (HDF5) file.

load_Neper(fname)

Load from Neper VTK file.

mirror(directions[, reflect])

Mirror grid along given directions.

renumber()

Renumber sorted material indices as 0,...,N-1.

rotate(R[, fill])

Rotate grid (pad if required).

save(fname[, compress])

Save as VTK image data file.

save_ASCII(fname)

Save as geom file.

scale(cells[, periodic])

Scale grid to new cells.

show()

Show on screen.

sort()

Sort material indices such that min(material) is located at (0,0,0).

substitute(from_material, to_material)

Substitute material indices.

vicinity_offset([vicinity, offset, trigger, ...])

Offset material index of points in the vicinity of xxx.

copy()

Create deep copy.

property material

Material indices.

property size

Physical size of grid in meter.

property origin

Coordinates of grid origin in meter.

property comments

Comments, e.g. history of operations.

property cells

Number of cells in x,y,z direction.

property N_materials

Number of (unique) material indices within grid.

static load(fname)[source]

Load from VTK image data file.

Parameters
fnamestr or or pathlib.Path

Grid file to read. Valid extension is .vti, which will be appended if not given.

Returns
loadeddamask.Grid

Grid-based geometry from file.

static load_ASCII(fname)[source]

Load from geom file.

Storing geometry files in ASCII format is deprecated. This function will be removed in a future version of DAMASK.

Parameters
fnamestr, pathlib.Path, or file handle

Geometry file to read.

Returns
loadeddamask.Grid

Grid-based geometry from file.

static load_Neper(fname)[source]

Load from Neper VTK file.

Parameters
fnamestr, pathlib.Path, or file handle

Geometry file to read.

Returns
loadeddamask.Grid

Grid-based geometry from file.

static load_DREAM3D(fname, feature_IDs=None, cell_data=None, phases='Phases', Euler_angles='EulerAngles', base_group=None)[source]

Load DREAM.3D (HDF5) file.

Data in DREAM.3D files can be stored per cell (‘CellData’) and/or per grain (‘Grain Data’). Per default, cell-wise data is assumed.

damask.ConfigMaterial.load_DREAM3D gives the corresponding material definition.

Parameters
fnamestr

Filename of the DREAM.3D (HDF5) file.

feature_IDsstr

Name of the dataset containing the mapping between cells and grain-wise data. Defaults to ‘None’, in which case cell-wise data is used.

cell_datastr

Name of the group (folder) containing cell-wise data. Defaults to None in wich case it is automatically detected.

phasesstr

Name of the dataset containing the phase ID. It is not used for grain-wise data, i.e. when feature_IDs is not None. Defaults to ‘Phases’.

Euler_anglesstr

Name of the dataset containing the crystallographic orientation as Euler angles in radians It is not used for grain-wise data, i.e. when feature_IDs is not None. Defaults to ‘EulerAngles’.

base_groupstr

Path to the group (folder) that contains geometry (_SIMPL_GEOMETRY), and grain- or cell-wise data. Defaults to None, in which case it is set as the path that contains _SIMPL_GEOMETRY/SPACING.

Returns
loadeddamask.Grid

Grid-based geometry from file.

static from_table(table, coordinates, labels)[source]

Create grid from ASCII table.

Parameters
tabledamask.Table

Table that contains material information.

coordinatesstr

Label of the vector column containing the spatial coordinates. Need to be ordered (1./x fast, 3./z slow).

labelsstr or list of str

Label(s) of the columns containing the material definition. Each unique combination of values results in one material ID.

Returns
newdamask.Grid

Grid-based geometry from values in table.

static from_Laguerre_tessellation(cells, size, seeds, weights, material=None, periodic=True)[source]

Create grid from Laguerre tessellation.

Parameters
cellsint numpy.ndarray of shape (3)

Number of cells in x,y,z direction.

sizelist or numpy.ndarray of shape (3)

Physical size of the grid in meter.

seedsnumpy.ndarray of shape (:,3)

Position of the seed points in meter. All points need to lay within the box.

weightsnumpy.ndarray of shape (seeds.shape[0])

Weights of the seeds. Setting all weights to 1.0 gives a standard Voronoi tessellation.

materialnumpy.ndarray of shape (seeds.shape[0]), optional

Material ID of the seeds. Defaults to None, in which case materials are consecutively numbered.

periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

Returns
newdamask.Grid

Grid-based geometry from tessellation.

static from_Voronoi_tessellation(cells, size, seeds, material=None, periodic=True)[source]

Create grid from Voronoi tessellation.

Parameters
cellsint numpy.ndarray of shape (3)

Number of cells in x,y,z direction.

sizelist or numpy.ndarray of shape (3)

Physical size of the grid in meter.

seedsnumpy.ndarray of shape (:,3)

Position of the seed points in meter. All points need to lay within the box.

materialnumpy.ndarray of shape (seeds.shape[0]), optional

Material ID of the seeds. Defaults to None, in which case materials are consecutively numbered.

periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

Returns
newdamask.Grid

Grid-based geometry from tessellation.

static from_minimal_surface(cells, size, surface, threshold=0.0, periods=1, materials=(0, 1))[source]

Create grid from definition of triply periodic minimal surface.

Parameters
cellsint numpy.ndarray of shape (3)

Number of cells in x,y,z direction.

sizelist or numpy.ndarray of shape (3)

Physical size of the grid in meter.

surfacestr

Type of the minimal surface. See notes for details.

thresholdfloat, optional.

Threshold of the minimal surface. Defaults to 0.0.

periodsinteger, optional.

Number of periods per unit cell. Defaults to 1.

materials(int, int), optional

Material IDs. Defaults to (0,1).

Returns
newdamask.Grid

Grid-based geometry from definition of minimal surface.

Notes

The following triply-periodic minimal surfaces are implemented:
  • Schwarz P

  • Double Primitive

  • Schwarz D

  • Complementary D

  • Double Diamond

  • Dprime

  • Gyroid

  • Gprime

  • Karcher K

  • Lidinoid

  • Neovius

  • Fisher-Koch S

References

S.B.G. Blanquer et al., Biofabrication 9(2):025001, 2017 https://doi.org/10.1088/1758-5090/aa6553

M. Wohlgemuth et al., Macromolecules 34(17):6083-6089, 2001 https://doi.org/10.1021/ma0019499

M.-T. Hsieh and L. Valdevit, Software Impacts 6:100026, 2020 https://doi.org/10.1016/j.simpa.2020.100026

Examples

Minimal surface of ‘Gyroid’ type.

>>> import numpy as np
>>> import damask
>>> damask.Grid.from_minimal_surface(np.array([64]*3,int),np.ones(3),
...                                  'Gyroid')
cells  a b c: 64 x 64 x 64
size   x y z: 1.0 x 1.0 x 1.0
origin x y z: 0.0   0.0   0.0
# materials: 2

Minimal surface of ‘Neovius’ type. non-default material IDs.

>>> import numpy as np
>>> import damask
>>> damask.Grid.from_minimal_surface(np.array([80]*3,int),np.ones(3),
...                                  'Neovius',materials=(1,5))
cells  a b c: 80 x 80 x 80
size   x y z: 1.0 x 1.0 x 1.0
origin x y z: 0.0   0.0   0.0
# materials: 2 (min: 1, max: 5)
save(fname, compress=True)[source]

Save as VTK image data file.

Parameters
fnamestr or pathlib.Path

Filename to write. Valid extension is .vti, it will be appended if not given.

compressbool, optional

Compress with zlib algorithm. Defaults to True.

save_ASCII(fname)[source]

Save as geom file.

Storing geometry files in ASCII format is deprecated. This function will be removed in a future version of DAMASK.

Parameters
fnamestr or file handle

Geometry file to write with extension ‘.geom’.

compressbool, optional

Compress geometry with ‘x of y’ and ‘a to b’.

show()[source]

Show on screen.

add_primitive(dimension, center, exponent, fill=None, R=Quaternion [1. 0. 0. 0.], inverse=False, periodic=True)[source]

Insert a primitive geometric object at a given position.

Parameters
dimensionint or float numpy.ndarray of shape (3)

Dimension (diameter/side length) of the primitive. If given as integers, cell centers are addressed. If given as floats, coordinates are addressed.

centerint or float numpy.ndarray of shape (3)

Center of the primitive. If given as integers, cell centers are addressed. If given as floats, coordinates in space are addressed.

exponentnumpy.ndarray of shape (3) or float

Exponents for the three axes. 0 gives octahedron (ǀxǀ^(2^0) + ǀyǀ^(2^0) + ǀzǀ^(2^0) < 1) 1 gives sphere (ǀxǀ^(2^1) + ǀyǀ^(2^1) + ǀzǀ^(2^1) < 1)

fillint, optional

Fill value for primitive. Defaults to material.max()+1.

Rdamask.Rotation, optional

Rotation of primitive. Defaults to no rotation.

inverseBoolean, optional

Retain original materials within primitive and fill outside. Defaults to False.

periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

Returns
updateddamask.Grid

Updated grid-based geometry.

Examples

Add a sphere at the center.

>>> import numpy as np
>>> import damask
>>> g = damask.Grid(np.zeros([64]*3,int), np.ones(3)*1e-4)
>>> g.add_primitive(np.ones(3)*5e-5,np.ones(3)*5e-5,1)
cells  a b c: 64 x 64 x 64
size   x y z: 0.0001 x 0.0001 x 0.0001
origin x y z: 0.0   0.0   0.0
# materials: 2

Add a cube at the origin.

>>> import numpy as np
>>> import damask
>>> g = damask.Grid(np.zeros([64]*3,int), np.ones(3)*1e-4)
>>> g.add_primitive(np.ones(3,int)*32,np.zeros(3),np.inf)
cells  a b c: 64 x 64 x 64
size   x y z: 0.0001 x 0.0001 x 0.0001
origin x y z: 0.0   0.0   0.0
# materials: 2
mirror(directions, reflect=False)[source]

Mirror grid along given directions.

Parameters
directionsiterable containing str

Direction(s) along which the grid is mirrored. Valid entries are ‘x’, ‘y’, ‘z’.

reflectbool, optional

Reflect (include) outermost layers. Defaults to False.

Returns
updateddamask.Grid

Updated grid-based geometry.

Examples

Mirror along x- and y-direction.

>>> import numpy as np
>>> import damask
>>> g = damask.Grid(np.zeros([32]*3,int), np.ones(3)*1e-4)
>>> g.mirror('xy',True)
cells  a b c: 64 x 64 x 32
size   x y z: 0.0002 x 0.0002 x 0.0001
origin x y z: 0.0   0.0   0.0
# materials: 1
flip(directions)[source]

Flip grid along given directions.

Parameters
directionsiterable containing str

Direction(s) along which the grid is flipped. Valid entries are ‘x’, ‘y’, ‘z’.

Returns
updateddamask.Grid

Updated grid-based geometry.

scale(cells, periodic=True)[source]

Scale grid to new cells.

Parameters
cellsnumpy.ndarray of shape (3)

Number of cells in x,y,z direction.

periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

Returns
updateddamask.Grid

Updated grid-based geometry.

Examples

Double resolution.

>>> import numpy as np
>>> import damask
>>> g = damask.Grid(np.zeros([32]*3,int),np.ones(3)*1e-4)
>>> g.scale(g.cells*2)
cells  a b c: 64 x 64 x 64
size   x y z: 0.0001 x 0.0001 x 0.0001
origin x y z: 0.0   0.0   0.0
# materials: 1
clean(stencil=3, selection=None, periodic=True)[source]

Smooth grid by selecting most frequent material index within given stencil at each location.

Parameters
stencilint, optional

Size of smoothing stencil.

selectionlist, optional

Field values that can be altered. Defaults to all.

periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

Returns
updateddamask.Grid

Updated grid-based geometry.

renumber()[source]

Renumber sorted material indices as 0,…,N-1.

Returns
updateddamask.Grid

Updated grid-based geometry.

rotate(R, fill=None)[source]

Rotate grid (pad if required).

Parameters
Rdamask.Rotation

Rotation to apply to the grid.

fillint or float, optional

Material index to fill the corners. Defaults to material.max() + 1.

Returns
updateddamask.Grid

Updated grid-based geometry.

canvas(cells=None, offset=None, fill=None)[source]

Crop or enlarge/pad grid.

Parameters
cellsnumpy.ndarray of shape (3)

Number of cells x,y,z direction.

offsetnumpy.ndarray of shape (3)

Offset (measured in cells) from old to new grid [0,0,0].

fillint or float, optional

Material index to fill the background. Defaults to material.max() + 1.

Returns
updateddamask.Grid

Updated grid-based geometry.

Examples

Remove 1/2 of the microstructure in z-direction.

>>> import numpy as np
>>> import damask
>>> g = damask.Grid(np.zeros([32]*3,int),np.ones(3)*1e-4)
>>> g.canvas(np.array([32,32,16],int))
cells  a b c: 33 x 32 x 16
size   x y z: 0.0001 x 0.0001 x 5e-05
origin x y z: 0.0   0.0   0.0
# materials: 1
substitute(from_material, to_material)[source]

Substitute material indices.

Parameters
from_materialiterable of ints

Material indices to be substituted.

to_materialiterable of ints

New material indices.

Returns
updateddamask.Grid

Updated grid-based geometry.

sort()[source]

Sort material indices such that min(material) is located at (0,0,0).

Returns
updateddamask.Grid

Updated grid-based geometry.

vicinity_offset(vicinity=1, offset=None, trigger=[], periodic=True)[source]

Offset material index of points in the vicinity of xxx.

Different from themselves (or listed as triggers) within a given (cubic) vicinity, i.e. within the region close to a grain/phase boundary. ToDo: use include/exclude as in seeds.from_grid

Parameters
vicinityint, optional

Voxel distance checked for presence of other materials. Defaults to 1.

offsetint, optional

Offset (positive or negative) to tag material indices, defaults to material.max()+1.

triggerlist of ints, optional

List of material indices that trigger a change. Defaults to [], meaning that any different neighbor triggers a change.

periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

Returns
updateddamask.Grid

Updated grid-based geometry.

get_grain_boundaries(periodic=True, directions='xyz')[source]

Create VTK unstructured grid containing grain boundaries.

Parameters
periodicBoolean, optional

Assume grid to be periodic. Defaults to True.

directionsiterable containing str, optional

Direction(s) along which the boundaries are determined. Valid entries are ‘x’, ‘y’, ‘z’. Defaults to ‘xyz’.

Returns
grain_boundariesdamask.VTK

VTK-based geometry of grain boundary network.


seeds

Functionality for generation of seed points for Voronoi or Laguerre tessellation.

damask.seeds.from_random(size, N_seeds, cells=None, rng_seed=None)[source]

Place seeds randomly in space.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the seeding domain.

N_seedsint

Number of seeds.

cellsnumpy.ndarray of shape (3), optional.

If given, ensures that each seed results in a grain when a standard Voronoi tessellation is performed using the given grid resolution (i.e. size/cells).

rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional

A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS.

Returns
coordsnumpy.ndarray of shape (N_seeds,3)

Seed coordinates in 3D space.

damask.seeds.from_Poisson_disc(size, N_seeds, N_candidates, distance, periodic=True, rng_seed=None)[source]

Place seeds according to a Poisson disc distribution.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the seeding domain.

N_seedsint

Number of seeds.

N_candidatesint

Number of candidates to consider for finding best candidate.

distancefloat

Minimum acceptable distance to other seeds.

periodicboolean, optional

Calculate minimum distance for periodically repeated grid.

rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional

A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS.

Returns
coordsnumpy.ndarray of shape (N_seeds,3)

Seed coordinates in 3D space.

damask.seeds.from_grid(grid, selection=None, invert=False, average=False, periodic=True)[source]

Create seeds from grid description.

Parameters
griddamask.Grid

Grid, from which the material IDs are used as seeds.

selectioniterable of integers, optional

Material IDs to consider.

invertboolean, false

Do not consider the material IDs given in selection. Defaults to False.

averageboolean, optional

Seed corresponds to center of gravity of material ID cloud.

periodicboolean, optional

Center of gravity with periodic boundaries.

Returns
coords, materialsnumpy.ndarray of shape (:,3), numpy.ndarray of shape (:)

Seed coordinates in 3D space, material IDs.


Rotation

class damask.Rotation(rotation=array([1., 0., 0., 0.]))[source]

Rotation with functionality for conversion between 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 angular ranges of [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).

References

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

Examples

Rotate vector ‘a’ (defined in coordinate system ‘A’) to coordinates ‘b’ expressed in system ‘B’:

>>> import damask
>>> import numpy as np
>>> Q = damask.Rotation.from_random()
>>> a = np.random.rand(3)
>>> b = Q @ a
>>> np.allclose(np.dot(Q.as_matrix(),a),b)
True

Compound rotations R1 (first) and R2 (second):

>>> import damask
>>> import numpy as np
>>> R1 = damask.Rotation.from_random()
>>> R2 = damask.Rotation.from_random()
>>> R = R2 * R1
>>> np.allclose(R.as_matrix(), np.dot(R2.as_matrix(),R1.as_matrix()))
True
Attributes
quaternion
shape
size

Methods

allclose(other[, rtol, atol, equal_nan])

Test whether all values are approximately equal to corresponding ones of other Rotation.

append(other)

Extend array along first dimension with other array(s).

apply(other)

Rotate vector, second order tensor, or fourth order tensor.

as_Euler_angles([degrees])

Represent as Bunge Euler angles.

as_Rodrigues_vector([compact])

Represent as Rodrigues–Frank vector with separate axis and angle argument.

as_axis_angle([degrees, pair])

Represent as axis–angle pair.

as_cubochoric()

Represent as cubochoric vector.

as_homochoric()

Represent as homochoric vector.

as_matrix()

Represent as rotation matrix.

as_quaternion()

Represent as unit quaternion.

average([weights])

Average along last array dimension.

broadcast_to(shape[, mode])

Broadcast array.

copy(**kwargs)

Create deep copy.

flatten([order])

Flatten array.

from_Euler_angles(phi[, degrees])

Initialize from Bunge Euler angles.

from_ODF(weights, phi[, N, degrees, ...])

Sample discrete values from a binned orientation distribution function (ODF).

from_Rodrigues_vector(rho[, normalize, P])

Initialize from Rodrigues–Frank vector (angle separated from axis).

from_axis_angle(axis_angle[, degrees, ...])

Initialize from Axis angle pair.

from_basis(basis[, orthonormal, reciprocal])

Initialize from lattice basis vectors.

from_cubochoric(x[, P])

Initialize from cubochoric vector.

from_fiber_component(alpha, beta[, sigma, ...])

Calculate set of rotations with Gaussian distribution around direction.

from_homochoric(h[, P])

Initialize from homochoric vector.

from_matrix(R)

Initialize from rotation matrix.

from_parallel(a, b, **kwargs)

Initialize from pairs of two orthogonal lattice basis vectors.

from_quaternion(q[, accept_homomorph, P])

Initialize from quaternion.

from_random([shape, rng_seed])

Initialize with random rotation.

from_spherical_component(center, sigma[, N, ...])

Calculate set of rotations with Gaussian distribution around center.

isclose(other[, rtol, atol, equal_nan])

Report where values are approximately equal to corresponding ones of other Rotation.

misorientation(other)

Calculate misorientation to other Rotation.

reshape(shape[, order])

Reshape array.

copy(**kwargs)

Create deep copy.

isclose(other, rtol=1e-05, atol=1e-08, equal_nan=True)[source]

Report where values are approximately equal to corresponding ones of other Rotation.

Parameters
otherRotation

Rotation to compare against.

rtolfloat, optional

Relative tolerance of equality.

atolfloat, optional

Absolute tolerance of equality.

equal_nanbool, optional

Consider matching NaN values as equal. Defaults to True.

Returns
masknumpy.ndarray bool

Mask indicating where corresponding rotations are close.

allclose(other, rtol=1e-05, atol=1e-08, equal_nan=True)[source]

Test whether all values are approximately equal to corresponding ones of other Rotation.

Parameters
otherRotation

Rotation to compare against.

rtolfloat, optional

Relative tolerance of equality.

atolfloat, optional

Absolute tolerance of equality.

equal_nanbool, optional

Consider matching NaN values as equal. Defaults to True.

Returns
answerbool

Whether all values are close between both rotations.

apply(other)

Rotate vector, second order tensor, or fourth order tensor.

Parameters
othernumpy.ndarray of shape (…,3), (…,3,3), or (…,3,3,3,3)

Vector or tensor on which to apply the rotation.

Returns
rotatednumpy.ndarray of shape (…,3), (…,3,3), or (…,3,3,3,3)

Rotated vector or tensor, i.e. transformed to frame defined by rotation.

append(other)[source]

Extend array along first dimension with other array(s).

Parameters
otherdamask.Rotation
flatten(order='C')[source]

Flatten array.

Returns
flatteneddamask.Rotation

Rotation flattened to single dimension.

reshape(shape, order='C')[source]

Reshape array.

Returns
reshapeddamask.Rotation

Rotation of given shape.

broadcast_to(shape, mode='right')[source]

Broadcast array.

Parameters
shapetuple

Shape of broadcasted array.

modestr, optional

Where to preferentially locate missing dimensions. Either ‘left’ or ‘right’ (default).

Returns
broadcasteddamask.Rotation

Rotation broadcasted to given shape.

average(weights=None)[source]

Average along last array dimension.

Parameters
weightslist of floats, optional

Relative weight of each rotation.

Returns
averagedamask.Rotation

Weighted average of original Rotation field.

References

F. Landis Markley et al., Journal of Guidance, Control, and Dynamics 30(4):1193-1197, 2007 https://doi.org/10.2514/1.28949

misorientation(other)[source]

Calculate misorientation to other Rotation.

Parameters
otherdamask.Rotation

Rotation to which the misorientation is computed.

Returns
gdamask.Rotation

Misorientation.

as_quaternion()[source]

Represent as unit quaternion.

Returns
qnumpy.ndarray of shape (…,4)

Unit quaternion (q_0, q_1, q_2, q_3) in positive real hemisphere, i.e. ǀqǀ = 1, q_0 ≥ 0.

as_Euler_angles(degrees=False)[source]

Represent as Bunge Euler angles.

Parameters
degreesbool, optional

Return angles in degrees.

Returns
phinumpy.ndarray of shape (…,3)

Bunge Euler angles (φ_1 ∈ [0,2π], ϕ ∈ [0,π], φ_2 ∈ [0,2π]) or (φ_1 ∈ [0,360], ϕ ∈ [0,180], φ_2 ∈ [0,360]) if degrees == True.

Notes

Bunge Euler angles correspond to a rotation axis sequence of z–x’–z’’.

Examples

Cube orientation as Bunge Euler angles.

>>> import damask
>>> import numpy as np
>>> damask.Rotation(np.array([1,0,0,0])).as_Euler_angles()
array([0., 0., 0.])
as_axis_angle(degrees=False, pair=False)[source]

Represent as axis–angle pair.

Parameters
degreesbool, optional

Return rotation angle in degrees. Defaults to False.

pairbool, optional

Return tuple of axis and angle. Defaults to False.

Returns
axis_anglenumpy.ndarray of shape (…,4) or tuple ((…,3), (…)) if pair == True

Axis and angle [n_1, n_2, n_3, ω] with ǀnǀ = 1 and ω ∈ [0,π] or ω ∈ [0,180] if degrees == True.

Examples

Cube orientation as axis–angle pair.

>>> import damask
>>> import numpy as np
>>> damask.Rotation(np.array([1,0,0,0])).as_axis_angle(pair=True)
(array([0., 0., 1.]), array(0.))
as_matrix()[source]

Represent as rotation matrix.

Returns
Rnumpy.ndarray of shape (…,3,3)

Rotation matrix R with det(R) = 1, R.T ∙ R = I.

Examples

Cube orientation as rotation matrix.

>>> import damask
>>> import numpy as np
>>> damask.Rotation(np.array([1,0,0,0])).as_matrix()
array([[1., 0., 0.],
       [0., 1., 0.],
       [0., 0., 1.]])
as_Rodrigues_vector(compact=False)[source]

Represent as Rodrigues–Frank vector with separate axis and angle argument.

Parameters
compactbool, optional

Return three-component Rodrigues–Frank vector, i.e. axis and angle argument are not separated.

Returns
rhonumpy.ndarray of shape (…,4) or (…,3) if compact == True

Rodrigues–Frank vector [n_1, n_2, n_3, tan(ω/2)] with ǀnǀ = 1 and ω ∈ [0,π] or [n_1, n_2, n_3] with ǀnǀ = tan(ω/2) and ω ∈ [0,π] if compact == True.

Examples

Cube orientation as three-component Rodrigues–Frank vector.

>>> import damask
>>> import numpy as np
>>> damask.Rotation(np.array([1,0,0,0])).as_Rodrigues_vector(compact=True)
array([ 0.,  0., 0.])
as_homochoric()[source]

Represent as homochoric vector.

Returns
hnumpy.ndarray of shape (…,3)

Homochoric vector (h_1, h_2, h_3) with ǀhǀ < (3/4*π)^(1/3).

Examples

Cube orientation as homochoric vector.

>>> import damask
>>> import numpy as np
>>> damask.Rotation(np.array([1,0,0,0])).as_homochoric()
array([0., 0., 0.])
as_cubochoric()[source]

Represent as cubochoric vector.

Returns
xnumpy.ndarray of shape (…,3)

Cubochoric vector (x_1, x_2, x_3) with max(x_i) < 1/2*π^(2/3).

Examples

Cube orientation as cubochoric vector.

>>> import damask
>>> import numpy as np
>>> damask.Rotation(np.array([1,0,0,0])).as_cubochoric()
array([0., 0., 0.])
static from_quaternion(q, accept_homomorph=False, P=- 1)[source]

Initialize from quaternion.

Parameters
qnumpy.ndarray of shape (…,4)

Unit quaternion (q_0, q_1, q_2, q_3) in positive real hemisphere, i.e. ǀqǀ = 1, q_0 ≥ 0.

accept_homomorphboolean, optional

Allow homomorphic variants, i.e. q_0 < 0 (negative real hemisphere). Defaults to False.

Pint ∈ {-1,1}, optional

Sign convention. Defaults to -1.

static from_Euler_angles(phi, degrees=False)[source]

Initialize from Bunge Euler angles.

Parameters
phinumpy.ndarray of shape (…,3)

Euler angles (φ_1 ∈ [0,2π], ϕ ∈ [0,π], φ_2 ∈ [0,2π]) or (φ_1 ∈ [0,360], ϕ ∈ [0,180], φ_2 ∈ [0,360]) if degrees == True.

degreesboolean, optional

Euler angles are given in degrees. Defaults to False.

Notes

Bunge Euler angles correspond to a rotation axis sequence of z–x’–z’’.

static from_axis_angle(axis_angle, degrees=False, normalize=False, P=- 1)[source]

Initialize from Axis angle pair.

Parameters
axis_anglenumpy.ndarray of shape (…,4)

Axis and angle (n_1, n_2, n_3, ω) with ǀnǀ = 1 and ω ∈ [0,π] or ω ∈ [0,180] if degrees == True.

degreesboolean, optional

Angle ω is given in degrees. Defaults to False.

normalize: boolean, optional

Allow ǀnǀ ≠ 1. Defaults to False.

Pint ∈ {-1,1}, optional

Sign convention. Defaults to -1.

static from_basis(basis, orthonormal=True, reciprocal=False)[source]

Initialize from lattice basis vectors.

Parameters
basisnumpy.ndarray of shape (…,3,3)

Three three-dimensional lattice basis vectors.

orthonormalboolean, optional

Basis is strictly orthonormal, i.e. is free of stretch components. Defaults to True.

reciprocalboolean, optional

Basis vectors are given in reciprocal (instead of real) space. Defaults to False.

static from_matrix(R)[source]

Initialize from rotation matrix.

Parameters
Rnumpy.ndarray of shape (…,3,3)

Rotation matrix with det(R) = 1, R.T ∙ R = I.

static from_parallel(a, b, **kwargs)[source]

Initialize from pairs of two orthogonal lattice basis vectors.

Parameters
anumpy.ndarray of shape (…,2,3)

Two three-dimensional lattice vectors of first orthogonal basis.

bnumpy.ndarray of shape (…,2,3)

Corresponding three-dimensional lattice vectors of second basis.

static from_Rodrigues_vector(rho, normalize=False, P=- 1)[source]

Initialize from Rodrigues–Frank vector (angle separated from axis).

Parameters
rhonumpy.ndarray of shape (…,4)

Rodrigues–Frank vector (n_1, n_2, n_3, tan(ω/2)) with ǀnǀ = 1 and ω ∈ [0,π].

normalizeboolean, optional

Allow ǀnǀ ≠ 1. Defaults to False.

Pint ∈ {-1,1}, optional

Sign convention. Defaults to -1.

static from_homochoric(h, P=- 1)[source]

Initialize from homochoric vector.

Parameters
hnumpy.ndarray of shape (…,3)

Homochoric vector (h_1, h_2, h_3) with ǀhǀ < (3/4*π)^(1/3).

Pint ∈ {-1,1}, optional

Sign convention. Defaults to -1.

static from_cubochoric(x, P=- 1)[source]

Initialize from cubochoric vector.

Parameters
xnumpy.ndarray of shape (…,3)

Cubochoric vector (x_1, x_2, x_3) with max(x_i) < 1/2*π^(2/3).

Pint ∈ {-1,1}, optional

Sign convention. Defaults to -1.

static from_random(shape=None, rng_seed=None)[source]

Initialize with random rotation.

Rotations are uniformly distributed.

Parameters
shapetuple of ints, optional

Shape of the sample. Defaults to None, which gives a single rotation.

rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional

A seed to initialize the BitGenerator. Defaults to None, i.e. unpredictable entropy will be pulled from the OS.

static from_ODF(weights, phi, N=500, degrees=True, fractions=True, rng_seed=None, **kwargs)[source]

Sample discrete values from a binned orientation distribution function (ODF).

Parameters
weightsnumpy.ndarray of shape (n)

Texture intensity values (probability density or volume fraction) at Euler space grid points.

phinumpy.ndarray of shape (n,3)

Grid coordinates in Euler space at which weights are defined.

Ninteger, optional

Number of discrete orientations to be sampled from the given ODF. Defaults to 500.

degreesboolean, optional

Euler space grid coordinates are in degrees. Defaults to True.

fractionsboolean, optional

ODF values correspond to volume fractions, not probability densities. Defaults to True.

rng_seed: {None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional

A seed to initialize the BitGenerator. Defaults to None, i.e. unpredictable entropy will be pulled from the OS.

Returns
samplesdamask.Rotation of shape (N)

Array of sampled rotations closely representing the input ODF.

Notes

Due to the distortion of Euler space in the vicinity of ϕ = 0, probability densities, p, defined on grid points with ϕ = 0 will never result in reconstructed orientations as their dV/V = p dγ = p × 0. Hence, it is recommended to transform any such dataset to cell centers that avoid grid points at ϕ = 0.

References

P. Eisenlohr and F. Roters, Computational Materials Science 42(4):670-678, 2008 https://doi.org/10.1016/j.commatsci.2007.09.015

static from_spherical_component(center, sigma, N=500, degrees=True, rng_seed=None)[source]

Calculate set of rotations with Gaussian distribution around center.

Parameters
centerRotation

Central Rotation.

sigmafloat

Standard deviation of (Gaussian) misorientation distribution.

Nint, optional

Number of samples. Defaults to 500.

degreesboolean, optional

sigma is given in degrees. Defaults to True.

rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional

A seed to initialize the BitGenerator. Defaults to None, i.e. unpredictable entropy will be pulled from the OS.

static from_fiber_component(alpha, beta, sigma=0.0, N=500, degrees=True, rng_seed=None)[source]

Calculate set of rotations with Gaussian distribution around direction.

Parameters
alphanumpy.ndarray of shape (2)

Polar coordinates (phi from x, theta from z) of fiber direction in crystal frame.

betanumpy.ndarray of shape (2)

Polar coordinates (phi from x, theta from z) of fiber direction in sample frame.

sigmafloat, optional

Standard deviation of (Gaussian) misorientation distribution. Defaults to 0.

Nint, optional

Number of samples. Defaults to 500.

degreesboolean, optional

sigma, alpha, and beta are given in degrees.

rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional

A seed to initialize the BitGenerator. Defaults to None, i.e. unpredictable entropy will be pulled from the OS.