Post-Processing

The following classes and modules support the analysis of results from DAMASK simulations.

Result

class damask.Result(fname)[source]

Add data to and export data from a DADF5 file.

A DADF5 (DAMASK HDF5) file contains DAMASK results. Its group/folder structure reflects the layout in material.yaml.

This class provides a customizable view on the DADF5 file. Upon initialization, all attributes are visible. Derived quantities are added to the file and existing data is exported based on the current view.

Examples

Open ‘my_file.hdf5’, which is assumed to contain deformation gradient ‘F’ and first Piola-Kirchhoff stress ‘P’, add the Mises equivalent of the Cauchy stress, and export it to VTK (file) and numpy.ndarray (memory).

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_Cauchy()
>>> r.add_equivalent_Mises('sigma')
>>> r.export_VTK()
>>> r_last = r.view('increments',-1)
>>> sigma_vM_last = r_last.get('sigma_vM')
Attributes
coordinates0_node

Initial/undeformed nodal coordinates.

coordinates0_point

Initial/undeformed cell center coordinates.

geometry0

Initial/undeformed geometry.

Methods

add_IPF_color(l[, q])

Add RGB color tuple of inverse pole figure (IPF) color.

add_absolute(x)

Add absolute value.

add_calculation(formula, name[, unit, …])

Add result of a general formula.

add_curl(f)

Add curl of a field.

add_determinant(T)

Add the determinant of a tensor.

add_deviator(T)

Add the deviatoric part of a tensor.

add_divergence(f)

Add divergence of a field.

add_eigenvalue(T_sym[, eigenvalue])

Add eigenvalues of symmetric tensor.

add_eigenvector(T_sym[, eigenvalue])

Add eigenvector of symmetric tensor.

add_equivalent_Mises(T_sym[, kind])

Add the equivalent Mises stress or strain of a symmetric tensor.

add_gradient(f)

Add gradient of a field.

add_maximum_shear(T_sym)

Add maximum shear components of symmetric tensor.

add_norm(x[, ord])

Add the norm of vector or tensor.

add_rotation(F)

Add rotational part of a deformation gradient.

add_spherical(T)

Add the spherical (hydrostatic) part of a tensor.

add_strain([F, t, m])

Add strain tensor of a deformation gradient.

add_stress_Cauchy([P, F])

Add Cauchy stress calculated from first Piola-Kirchhoff stress and deformation gradient.

add_stress_second_Piola_Kirchhoff([P, F])

Add second Piola-Kirchhoff stress calculated from first Piola-Kirchhoff stress and deformation gradient.

add_stretch_tensor([F, t])

Add stretch tensor of a deformation gradient.

copy()

Create deep copy.

export_VTK([output, mode, constituents, …])

Export to VTK cell/point data.

export_XDMF([output])

Write XDMF file to directly visualize data in DADF5 file.

get([output, flatten, prune])

Collect data per phase/homogenization reflecting the group/folder structure in the DADF5 file.

increments_in_range(start, end)

Get all increments within a given range.

list_data()

Return information on all active datasets in the file.

modification_disable()

Prevent modification of existing data (default case).

modification_enable()

Allow modification of existing data.

place([output, flatten, prune, …])

Merge data into spatial order that is compatible with the damask.VTK geometry representation.

remove(name)

Remove/delete datasets.

rename(name_src, name_dst)

Rename/move datasets (within the same group/folder).

save_VTK([output, mode, constituents, …])

Export to VTK cell/point data.

save_XDMF([output])

Write XDMF file to directly visualize data in DADF5 file.

times_in_range(start, end)

Get all increments within a given time range.

view(what, datasets)

Set view.

view_less(what, datasets)

Remove from view.

view_more(what, datasets)

Add to view.

enable_user_function

copy()

Create deep copy.

modification_enable()[source]

Allow modification of existing data.

Returns
modified_viewdamask.Result

View without write-protection of existing data.

modification_disable()[source]

Prevent modification of existing data (default case).

Returns
modified_viewdamask.Result

View with write-protection of existing data.

increments_in_range(start, end)[source]

Get all increments within a given range.

Parameters
startint or str

Start increment.

endint or str

End increment.

Returns
incrementslist of ints

Increment number of all increments within the given bounds.

times_in_range(start, end)[source]

Get all increments within a given time range.

Parameters
startfloat

Time of start increment.

endfloat

Time of end increment.

Returns
timeslist of float

Simulation time of all increments within the given bounds.

view(what, datasets)[source]

Set view.

Parameters
what{‘increments’, ‘times’, ‘phases’, ‘homogenizations’, ‘fields’}

Attribute to change.

datasets(list of) int (for increments), (list of) float (for times), (list of) str, or bool

Name of datasets; supports ‘?’ and ‘*’ wildcards. True is equivalent to ‘*’, False is equivalent to [].

Returns
viewdamask.Result

View with only the selected attributes being visible.

Examples

Get a view that shows only results from the initial configuration:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r_first = r.view('increment',0)

Get a view that shows all results between simulation times of 10 to 40:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r_t10to40 = r.view('times',r.times_in_range(10.0,40.0))
view_more(what, datasets)[source]

Add to view.

Parameters
what{‘increments’, ‘times’, ‘phases’, ‘homogenizations’, ‘fields’}

Attribute to change.

datasets(list of) int (for increments), (list of) float (for times), (list of) str, or bool

Name of datasets; supports ‘?’ and ‘*’ wildcards. True is equivalent to ‘*’, False is equivalent to [].

Returns
modified_viewdamask.Result

View with additional visible attributes.

Examples

Get a view that shows only results from first and last increment:

>>> import damask
>>> r_empty = damask.Result('my_file.hdf5').view('increments',False)
>>> r_first = r_empty.view_more('increments',0)
>>> r_first_and_last = r.first.view_more('increments',-1)
view_less(what, datasets)[source]

Remove from view.

Parameters
what{‘increments’, ‘times’, ‘phases’, ‘homogenizations’, ‘fields’}

Attribute to change.

datasets(list of) int (for increments), (list of) float (for times), (list of) str, or bool

Name of datasets; supports ‘?’ and ‘*’ wildcards. True is equivalent to ‘*’, False is equivalent to [].

Returns
modified_viewdamask.Result

View with fewer visible attributes.

Examples

Get a view that omits the undeformed configuration:

>>> import damask
>>> r_all = damask.Result('my_file.hdf5')
>>> r_deformed = r_all.view_less('increments',0)
rename(name_src, name_dst)[source]

Rename/move datasets (within the same group/folder).

This operation is discouraged because the history of the data becomes untraceable and data integrity is not ensured.

Parameters
name_srcstr

Name of the datasets to be renamed.

name_dststr

New name of the datasets.

Examples

Rename datasets containing the deformation gradient from ‘F’ to ‘def_grad’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r_unprotected = r.modification_enable()
>>> r_unprotected.rename('F','def_grad')
remove(name)[source]

Remove/delete datasets.

This operation is discouraged because the history of the data becomes untraceable and data integrity is not ensured.

Parameters
namestr

Name of the datasets to be deleted.

Examples

Delete the deformation gradient ‘F’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r_unprotected = r.modification_enable()
>>> r_unprotected.remove('F')
list_data()[source]

Return information on all active datasets in the file.

property coordinates0_point

Initial/undeformed cell center coordinates.

property coordinates0_node

Initial/undeformed nodal coordinates.

property geometry0

Initial/undeformed geometry.

add_absolute(x)[source]

Add absolute value.

Parameters
xstr

Name of scalar, vector, or tensor dataset to take absolute value of.

add_calculation(formula, name, unit='n/a', description=None)[source]

Add result of a general formula.

Parameters
formulastr

Formula to calculate resulting dataset. Existing datasets are referenced by ‘#TheirName#’.

namestr

Name of resulting dataset.

unitstr, optional

Physical unit of the result.

descriptionstr, optional

Human-readable description of the result.

Examples

Add total dislocation density, i.e. the sum of mobile dislocation density ‘rho_mob’ and dislocation dipole density ‘rho_dip’ over all slip systems:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_calculation('np.sum(#rho_mob#,axis=1)','rho_mob_total',
...                    '1/m²','total mobile dislocation density')
>>> r.add_calculation('np.sum(#rho_dip#,axis=1)','rho_dip_total',
...                    '1/m²','total dislocation dipole density')
>>> r.add_calculation('#rho_dip_total#+#rho_mob_total','rho_total',
...                    '1/m²','total dislocation density')

Add Mises equivalent of the Cauchy stress without storage of intermediate results. Define a user function for better readability:

>>> import damask
>>> def equivalent_stress(F,P):
...     sigma = damask.mechanics.stress_Cauchy(F=F,P=P)
...     return damask.mechanics.equivalent_stress_Mises(sigma)
>>> r = damask.Result('my_file.hdf5')
>>> r.enable_user_function(equivalent_stress)
>>> r.add_calculation('equivalent_stress(#F#,#P#)','sigma_vM','Pa',
...                   'Mises equivalent of the Cauchy stress')
add_stress_Cauchy(P='P', F='F')[source]

Add Cauchy stress calculated from first Piola-Kirchhoff stress and deformation gradient.

Parameters
Pstr, optional

Name of the dataset containing the first Piola-Kirchhoff stress. Defaults to ‘P’.

Fstr, optional

Name of the dataset containing the deformation gradient. Defaults to ‘F’.

add_determinant(T)[source]

Add the determinant of a tensor.

Parameters
Tstr

Name of tensor dataset.

Examples

Add the determinant of plastic deformation gradient ‘F_p’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_determinant('F_p')
add_deviator(T)[source]

Add the deviatoric part of a tensor.

Parameters
Tstr

Name of tensor dataset.

add_eigenvalue(T_sym, eigenvalue='max')[source]

Add eigenvalues of symmetric tensor.

Parameters
T_symstr

Name of symmetric tensor dataset.

eigenvalue{‘max’, ‘mid’, ‘min’}

Eigenvalue. Defaults to ‘max’.

add_eigenvector(T_sym, eigenvalue='max')[source]

Add eigenvector of symmetric tensor.

Parameters
T_symstr

Name of symmetric tensor dataset.

eigenvalue{‘max’, ‘mid’, ‘min’}

Eigenvalue to which the eigenvector corresponds. Defaults to ‘max’.

add_IPF_color(l, q='O')[source]

Add RGB color tuple of inverse pole figure (IPF) color.

Parameters
lnumpy.array of shape (3)

Lab frame direction for inverse pole figure.

qstr

Name of the dataset containing the crystallographic orientation as quaternions. Defaults to ‘O’.

Examples

Add the IPF color along [0,1,1] for orientation ‘O’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_IPF_color(np.array([0,1,1]))
add_maximum_shear(T_sym)[source]

Add maximum shear components of symmetric tensor.

Parameters
T_symstr

Name of symmetric tensor dataset.

add_equivalent_Mises(T_sym, kind=None)[source]

Add the equivalent Mises stress or strain of a symmetric tensor.

Parameters
T_symstr

Name of symmetric tensorial stress or strain dataset.

kind{‘stress’, ‘strain’, None}, optional

Kind of the von Mises equivalent. Defaults to None, in which case it is selected based on the unit of the dataset (‘1’ -> strain, ‘Pa’ -> stress).

Examples

Add the Mises equivalent of the Cauchy stress ‘sigma’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_equivalent_Mises('sigma')

Add the Mises equivalent of the spatial logarithmic strain ‘epsilon_V^0.0(F)’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_equivalent_Mises('epsilon_V^0.0(F)')
add_norm(x, ord=None)[source]

Add the norm of vector or tensor.

Parameters
xstr

Name of vector or tensor dataset.

ord{non-zero int, inf, -inf, ‘fro’, ‘nuc’}, optional

Order of the norm. inf means NumPy’s inf object. For details refer to numpy.linalg.norm.

add_stress_second_Piola_Kirchhoff(P='P', F='F')[source]

Add second Piola-Kirchhoff stress calculated from first Piola-Kirchhoff stress and deformation gradient.

Parameters
Pstr, optional

Name of first Piola-Kirchhoff stress dataset. Defaults to ‘P’.

Fstr, optional

Name of deformation gradient dataset. Defaults to ‘F’.

Notes

The definition of the second Piola-Kirchhoff stress (S = [F^-1 P]_sym) follows the standard definition in nonlinear continuum mechanics. As such, no intermediate configuration, for instance that reached by F_p, is taken into account.

add_rotation(F)[source]

Add rotational part of a deformation gradient.

Parameters
Fstr

Name of deformation gradient dataset.

Examples

Add the rotational part of deformation gradient ‘F’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_rotation('F')
add_spherical(T)[source]

Add the spherical (hydrostatic) part of a tensor.

Parameters
Tstr

Name of tensor dataset.

Examples

Add the hydrostatic part of the Cauchy stress ‘sigma’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.add_spherical('sigma')
add_strain(F='F', t='V', m=0.0)[source]

Add strain tensor of a deformation gradient.

For details, see damask.mechanics.strain.

Parameters
Fstr, optional

Name of deformation gradient dataset. Defaults to ‘F’.

t{‘V’, ‘U’}, optional

Type of the polar decomposition, ‘V’ for left stretch tensor and ‘U’ for right stretch tensor. Defaults to ‘V’.

mfloat, optional

Order of the strain calculation. Defaults to 0.0.

Examples

Add the Biot strain based on the deformation gradient ‘F’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.strain(t='U',m=0.5)

Add the plastic Euler-Almansi strain based on the plastic deformation gradient ‘F_p’:

>>> import damask
>>> r = damask.Result('my_file.hdf5')
>>> r.strain('F_p','V',-1)
add_stretch_tensor(F='F', t='V')[source]

Add stretch tensor of a deformation gradient.

Parameters
Fstr, optional

Name of deformation gradient dataset. Defaults to ‘F’.

t{‘V’, ‘U’}, optional

Type of the polar decomposition, ‘V’ for left stretch tensor and ‘U’ for right stretch tensor. Defaults to ‘V’.

add_curl(f)[source]

Add curl of a field.

Parameters
fstr

Name of vector or tensor field dataset.

Notes

This function is only available for structured grids, i.e. results from the grid solver.

add_divergence(f)[source]

Add divergence of a field.

Parameters
fstr

Name of vector or tensor field dataset.

Notes

This function is only available for structured grids, i.e. results from the grid solver.

add_gradient(f)[source]

Add gradient of a field.

Parameters
fstr

Name of scalar or vector field dataset.

Notes

This function is only available for structured grids, i.e. results from the grid solver.

export_XDMF(output='*')[source]

Write XDMF file to directly visualize data in DADF5 file.

The XDMF format is only supported for structured grids with single phase and single constituent. For other cases use export_VTK.

Parameters
output(list of) str

Names of the datasets included in the XDMF file. Defaults to ‘*’, in which case all datasets are considered.

export_VTK(output='*', mode='cell', constituents=None, fill_float=nan, fill_int=0, parallel=True)[source]

Export to VTK cell/point data.

One VTK file per visible increment is created. For point data, the VTK format is poly data (.vtp). For cell data, either an image (.vti) or unstructured (.vtu) dataset is written for grid-based or mesh-based simulations, respectively.

Parameters
output(list of) str, optional

Names of the datasets to export to the VTK file. Defaults to ‘*’, in which case all datasets are exported.

mode{‘cell’, ‘point’}

Export in cell format or point format. Defaults to ‘cell’.

constituents(list of) int, optional

Constituents to consider. Defaults to None, in which case all constituents are considered.

fill_floatfloat

Fill value for non-existent entries of floating point type. Defaults to NaN.

fill_intint

Fill value for non-existent entries of integer type. Defaults to 0.

parallelbool

Write VTK files in parallel in a separate background process. Defaults to True.

get(output='*', flatten=True, prune=True)[source]

Collect data per phase/homogenization reflecting the group/folder structure in the DADF5 file.

Parameters
output(list of) str

Names of the datasets to read. Defaults to ‘*’, in which case all datasets are read.

flattenbool

Remove singular levels of the folder hierarchy. This might be beneficial in case of single increment, phase/homogenization, or field. Defaults to True.

prunebool

Remove branches with no data. Defaults to True.

Returns
datadict of numpy.ndarray

Datasets structured by phase/homogenization and according to selected view.

place(output='*', flatten=True, prune=True, constituents=None, fill_float=nan, fill_int=0)[source]

Merge data into spatial order that is compatible with the damask.VTK geometry representation.

The returned data structure reflects the group/folder structure in the DADF5 file.

Multi-phase data is fused into a single output. place is equivalent to get if only one phase/homogenization and one constituent is present.

Parameters
output(list of) str, optional

Names of the datasets to read. Defaults to ‘*’, in which case all datasets are placed.

flattenbool

Remove singular levels of the folder hierarchy. This might be beneficial in case of single increment or field. Defaults to True.

prunebool

Remove branches with no data. Defaults to True.

constituents(list of) int, optional

Constituents to consider. Defaults to None, in which case all constituents are considered.

fill_floatfloat

Fill value for non-existent entries of floating point type. Defaults to NaN.

fill_intint

Fill value for non-existent entries of integer type. Defaults to 0.

Returns
datadict of numpy.ma.MaskedArray

Datasets structured by spatial position and according to selected view.

save_VTK(output='*', mode='cell', constituents=None, fill_float=nan, fill_int=0, parallel=True)

Export to VTK cell/point data.

One VTK file per visible increment is created. For point data, the VTK format is poly data (.vtp). For cell data, either an image (.vti) or unstructured (.vtu) dataset is written for grid-based or mesh-based simulations, respectively.

Parameters
output(list of) str, optional

Names of the datasets to export to the VTK file. Defaults to ‘*’, in which case all datasets are exported.

mode{‘cell’, ‘point’}

Export in cell format or point format. Defaults to ‘cell’.

constituents(list of) int, optional

Constituents to consider. Defaults to None, in which case all constituents are considered.

fill_floatfloat

Fill value for non-existent entries of floating point type. Defaults to NaN.

fill_intint

Fill value for non-existent entries of integer type. Defaults to 0.

parallelbool

Write VTK files in parallel in a separate background process. Defaults to True.

save_XDMF(output='*')

Write XDMF file to directly visualize data in DADF5 file.

The XDMF format is only supported for structured grids with single phase and single constituent. For other cases use export_VTK.

Parameters
output(list of) str

Names of the datasets included in the XDMF file. Defaults to ‘*’, in which case all datasets are considered.


Colormap

class damask.Colormap(colors, name='from_list', N=None)[source]

Enhance matplotlib colormap functionality to be used within DAMASK.

References

K. Moreland, Proceedings of the 5th International Symposium on Advances in Visual Computing, 2009 https://doi.org/10.1007/978-3-642-10520-3_9

P. Eisenlohr et al., International Journal of Plasticity 46:37–53, 2013 https://doi.org/10.1016/j.ijplas.2012.09.012

Matplotlib colormaps overview https://matplotlib.org/tutorials/colors/colormaps.html

Methods

__call__(X[, alpha, bytes])

Parameters

copy()

Return a copy of the colormap.

from_predefined(name[, N])

Select from a set of predefined colormaps.

from_range(low, high[, name, N, model])

Create a perceptually uniform colormap between given (inclusive) bounds.

get_bad()

Get the color for masked values.

get_over()

Get the color for high out-of-range values.

get_under()

Get the color for low out-of-range values.

is_gray()

Return whether the colormap is grayscale.

reversed([name])

Reverse.

save_ASCII([fname])

Save as ASCII file.

save_GOM([fname])

Save as ASCII file for use in GOM Aramis.

save_gmsh([fname])

Save as ASCII file for use in gmsh.

save_paraview([fname])

Save as JSON file for use in Paraview.

set_bad([color, alpha])

Set the color for masked values.

set_extremes(*[, bad, under, over])

Set the colors for masked (bad) values and, when norm.clip = False, low (under) and high (over) out-of-range values.

set_over([color, alpha])

Set the color for high out-of-range values.

set_under([color, alpha])

Set the color for low out-of-range values.

shade(field[, bounds, gap])

Generate PIL image of 2D field using colormap.

with_extremes(*[, bad, under, over])

Return a copy of the colormap, for which the colors for masked (bad) values and, when norm.clip = False, low (under) and high (over) out-of-range values, have been set accordingly.

static from_range(low, high, name='DAMASK colormap', N=256, model='rgb')[source]

Create a perceptually uniform colormap between given (inclusive) bounds.

Colors are internally stored as R(ed) G(green) B(lue) values. The colormap can be used in matplotlib/seaborn or exported to file for external use.

Parameters
lownumpy.ndarray of shape (3)

Color definition for minimum value.

highnumpy.ndarray of shape (3)

Color definition for maximum value.

Nint, optional

The number of color quantization levels. Defaults to 256.

namestr, optional

The name of the colormap. Defaults to DAMASK colormap.

model{‘rgb’, ‘hsv’, ‘hsl’, ‘xyz’, ‘lab’, ‘msh’}

Colormodel used for input color definitions. Defaults to rgb. The available color models are: - ‘rgb’: R(ed) G(green) B(lue). - ‘hsv’: H(ue) S(aturation) V(alue). - ‘hsl’: H(ue) S(aturation) L(uminance). - ‘xyz’: CIE Xyz. - ‘lab’: CIE Lab. - ‘msh’: Msh (for perceptual uniform interpolation).

Returns
newdamask.Colormap

Colormap within given bounds.

Examples

>>> import damask
>>> damask.Colormap.from_range((0,0,1),(0,0,0),'blue_to_black')
static from_predefined(name, N=256)[source]

Select from a set of predefined colormaps.

Predefined colormaps include native matplotlib colormaps and common DAMASK colormaps.

Parameters
namestr

The name of the colormap.

Nint, optional

The number of color quantization levels. Defaults to 256. This parameter is not used for matplotlib colormaps that are of type ListedColormap.

Returns
newdamask.Colormap

Predefined colormap.

Examples

>>> import damask
>>> damask.Colormap.from_predefined('strain')
shade(field, bounds=None, gap=None)[source]

Generate PIL image of 2D field using colormap.

Parameters
fieldnumpy.array of shape (:,:)

Data to be shaded.

boundsiterable of len (2), optional

Colormap value range (low,high).

gapfield.dtype, optional

Transparent value. NaN will always be rendered transparent.

Returns
PIL.Image

RGBA image of shaded data.

reversed(name=None)[source]

Reverse.

Parameters
namestr, optional

The name for the reversed colormap. A name of None will be replaced by the name of the parent colormap + “_r”.

Returns
damask.Colormap

The reversed colormap.

Examples

>>> import damask
>>> damask.Colormap.from_predefined('stress').reversed()
save_paraview(fname=None)[source]

Save as JSON file for use in Paraview.

Parameters
fnamefile, str, or pathlib.Path, optional

Filename to store results. If not given, the filename will consist of the name of the colormap and extension ‘.json’.

save_ASCII(fname=None)[source]

Save as ASCII file.

Parameters
fnamefile, str, or pathlib.Path, optional

Filename to store results. If not given, the filename will consist of the name of the colormap and extension ‘.txt’.

save_GOM(fname=None)[source]

Save as ASCII file for use in GOM Aramis.

Parameters
fnamefile, str, or pathlib.Path, optional

Filename to store results. If not given, the filename will consist of the name of the colormap and extension ‘.legend’.

save_gmsh(fname=None)[source]

Save as ASCII file for use in gmsh.

Parameters
fnamefile, str, or pathlib.Path, optional

Filename to store results. If not given, the filename will consist of the name of the colormap and extension ‘.msh’.

copy()[source]

Return a copy of the colormap.

get_bad()[source]

Get the color for masked values.

get_over()[source]

Get the color for high out-of-range values.

get_under()[source]

Get the color for low out-of-range values.

is_gray()[source]

Return whether the colormap is grayscale.

set_bad(color='k', alpha=None)[source]

Set the color for masked values.

set_extremes(*, bad=None, under=None, over=None)[source]

Set the colors for masked (bad) values and, when norm.clip = False, low (under) and high (over) out-of-range values.

set_over(color='k', alpha=None)[source]

Set the color for high out-of-range values.

set_under(color='k', alpha=None)[source]

Set the color for low out-of-range values.

with_extremes(*, bad=None, under=None, over=None)[source]

Return a copy of the colormap, for which the colors for masked (bad) values and, when norm.clip = False, low (under) and high (over) out-of-range values, have been set accordingly.

colorbar_extend

When this colormap exists on a scalar mappable and colorbar_extend is not False, colorbar creation will pick up colorbar_extend as the default value for the extend keyword in the matplotlib.colorbar.Colorbar constructor.


Orientation

class damask.Orientation(rotation=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]

Representation of crystallographic orientation as combination of rotation and either crystal family or Bravais lattice.

The crystal family is one of:

  • triclinic

  • monoclinic

  • orthorhombic

  • tetragonal

  • hexagonal

  • cubic

and enables symmetry-related operations such as “equivalent”, “reduced”, “disorientation”, “IPF_color”, or “to_SST”.

The Bravais lattice is given in the Pearson notation:

  • triclinic
    • aP : primitive

  • monoclininic
    • mP : primitive

    • mS : base-centered

  • orthorhombic
    • oP : primitive

    • oS : base-centered

    • oI : body-centered

    • oF : face-centered

  • tetragonal
    • tP : primitive

    • tI : body-centered

  • hexagonal
    • hP : primitive

  • cubic
    • cP : primitive

    • cI : body-centered

    • cF : face-centered

and inherits the corresponding crystal family. Specifying a Bravais lattice, compared to just the crystal family, extends the functionality of Orientation objects to include operations such as “Schmid”, “related”, or “to_pole” that require a lattice type and its parameters.

Examples

An array of 3 x 5 random orientations reduced to the fundamental zone of tetragonal symmetry:

>>> damask.Orientation.from_random(shape=(3,5),lattice='tetragonal').reduced
Attributes
basis_real

Calculate orthogonal real space crystal basis.

basis_reciprocal

Calculate reciprocal (dual) crystal basis.

equivalent

Orientations that are symmetrically equivalent.

immutable

Return immutable parameters of own lattice.

in_FZ

Check whether orientation falls into fundamental zone of own symmetry.

in_disorientation_FZ

Check whether orientation falls into fundamental zone of disorientations.

parameters

Return lattice parameters a, b, c, alpha, beta, gamma.

quaternion
ratio

Return axes ratios of own lattice.

reduced

Select symmetrically equivalent orientation that falls into fundamental zone according to symmetry.

shape
size
symmetry_operations

Symmetry operations as Rotations.

Methods

IPF_color(vector[, in_SST, proper])

Map vector to RGB color within standard stereographic triangle of own symmetry.

Schmid(mode)

Calculate Schmid matrix P = d ⨂ n in the lab frame for given lattice shear kinematics.

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

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

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 separated 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, return_cloud])

Return orientation average over last dimension.

broadcast_to(shape[, mode])

Broadcast array.

copy(**kwargs)

Create deep copy.

disorientation(other[, return_operators])

Calculate disorientation between myself and given other orientation.

flatten([order])

Flatten array.

from_Euler_angles(**kwargs)

Initialize from Bunge-Euler angles.

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

Sample discrete values from a binned ODF.

from_Rodrigues_vector(**kwargs)

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

from_axis_angle(**kwargs)

Initialize from Axis angle pair.

from_basis(**kwargs)

Initialize from lattice basis vectors.

from_cubochoric(**kwargs)

Initialize from cubochoric vector.

from_directions(uvw, hkl, **kwargs)

Initialize orientation object from two crystallographic directions.

from_fiber_component(**kwargs)

Calculate set of rotations with Gaussian distribution around direction.

from_homochoric(**kwargs)

Initialize from homochoric vector.

from_matrix(**kwargs)

Initialize from rotation matrix.

from_parallel(a, b, **kwargs)

Initialize from pairs of two orthogonal lattice basis vectors.

from_quaternion(**kwargs)

Initialize from quaternion.

from_random(**kwargs)

Draw random rotation.

from_spherical_component(**kwargs)

Calculate set of rotations with Gaussian distribution around center.

in_SST(vector[, proper])

Check whether given crystal frame vector falls into standard stereographic triangle of own symmetry.

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

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

misorientation(other)

Calculate misorientation to other Rotation.

related(model)

Orientations derived from the given relationship.

relation_operations(model[, return_lattice])

Crystallographic orientation relationships for phase transformations.

reshape(shape[, order])

Reshape array.

to_SST(vector[, proper, return_operators])

Rotate vector to ensure it falls into (improper or proper) standard stereographic triangle of crystal symmetry.

to_frame(*[, uvw, hkl, with_symmetry])

Calculate crystal frame vector along lattice direction [uvw] or plane normal (hkl).

to_lattice(*[, direction, plane])

Calculate lattice vector corresponding to crystal frame direction or plane normal.

to_pole(*[, uvw, hkl, with_symmetry])

Calculate lab frame vector along lattice direction [uvw] or plane normal (hkl).

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 Orientation.

Parameters
otherOrientation

Orientation 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 orientations 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 Orientation.

Parameters
otherOrientation

Orientation 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 orientations.

classmethod from_random(**kwargs)[source]

Draw 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. If None, then fresh, unpredictable entropy will be pulled from the OS.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_quaternion(**kwargs)[source]

Initialize from quaternion.

Parameters
qnumpy.ndarray of shape (…,4)

Unit quaternion in positive real hemisphere: (q_0, q_1, q_2, q_3), ǀ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

Convention used. Defaults to -1.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_Euler_angles(**kwargs)[source]

Initialize from Bunge-Euler angles.

Parameters
phinumpy.ndarray of shape (…,3)

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

degreesboolean, optional

Bunge-Euler angles are given in degrees. Defaults to False.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_axis_angle(**kwargs)[source]

Initialize from Axis angle pair.

Parameters
axis_anglenumpy.ndarray of shape (…,4)

Axis angle pair: [n_1, n_2, n_3, ω], ǀnǀ = 1 and ω ∈ [0,π] unless degrees = True: ω ∈ [0,180].

degreesboolean, optional

Angle ω is given in degrees. Defaults to False.

normalize: boolean, optional

Allow ǀnǀ ≠ 1. Defaults to False.

Pint ∈ {-1,1}, optional

Convention used. Defaults to -1.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_basis(**kwargs)[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.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_matrix(**kwargs)[source]

Initialize from rotation matrix.

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

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

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_Rodrigues_vector(**kwargs)[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)), ǀnǀ = 1 and ω ∈ [0,π].

normalizeboolean, optional

Allow ǀnǀ ≠ 1. Defaults to False.

Pint ∈ {-1,1}, optional

Convention used. Defaults to -1.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_homochoric(**kwargs)[source]

Initialize from homochoric vector.

Parameters
hnumpy.ndarray of shape (…,3)

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

Pint ∈ {-1,1}, optional

Convention used. Defaults to -1.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_cubochoric(**kwargs)[source]

Initialize from cubochoric vector.

Parameters
xnumpy.ndarray of shape (…,3)

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

Pint ∈ {-1,1}, optional

Convention used. Defaults to -1.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_spherical_component(**kwargs)[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.

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.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_fiber_component(**kwargs)[source]

Calculate set of rotations with Gaussian distribution around direction.

Parameters
alphanumpy.ndarray of size 2

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

betanumpy.ndarray of size 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.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

classmethod from_directions(uvw, hkl, **kwargs)[source]

Initialize orientation object from two crystallographic directions.

Parameters
uvwlist, numpy.ndarray of shape (…,3)

lattice direction aligned with lab frame x-direction.

hkllist, numpy.ndarray of shape (…,3)

lattice plane normal aligned with lab frame z-direction.

latticestr

Either a crystal family out of {triclinic, monoclinic, orthorhombic, tetragonal, hexagonal, cubic} or a Bravais lattice out of {aP, mP, mS, oP, oS, oI, oF, tP, tI, hP, cP, cI, cF}. When specifying a Bravais lattice, additional lattice parameters might be required.

afloat, optional

Length of lattice parameter ‘a’.

bfloat, optional

Length of lattice parameter ‘b’.

cfloat, optional

Length of lattice parameter ‘c’.

alphafloat, optional

Angle between b and c lattice basis.

betafloat, optional

Angle between c and a lattice basis.

gammafloat, optional

Angle between a and b lattice basis.

degreesbool, optional

Angles are given in degrees. Defaults to False.

property symmetry_operations

Symmetry operations as Rotations.

property equivalent

Orientations that are symmetrically equivalent.

One dimension (length corresponds to number of symmetrically equivalent orientations) is added to the left of the Rotation array.

property reduced

Select symmetrically equivalent orientation that falls into fundamental zone according to symmetry.

property in_FZ

Check whether orientation falls into fundamental zone of own symmetry.

Returns
innumpy.ndarray of quaternion.shape

Boolean array indicating whether Rodrigues-Frank vector falls into fundamental zone.

Notes

Fundamental zones in Rodrigues space are point-symmetric around origin.

References

A. Heinz and P. Neumann, Acta Crystallographica Section A 47:780-789, 1991 https://doi.org/10.1107/S0108767391006864

property in_disorientation_FZ

Check whether orientation falls into fundamental zone of disorientations.

Returns
innumpy.ndarray of quaternion.shape

Boolean array indicating whether Rodrigues-Frank vector falls into disorientation FZ.

References

A. Heinz and P. Neumann, Acta Crystallographica Section A 47:780-789, 1991 https://doi.org/10.1107/S0108767391006864

relation_operations(model, return_lattice=False)[source]

Crystallographic orientation relationships for phase transformations.

Parameters
modelstr

Name of orientation relationship.

return_latticebool, optional

Return the target lattice in addition.

Returns
operationsRotations

Rotations characterizing the orientation relationship.

References

S. Morito et al., Journal of Alloys and Compounds 577:s587-s592, 2013 https://doi.org/10.1016/j.jallcom.2012.02.004

K. Kitahara et al., Acta Materialia 54(5):1279-1288, 2006 https://doi.org/10.1016/j.actamat.2005.11.001

Y. He et al., Journal of Applied Crystallography 39:72-81, 2006 https://doi.org/10.1107/S0021889805038276

H. Kitahara et al., Materials Characterization 54(4-5):378-386, 2005 https://doi.org/10.1016/j.matchar.2004.12.015

Y. He et al., Acta Materialia 53(4):1179-1190, 2005 https://doi.org/10.1016/j.actamat.2004.11.021

related(model)[source]

Orientations derived from the given relationship.

One dimension (length according to number of related orientations) is added to the left of the Rotation array.

append(other)

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

Parameters
otherdamask.Rotation
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.

as_Euler_angles(degrees=False)

Represent as Bunge-Euler angles.

Parameters
degreesbool, optional

Return angles in degrees.

Returns
phinumpy.ndarray of shape (…,3)

Bunge-Euler angles: (φ_1, ϕ, φ_2), φ_1 ∈ [0,2π], ϕ ∈ [0,π], φ_2 ∈ [0,2π] unless degrees == True: φ_1 ∈ [0,360], ϕ ∈ [0,180], φ_2 ∈ [0,360]

as_Rodrigues_vector(compact=False)

Represent as Rodrigues-Frank vector with separated axis and angle argument.

Parameters
compactbool, optional

Return as actual Rodrigues-Frank vector, i.e. axis and angle argument are not separated.

Returns
rhonumpy.ndarray of shape (…,4) containing

[n_1, n_2, n_3, tan(ω/2)], ǀnǀ = 1 and ω ∈ [0,π] unless compact == True: numpy.ndarray of shape (…,3) containing tan(ω/2) [n_1, n_2, n_3], ω ∈ [0,π].

as_axis_angle(degrees=False, pair=False)

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) unless pair == True:

tuple containing numpy.ndarray of shapes (…,3) and (…) Axis angle pair: (n_1, n_2, n_3, ω), ǀnǀ = 1 and ω ∈ [0,π] unless degrees = True: ω ∈ [0,180].

as_cubochoric()

Represent as cubochoric vector.

Returns
xnumpy.ndarray of shape (…,3)

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

as_homochoric()

Represent as homochoric vector.

Returns
hnumpy.ndarray of shape (…,3)

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

as_matrix()

Represent as rotation matrix.

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

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

as_quaternion()

Represent as unit quaternion.

Returns
qnumpy.ndarray of shape (…,4)

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

broadcast_to(shape, mode='right')

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.

flatten(order='C')

Flatten array.

Returns
flatteneddamask.Rotation

Rotation flattened to single dimension.

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

Sample discrete values from a binned ODF.

Parameters
weightsnumpy.ndarray of shape (n)

Texture intensity values (probability density or volume fraction) at Euler 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 grid values are in degrees. Defaults to True.

fractionsboolean, optional

ODF values correspond to volume fractions, not probability density. 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_parallel(a, b, **kwargs)

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.

misorientation(other)

Calculate misorientation to other Rotation.

Parameters
otherdamask.Rotation

Rotation to which the misorientation is computed.

Returns
gdamask.Rotation

Misorientation.

property parameters

Return lattice parameters a, b, c, alpha, beta, gamma.

reshape(shape, order='C')

Reshape array.

Returns
reshapeddamask.Rotation

Rotation of given shape.

property immutable

Return immutable parameters of own lattice.

property ratio

Return axes ratios of own lattice.

property basis_real

Calculate orthogonal real space crystal basis.

References

C.T. Young and J.L. Lytton, Journal of Applied Physics 43:1408–1417, 1972 https://doi.org/10.1063/1.1661333

property basis_reciprocal

Calculate reciprocal (dual) crystal basis.

in_SST(vector, proper=False)[source]

Check whether given crystal frame vector falls into standard stereographic triangle of own symmetry.

Parameters
vectornumpy.ndarray of shape (…,3)

Vector to check.

properbool, optional

Consider only vectors with z >= 0, hence combine two neighboring SSTs. Defaults to False.

Returns
innumpy.ndarray of shape (…)

Boolean array indicating whether vector falls into SST.

References

Bases are computed from

>>> basis = {
...    'cubic' :       np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,1.]/np.sqrt(2.),                #              green
...                                            [1.,1.,1.]/np.sqrt(3.)]).T),           #              blue
...    'hexagonal' :   np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,0.],                            #              green
...                                            [np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T),  #              blue
...    'tetragonal' :  np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,0.],                            #              green
...                                            [1.,1.,0.]/np.sqrt(2.)]).T),           #              blue
...    'orthorhombic': np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,0.],                            #              green
...                                            [0.,1.,0.]]).T),                       #              blue
...    }
IPF_color(vector, in_SST=True, proper=False)[source]

Map vector to RGB color within standard stereographic triangle of own symmetry.

Parameters
vectornumpy.ndarray of shape (…,3)

Vector to colorize.

in_SSTbool, optional

Consider symmetrically equivalent orientations such that poles are located in SST. Defaults to True.

properbool, optional

Consider only vectors with z >= 0, hence combine two neighboring SSTs (with mirrored colors). Defaults to False.

Returns
rgbnumpy.ndarray of shape (…,3)

RGB array of IPF colors.

References

Bases are computed from

>>> basis = {
...    'cubic' :       np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,1.]/np.sqrt(2.),                #              green
...                                            [1.,1.,1.]/np.sqrt(3.)]).T),           #              blue
...    'hexagonal' :   np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,0.],                            #              green
...                                            [np.sqrt(3.),1.,0.]/np.sqrt(4.)]).T),  #              blue
...    'tetragonal' :  np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,0.],                            #              green
...                                            [1.,1.,0.]/np.sqrt(2.)]).T),           #              blue
...    'orthorhombic': np.linalg.inv(np.array([[0.,0.,1.],                            # direction of red
...                                            [1.,0.,0.],                            #              green
...                                            [0.,1.,0.]]).T),                       #              blue
...    }

Examples

Inverse pole figure color of the e_3 direction for a crystal in “Cube” orientation with cubic symmetry:

>>> o = damask.Orientation(lattice='cubic')
>>> o.IPF_color([0,0,1])
array([1., 0., 0.])
disorientation(other, return_operators=False)[source]

Calculate disorientation between myself and given other orientation.

Parameters
otherOrientation

Orientation to calculate disorientation for. Shape of other blends with shape of own rotation array. For example, shapes of (2,3) for own rotations and (3,2) for other’s result in (2,3,2) disorientations.

return_operatorsbool, optional

Return index pair of symmetrically equivalent orientations that result in disorientation axis falling into FZ. Defaults to False.

Returns
disorientationOrientation

Disorientation between self and other.

operatorsnumpy.ndarray int of shape (…,2), conditional

Index of symmetrically equivalent orientation that rotated vector to the SST.

Notes

Currently requires same crystal family for both orientations. For extension to cases with differing symmetry see A. Heinz and P. Neumann 1991 and 10.1107/S0021889808016373.

Examples

Disorientation between two specific orientations of hexagonal symmetry:

>>> import damask
>>> a = damask.Orientation.from_Eulers(phi=[123,32,21],degrees=True,lattice='hexagonal')
>>> b = damask.Orientation.from_Eulers(phi=[104,11,87],degrees=True,lattice='hexagonal')
>>> a.disorientation(b)
Crystal family hexagonal
Quaternion: (real=0.976, imag=<+0.189, +0.018, +0.103>)
Matrix:
[[ 0.97831006  0.20710935  0.00389135]
 [-0.19363288  0.90765544  0.37238141]
 [ 0.07359167 -0.36505797  0.92807163]]
Bunge Eulers / deg: (11.40, 21.86, 0.60)
average(weights=None, return_cloud=False)[source]

Return orientation average over last dimension.

Parameters
weightsnumpy.ndarray, optional

Relative weights of orientations.

return_cloudbool, optional

Return the set of symmetrically equivalent orientations that was used in averaging. Defaults to False.

Returns
averageOrientation

Weighted average of original Orientation field.

cloudOrientations, conditional

Set of symmetrically equivalent orientations that were used in averaging.

References

J.C. Glez and J. Driver, Journal of Applied Crystallography 34:280-288, 2001 https://doi.org/10.1107/S0021889801003077

to_SST(vector, proper=False, return_operators=False)[source]

Rotate vector to ensure it falls into (improper or proper) standard stereographic triangle of crystal symmetry.

Parameters
vectornumpy.ndarray of shape (…,3)

Lab frame vector to align with crystal frame direction. Shape of other blends with shape of own rotation array. For example, a rotation array of shape (3,2) and a (2,4) vector array result in (3,2,4) outputs.

properbool, optional

Consider only vectors with z >= 0, hence combine two neighboring SSTs. Defaults to False.

return_operatorsbool, optional

Return the symmetrically equivalent orientation that rotated vector to SST. Defaults to False.

Returns
vector_SSTnumpy.ndarray of shape (…,3)

Rotated vector falling into SST.

operatorsnumpy.ndarray int of shape (…), conditional

Index of symmetrically equivalent orientation that rotated vector to SST.

to_lattice(*, direction=None, plane=None)[source]

Calculate lattice vector corresponding to crystal frame direction or plane normal.

Parameters
direction|normalnumpy.ndarray of shape (…,3)

Vector along direction or plane normal.

Returns
Millernumpy.ndarray of shape (…,3)

lattice vector of direction or plane. Use util.scale_to_coprime to convert to (integer) Miller indices.

to_frame(*, uvw=None, hkl=None, with_symmetry=False)[source]

Calculate crystal frame vector along lattice direction [uvw] or plane normal (hkl).

Parameters
uvw|hklnumpy.ndarray of shape (…,3)

Miller indices of crystallographic direction or plane normal.

with_symmetrybool, optional

Calculate all N symmetrically equivalent vectors.

Returns
vectornumpy.ndarray of shape (…,3) or (N,…,3)

Crystal frame vector (or vectors if with_symmetry) along [uvw] direction or (hkl) plane normal.

to_pole(*, uvw=None, hkl=None, with_symmetry=False)[source]

Calculate lab frame vector along lattice direction [uvw] or plane normal (hkl).

Parameters
uvw|hklnumpy.ndarray of shape (…,3)

Miller indices of crystallographic direction or plane normal.

with_symmetrybool, optional

Calculate all N symmetrically equivalent vectors.

Returns
vectornumpy.ndarray of shape (…,3) or (N,…,3)

Lab frame vector (or vectors if with_symmetry) along [uvw] direction or (hkl) plane normal.

Schmid(mode)[source]

Calculate Schmid matrix P = d ⨂ n in the lab frame for given lattice shear kinematics.

Parameters
mode{‘slip’, ‘twin’}

Type of kinematics.

Returns
Pnumpy.ndarray of shape (…,N,3,3)

Schmid matrix for each of the N deformation systems.

Examples

Schmid matrix (in lab frame) of slip systems of a face-centered cubic crystal in “Goss” orientation.

>>> import damask
>>> import numpy as np
>>> np.set_printoptions(3,suppress=True,floatmode='fixed')
>>> damask.Orientation.from_Eulers(phi=[0,45,0],degrees=True,lattice='cF').Schmid('slip')[0]
array([[ 0.000,  0.000,  0.000],
       [ 0.577, -0.000,  0.816],
       [ 0.000,  0.000,  0.000]])