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 (DAMASK HDF5) file.
A DADF5 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 \x1b[2mCreated by DAMASK_grid ... on ... executing "..."\x1b[0m ...
- Attributes:
coordinates0_node
Initial/undeformed nodal coordinates.
coordinates0_point
Initial/undeformed cell center coordinates.
- fields
geometry0
Initial/undeformed geometry.
- homogenizations
- increments
- phases
simulation_setup_files
Simulation setup files used to generate the Result object.
- times
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 the determinant of a tensor.
add_deviator
(T)Add the deviatoric part of a tensor.
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 a vector or tensor.
add_pole
([q, uvw, hkl, with_symmetry, normalize])Add lab frame vector along lattice direction [uvw] or plane normal (hkl).
add_rotation
(F)Add rotational part of a deformation gradient.
Add the spherical (hydrostatic) part of a tensor.
add_strain
([F, t, m])Add strain tensor (Seth-Hill family) of a deformation gradient.
add_stress_Cauchy
([P, F])Add Cauchy stress calculated from first Piola-Kirchhoff stress and deformation gradient.
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
()Return deepcopy(self).
export_DADF5
(fname[, output, mapping])Export visible components into a new DADF5 file.
export_DREAM3D
([q, target_dir])Export the visible components to DREAM3D compatible files.
export_VTK
([output, mode, constituents, ...])Export to VTK cell/point data.
export_XDMF
([output, target_dir, absolute_path])Write XDMF file to directly visualize data from DADF5 file.
export_simulation_setup
([output, ...])Export original simulation setup of the Result object.
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.
Collect information on all active datasets in the file.
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).
times_in_range
([start, end])Get times of all increments within a given time range.
view
(*[, increments, times, phases, ...])Set view.
view_all
()Make all attributes visible.
view_less
(*[, increments, times, phases, ...])Remove from view.
view_more
(*[, increments, times, phases, ...])Add to view.
enable_user_function
- copy()#
Return deepcopy(self).
Create deep copy.
- increments_in_range(start=None, end=None)[source]#
Get all increments within a given range.
- Parameters:
- startint or str, optional
Start increment. Defaults to first.
- endint or str, optional
End increment. Defaults to last.
- Returns:
- incrementslist of ints
Increment number of all increments within the given bounds.
- times_in_range(start=None, end=None)[source]#
Get times of all increments within a given time range.
- Parameters:
- startfloat, optional
Time of start increment. Defaults to time of first.
- endfloat, optional
Time of end increment. Defaults to time of last.
- Returns:
- timeslist of float
Time of each increment within the given bounds.
- view(*, increments=None, times=None, phases=None, homogenizations=None, fields=None, protected=None)[source]#
Set view.
Wildcard matching with ‘?’ and ‘*’ is supported. True is equivalent to ‘*’, False is equivalent to [].
- Parameters:
- increments: (list of) int, (list of) str, or bool, optional.
Numbers of increments to select.
- times: (list of) float, (list of) str, or bool, optional.
Simulation times of increments to select.
- phases: (list of) str, or bool, optional.
Names of phases to select.
- homogenizations: (list of) str, or bool, optional.
Names of homogenizations to select.
- fields: (list of) str, or bool, optional.
Names of fields to select.
- protected: bool, optional.
Protection status of existing data.
- 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(increments=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(*, increments=None, times=None, phases=None, homogenizations=None, fields=None)[source]#
Add to view.
Wildcard matching with ‘?’ and ‘*’ is supported. True is equivalent to ‘*’, False is equivalent to [].
- Parameters:
- increments: (list of) int, (list of) str, or bool, optional.
Numbers of increments to select.
- times: (list of) float, (list of) str, or bool, optional.
Simulation times of increments to select.
- phases: (list of) str, or bool, optional.
Names of phases to select.
- homogenizations: (list of) str, or bool, optional.
Names of homogenizations to select.
- fields: (list of) str, or bool, optional.
Names of fields to select.
- 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(*, increments=None, times=None, phases=None, homogenizations=None, fields=None)[source]#
Remove from view.
Wildcard matching with ‘?’ and ‘*’ is supported. True is equivalent to ‘*’, False is equivalent to [].
- Parameters:
- increments: (list of) int, (list of) str, or bool, optional.
Numbers of increments to select.
- times: (list of) float, (list of) str, or bool, optional.
Simulation times of increments to select.
- phases: (list of) str, or bool, optional.
Names of phases to select.
- homogenizations: (list of) str, or bool, optional.
Names of homogenizations to select.
- fields: (list of) str, or bool, optional.
Names of fields to select.
- 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)
- view_all()[source]#
Make all attributes visible.
- Returns:
- modified_viewdamask.Result
View with all attributes visible.
- 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.view(protected=False) \x1b[93m\x1b[1mWarning: Modification of existing datasets allowed!\x1b[0m\x1b[0m >>> 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.view(protected=False) \x1b[93m\x1b[1mWarning: Modification of existing datasets allowed!\x1b[0m\x1b[0m >>> r_unprotected.remove('F')
- list_data()[source]#
Collect information on all active datasets in the file.
- Returns:
- datalist of str
Line-formatted information about active datasets.
- property simulation_setup_files#
Simulation setup files used to generate the Result object.
- 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) Function equivalent_stress enabled in add_calculation. >>> 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.
Examples
Add the deviatoric part of Cauchy stress ‘sigma’:
>>> import damask >>> r = damask.Result('my_file.hdf5') >>> r.add_deviator('sigma') [...]
- add_eigenvalue(T_sym, eigenvalue='max')[source]#
Add eigenvalues of symmetric tensor.
- Parameters:
- T_symstr
Name of symmetric tensor dataset.
- eigenvalue{‘max’, ‘mid’, ‘min’}, optional
Eigenvalue. Defaults to ‘max’.
Examples
Add the minimum eigenvalue of Cauchy stress ‘sigma’:
>>> import damask >>> r = damask.Result('my_file.hdf5') >>> r.add_eigenvalue('sigma','min') [...]
- add_eigenvector(T_sym, eigenvalue='max')[source]#
Add eigenvector of symmetric tensor.
- Parameters:
- T_symstr
Name of symmetric tensor dataset.
- eigenvalue{‘max’, ‘mid’, ‘min’}, optional
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) or compatible
Lab frame direction for inverse pole figure.
- qstr, optional
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(l = [0,1,1], q = 'O') [...]
- 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 a 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 \(\vb{S} = \left(\vb{F}^{-1} \vb{P}\right)_\text{sym}\) follows the standard definition in nonlinear continuum mechanics. As such, no intermediate configuration, for instance that reached by \(\vb{F}_\text{p}\), is taken into account.
- add_pole(q='O', *, uvw=None, hkl=None, with_symmetry=False, normalize=True)[source]#
Add lab frame vector along lattice direction [uvw] or plane normal (hkl).
- Parameters:
- qstr, optional
Name of the dataset containing the crystallographic orientation as quaternions. Defaults to ‘O’.
- uvw|hklnumpy.ndarray of shape (3) or compatible
Miller indices of crystallographic direction or plane normal.
- with_symmetrybool, optional
Calculate all N symmetrically equivalent vectors. Defaults to False.
- normalizebool, optional
Normalize output vector. Defaults to True.
- 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 (Seth-Hill family) of a deformation gradient.
By default, the logarithmic strain based on the left stretch tensor is added.
- 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.
Notes
The presence of rotational parts in the elastic and plastic deformation gradient calls for the use of material/Lagragian strain measures (based on ‘U’) for plastic strains and spatial/Eulerian strain measures (based on ‘V’) for elastic strains when calculating averages.
The strain is defined as:
\[ \begin{align}\begin{aligned}\begin{split}m = 0 \\\\ \vb*{\epsilon}_V^{(0)} = \ln (\vb{V}) \\\\ \vb*{\epsilon}_U^{(0)} = \ln (\vb{U}) \\\\\end{split}\\\begin{split}m \neq 0 \\\\ \vb*{\epsilon}_V^{(m)} = \frac{1}{2m} (\vb{V}^{2m} - \vb{I}) \\\\ \vb*{\epsilon}_U^{(m)} = \frac{1}{2m} (\vb{U}^{2m} - \vb{I})\end{split}\end{aligned}\end{align} \]References
Examples
Add the Euler-Almansi strain:
>>> import damask >>> r = damask.Result('my_file.hdf5') >>> r.add_strain(t='V',m=-1.0) [...]
Add the plastic Biot strain:
>>> import damask >>> r = damask.Result('my_file.hdf5') >>> r.add_strain('F_p','U',0.5) [...]
- 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 implemented only for structured grids with one constituent and a single phase.
- add_divergence(f)[source]#
Add divergence of a field.
- Parameters:
- fstr
Name of vector or tensor field dataset.
Notes
This function is implemented only for structured grids with one constituent and a single phase.
- add_gradient(f)[source]#
Add gradient of a field.
- Parameters:
- fstr
Name of scalar or vector field dataset.
Notes
This function is implemented only for structured grids with one constituent and a single phase.
- 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, optional
Names of the datasets to read. Defaults to ‘*’, in which case all datasets are read.
- flattenbool, optional
Remove singular levels of the folder hierarchy. This might be beneficial in case of single increment, phase/homogenization, or field. Defaults to True.
- prunebool, optional
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 visible datasets are placed.
- flattenbool, optional
Remove singular levels of the folder hierarchy. This might be beneficial in case of single increment or field. Defaults to True.
- prunebool, optional
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, optional
Fill value for non-existent entries of floating point type. Defaults to NaN.
- fill_intint, optional
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.
- export_XDMF(output='*', target_dir=None, absolute_path=False)[source]#
Write XDMF file to directly visualize data from 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, optional
Names of the datasets included in the XDMF file. Defaults to ‘*’, in which case all datasets are considered.
- target_dirstr or pathlib.Path, optional
Directory to save XDMF file. Will be created if non-existent.
- absolute_pathbool, optional
Store absolute (instead of relative) path to DADF5 file. Defaults to False, i.e. the XDMF file expects the DADF5 file at a stable relative path.
Notes
This function is implemented only for structured grids with one constituent and a single phase.
- export_VTK(output='*', mode='cell', constituents=None, target_dir=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 PolyData (.vtp). For cell data, the file format is either ImageData (.vti) or UnstructuredGrid (.vtu) 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 visible datasets are exported.
- mode{‘cell’, ‘point’}, optional
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.
- target_dirstr or pathlib.Path, optional
Directory to save VTK files. Will be created if non-existent.
- fill_floatfloat, optional
Fill value for non-existent entries of floating point type. Defaults to NaN.
- fill_intint, optional
Fill value for non-existent entries of integer type. Defaults to 0.
- parallelbool, optional
Write VTK files in parallel in a separate background process. Defaults to True.
- export_DREAM3D(q='O', target_dir=None)[source]#
Export the visible components to DREAM3D compatible files.
One DREAM3D file per visible increment is created. The geometry is based on the undeformed configuration.
- Parameters:
- qstr, optional
Name of the dataset containing the crystallographic orientation as quaternions. Defaults to ‘O’.
- target_dirstr or pathlib.Path, optional
Directory to save DREAM3D files. Will be created if non-existent.
Notes
This function is implemented only for structured grids with one constituent.
- export_DADF5(fname, output='*', mapping=None)[source]#
Export visible components into a new DADF5 file.
A DADF5 (DAMASK HDF5) file contains DAMASK results. Its group/folder structure reflects the layout in material.yaml.
- Parameters:
- fnamestr or pathlib.Path
Name of the DADF5 file to be created.
- output(list of) str, optional
Names of the datasets to export. Defaults to ‘*’, in which case all visible datasets are exported.
- mappingnumpy.ndarray of int, shape (:,:,:), optional
Indices for regridding. Only applicable for grid solver results.
- export_simulation_setup(output='*', target_dir=None, overwrite=False)[source]#
Export original simulation setup of the Result object.
- Parameters:
- output(list of) str, optional
Names of the datasets to export to the file. Defaults to ‘*’, in which case all setup files are exported.
- target_dirstr or pathlib.Path, optional
Directory to save setup files. Will be created if non-existent.
- overwritebool, optional
Overwrite any existing setup files. Defaults to False.
Colormap#
- class damask.Colormap(colors, name)[source]#
Enhance matplotlib colormap functionality for use within DAMASK.
Colors are internally stored as R(ed) G(green) B(lue) values. A colormap can be used in matplotlib, seaborn, etc., or can be exported to file for external use.
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/stable/tutorials/colors/colormaps.html
Methods
__call__
(X[, alpha, bytes])- Parameters:
at
(fraction)Interpolate color at fraction.
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 the color for low out-of-range values.
is_gray
()Return whether the colormap is grayscale.
resampled
(lutsize)Return a new colormap with lutsize entries.
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.
- Parameters:
- lowsequence of float, len (3)
Color definition for minimum value.
- highsequence of float, len (3)
Color definition for maximum value.
- namestr, optional
Name of the colormap. Defaults to ‘DAMASK colormap’.
- Nint, optional
Number of color quantization levels. Defaults to 256.
- model{‘rgb’, ‘hsv’, ‘hsl’, ‘xyz’, ‘lab’, ‘msh’}
Color model used for input color definitions. Defaults to ‘rgb’. The available color models are: - ‘rgb’: Red Green Blue. - ‘hsv’: Hue Saturation Value. - ‘hsl’: Hue Saturation Luminance. - ‘xyz’: CIE Xyz. - ‘lab’: CIE Lab. - ‘msh’: Msh (for perceptually uniform interpolation).
- Returns:
- newdamask.Colormap
Colormap spanning given bounds.
Examples
>>> import damask >>> damask.Colormap.from_range((0,0,1),(0,0,0),'blue_to_black') Colormap: blue_to_black
- static from_predefined(name, N=256)[source]#
Select from a set of predefined colormaps.
Predefined colormaps (Colormap.predefined) include native matplotlib colormaps and common DAMASK colormaps.
- Parameters:
- namestr
Name of the colormap.
- Nint, optional
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') Colormap: strain
- at(fraction)[source]#
Interpolate color at fraction.
- Parameters:
- fraction(sequence of) float
Fractional coordinate(s) to evaluate Colormap at.
- Returns:
- colornumpy.ndarray, shape(…,4)
RGBA values of interpolated color(s).
Examples
>>> import damask >>> cmap = damask.Colormap.from_predefined('gray') >>> cmap.at(0.5) array([0.5, 0.5, 0.5, 1. ]) >>> 'rgb({},{},{})'.format(*cmap.at(0.5)) 'rgb(0.5,0.5,0.5)'
- shade(field, bounds=None, gap=None)[source]#
Generate PIL image of 2D field using colormap.
- Parameters:
- fieldnumpy.ndarray, shape (:,:)
Data to be shaded.
- boundssequence of float, len (2), optional
Value range (left,right) spanned by colormap.
- gapfield.dtype, optional
Transparent value. NaN will always be rendered transparent. Defaults to None.
- Returns:
- PIL.Image
RGBA image of shaded data.
- reversed(name=None)[source]#
Reverse.
- Parameters:
- namestr, optional
Name of the reversed colormap. Defaults to parent colormap name + ‘_r’.
- Returns:
- damask.Colormap
Reversed colormap.
Examples
>>> import damask >>> damask.Colormap.from_predefined('stress').reversed() Colormap: stress_r
- save_paraview(fname=None)[source]#
Save as JSON file for use in Paraview.
- Parameters:
- fnamefile, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + ‘.json’.
- save_ASCII(fname=None)[source]#
Save as ASCII file.
- Parameters:
- fnamefile, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + ‘.txt’.
- save_GOM(fname=None)[source]#
Save as ASCII file for use in GOM Aramis.
- Parameters:
- fnamefile, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + ‘.legend’.
- save_gmsh(fname=None)[source]#
Save as ASCII file for use in gmsh.
- Parameters:
- fnamefile, str, or pathlib.Path, optional
File to store results. Defaults to colormap name + ‘.msh’.
- 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.
- 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 theextend
keyword in the matplotlib.colorbar.Colorbar constructor.
Orientation#
- class damask.Orientation(rotation=array([1., 0., 0., 0.]), *, family=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
- monoclinic
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_frame” 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:
>>> import damask >>> o=damask.Orientation.from_random(shape=(3,5),family='tetragonal').reduced
- Attributes:
basis_real
Return orthogonal real space crystal basis.
basis_reciprocal
Return reciprocal (dual) crystal basis.
equivalent
Orientations that are symmetrically equivalent.
immutable
Return immutable lattice parameters.
in_FZ
Check whether orientation falls into fundamental zone of own symmetry.
in_disorientation_FZ
Check whether orientation falls into fundamental zone of disorientations.
lattice_points
Return lattice points.
orientation_relationships
Return labels of orientation relationships.
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
standard_triangle
Corners of the standard triangle.
symmetry_operations
Return symmetry operations.
Methods
IPF_color
(vector[, in_SST, proper])Map lab frame vector to RGB color within standard stereographic triangle of own symmetry.
Schmid
(*[, N_slip, N_twin])Calculate Schmid matrix P = d ⨂ n in the lab frame for selected deformation systems.
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)Return self@other.
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.
Represent as cubochoric vector.
Represent as homochoric vector.
Represent as rotation matrix.
Represent as unit quaternion.
average
([weights, return_cloud])Return orientation average over last dimension.
broadcast_to
(shape[, mode])Broadcast array.
copy
([rotation])Return deepcopy(self).
disorientation
(other[, return_operators])Calculate disorientation between self and given other orientation.
flatten
([order])Flatten array.
from_Euler_angles
(*, phi[, degrees, family, ...])Initialize from Bunge Euler angles.
from_ODF
(*, weights, phi[, shape, degrees, ...])Initialize with samples from a binned orientation distribution function (ODF).
from_Rodrigues_vector
(*, rho[, normalize, ...])Initialize from Rodrigues–Frank vector (with angle separated from axis).
from_axis_angle
(*, n_omega[, degrees, ...])Initialize from axis–angle pair.
from_basis
(*, basis[, orthonormal, ...])Initialize from basis vector triplet.
from_cubochoric
(*, x[, P, family, lattice, ...])Initialize from cubochoric vector.
from_directions
(uvw, hkl, **kwargs)Initialize orientation object from the crystallographic direction and plane parallel to lab x and z, respectively.
from_fiber_component
(*, crystal, sample[, ...])Initialize with samples from a Gaussian distribution around a given direction.
from_homochoric
(*, h[, P, family, lattice, ...])Initialize from homochoric vector.
from_matrix
(*, R[, normalize, family, ...])Initialize from rotation matrix.
from_parallel
(*, a, b[, family, lattice, c, ...])Initialize from pairs of two orthogonal basis vectors.
from_quaternion
(*, q[, accept_homomorph, ...])Initialize from quaternion.
from_random
(*[, shape, rng_seed, family, ...])Initialize with samples from a uniform distribution.
from_spherical_component
(*, center, sigma[, ...])Initialize with samples from a Gaussian distribution around a given 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.
kinematics
(mode)Return crystal kinematics systems.
misorientation
(other)Calculate misorientation to other Rotation.
related
(model[, target])All orientations related to self by given relationship model.
relation_operations
(model[, target])Crystallographic orientation relationships for phase transformations.
reshape
(newshape[, order])Reshape array.
to_SST
(vector[, proper, return_operators])Rotate lab frame vector to ensure it falls into (improper or proper) standard stereographic triangle of crystal symmetry.
to_frame
(*[, uvw, hkl, uvtw, hkil, ...])Calculate lab frame vector along lattice direction [uvw]/[uvtw] or plane normal (hkl)/(hkil).
to_lattice
(*[, direction, plane])Calculate lattice vector corresponding to lab frame direction or plane normal.
- copy(rotation=None)#
Return deepcopy(self).
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 of bool, shape (self.shape)
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_quaternion(*, q, accept_homomorph=False, normalize=False, P=-1, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from quaternion.
- Parameters:
- qnumpy.ndarray, shape (…,4)
Unit quaternion (q_0, q_1, q_2, q_3) in positive real hemisphere, i.e. ǀqǀ = 1 and q_0 ≥ 0.
- accept_homomorphbool, optional
Allow homomorphic variants, i.e. q_0 < 0 (negative real hemisphere). Defaults to False.
- normalize: bool, optional
Allow ǀqǀ ≠ 1. Defaults to False.
- Pint ∈ {-1,1}, optional
Sign convention. Defaults to -1.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
Examples
>>> import damask >>> damask.Rotation.from_quaternion([[1,0,0,0],[0,1,0,0]]) Quaternions of shape (2,) [[1. 0. 0. 0.] [0. 1. 0. 0.]]
- classmethod from_Euler_angles(*, phi, degrees=False, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None)[source]#
Initialize from Bunge Euler angles.
- Parameters:
- phinumpy.ndarray, shape (…,3)
Euler angles (φ_1 ∈ [0,2π], ϕ ∈ [0,π], φ_2 ∈ [0,2π]) or (φ_1 ∈ [0,360], ϕ ∈ [0,180], φ_2 ∈ [0,360]) if degrees == True.
- degreesbool, optional
Euler angles are given in degrees. Defaults to False.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
Notes
Bunge Euler angles correspond to a rotation axis sequence of z–x’–z’’.
Examples
>>> import damask >>> damask.Rotation.from_Euler_angles([180,0,0],degrees=True) Quaternion [0. 0. 0. 1.]
- classmethod from_axis_angle(*, n_omega, degrees=False, normalize=False, P=-1, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None)[source]#
Initialize from axis–angle pair.
- Parameters:
- n_omeganumpy.ndarray, shape (…,4)
Axis and angle (n_1, n_2, n_3, ω) with ǀnǀ = 1 and ω ∈ [0,π] or ω ∈ [0,180] if degrees == True.
- degreesbool, optional
Angle ω is given in degrees. Defaults to False.
- normalize: bool, optional
Allow ǀnǀ ≠ 1. Defaults to False.
- Pint ∈ {-1,1}, optional
Sign convention. Defaults to -1.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
Examples
>>> import damask >>> damask.Rotation.from_axis_angle([[0,0,1,90],[1,0,0,90]],degrees=True) Quaternions of shape (2,) [[0.707 0. 0. 0.707] [0.707 0.707 0. 0. ]]
- classmethod from_basis(*, basis, orthonormal=True, reciprocal=False, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from basis vector triplet.
- Parameters:
- basisnumpy.ndarray, shape (…,3,3)
Three three-dimensional basis vectors.
- orthonormalbool, optional
Basis is strictly orthonormal, i.e. is free of stretch components. Defaults to True.
- reciprocalbool, optional
Basis vectors are given in reciprocal (instead of real) space. Defaults to False.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
- classmethod from_matrix(*, R, normalize=False, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from rotation matrix.
- Parameters:
- Rnumpy.ndarray, shape (…,3,3)
Rotation matrix with det(R) = 1 and R.T ∙ R = I.
- normalizebool, optional
Rescales rotation matrix to unit determinant. Defaults to False.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
Examples
>>> import damask >>> damask.Rotation.from_matrix([[1,0,0],[0,0,-1],[0,1,0]]) Quaternion [ 0.707 -0.707 0. 0. ]
- classmethod from_parallel(*, a, b, family=None, lattice=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from pairs of two orthogonal basis vectors.
- Parameters:
- anumpy.ndarray, shape (…,2,3)
Two three-dimensional vectors of first orthogonal basis.
- bnumpy.ndarray, shape (…,2,3)
Corresponding three-dimensional vectors of second basis.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- cfloat, optional
Length of lattice parameter ‘c’.
- gammafloat, optional
Angle between a and b lattice basis.
- degreesbool, optional
Angles are given in degrees. Defaults to False.
- Returns:
- newdamask.Orientation
Examples
>>> import damask >>> damask.Rotation.from_parallel([[2,0,0],[0,1,0]],[[1,0,0],[0,2,0]]) Quaternion [1. 0. 0. 0.]
- classmethod from_Rodrigues_vector(*, rho, normalize=False, P=-1, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from Rodrigues–Frank vector (with angle separated from axis).
- Parameters:
- rhonumpy.ndarray, shape (…,4)
Rodrigues–Frank vector (n_1, n_2, n_3, tan(ω/2)) with ǀnǀ = 1 and ω ∈ [0,π].
- normalizebool, optional
Allow ǀnǀ ≠ 1. Defaults to False.
- Pint ∈ {-1,1}, optional
Sign convention. Defaults to -1.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
Examples
>>> import damask >>> damask.Rotation.from_Rodrigues_vector([0,0,1,1]) Quaternion [0.707 0. 0. 0.707]
- classmethod from_homochoric(*, h, P=-1, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from homochoric vector.
- Parameters:
- hnumpy.ndarray, 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.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
- classmethod from_cubochoric(*, x, P=-1, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize from cubochoric vector.
- Parameters:
- xnumpy.ndarray, 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.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
- classmethod from_random(*, shape=None, rng_seed=None, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Initialize with samples from a uniform distribution.
- Parameters:
- shape(sequence of) int, optional
Output shape. Defaults to None, which gives a scalar.
- 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.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
- classmethod from_ODF(*, weights, phi, shape=None, degrees=False, fractions=True, rng_seed=None, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None)[source]#
Initialize with samples from a binned orientation distribution function (ODF).
- Parameters:
- weightsnumpy.ndarray, shape (n)
Texture intensity values (probability density or volume fraction) at Euler space grid points.
- phinumpy.ndarray, shape (n,3)
Grid coordinates in Euler space at which weights are defined.
- shape(sequence of) int, optional
Output shape. Defaults to None, which gives a scalar.
- degreesbool, optional
Euler space grid coordinates are in degrees. Defaults to True.
- fractionsbool, 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.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
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 a cell-centered version, which avoids 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
- classmethod from_spherical_component(*, center, sigma, shape=None, degrees=False, rng_seed=None, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None)[source]#
Initialize with samples from a Gaussian distribution around a given center.
- Parameters:
- centerRotation or Orientation
Central rotation.
- sigmafloat
Standard deviation of (Gaussian) misorientation distribution.
- shape(sequence of) int, optional
Output shape. Defaults to None, which gives a scalar.
- degreesbool, optional
sigma is given in degrees. Defaults to False.
- 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.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
Examples
Create a brass texture consisting of 200 orientations:
>>> import damask >>> center = damask.Rotation.from_Euler_angles([35.,45.,0.],degrees=True) >>> brass = damask.Rotation.from_spherical_component(center=center,sigma=3.,shape=200,degrees=True)
Create a Goss texture consisting of 100 orientations:
>>> import damask >>> center = damask.Rotation.from_Euler_angles([0.,45.,0.],degrees=True) >>> goss = damask.Rotation.from_spherical_component(center=center,sigma=3.,shape=100,degrees=True)
- classmethod from_fiber_component(*, crystal, sample, sigma=0.0, shape=None, degrees=False, rng_seed=None, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None)[source]#
Initialize with samples from a Gaussian distribution around a given direction.
- Parameters:
- crystalnumpy.ndarray, shape (2)
Polar coordinates (polar angle θ from [0 0 1], azimuthal angle φ from [1 0 0]) of fiber direction in crystal frame.
- samplenumpy.ndarray, shape (2)
Polar coordinates (polar angle θ from z, azimuthal angle φ from x) of fiber direction in sample frame.
- sigmafloat, optional
Standard deviation of (Gaussian) misorientation distribution. Defaults to 0.
- shape(sequence of) int, optional
Output shape. Defaults to None, which gives a scalar.
- degreesbool, optional
sigma and polar coordinates are given in degrees. Defaults to False.
- 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.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
Notes
The crystal direction for (θ=0,φ=0) is [0 0 1], the sample direction for (θ=0,φ=0) is z.
Polar coordinates follow the ISO 80000-2:2019 convention typically used in physics. See https://en.wikipedia.org/wiki/Spherical_coordinate_system.
Ranges 0≤θ≤π and 0≤φ≤2π give a unique set of coordinates.
Examples
Create an ideal α-fiber texture (<1 1 0> ǀǀ RD=x) consisting of 200 orientations:
>>> import damask >>> import numpy as np >>> alpha = damask.Rotation.from_fiber_component([np.pi/4.,0.],[np.pi/2.,0.],shape=200)
Create an ideal γ-fiber texture (<1 1 1> ǀǀ ND=z) consisting of 100 orientations:
>>> import damask >>> gamma = damask.Rotation.from_fiber_component([54.7,45.],[0.,0.],shape=100,degrees=True)
- classmethod from_directions(uvw, hkl, **kwargs)[source]#
Initialize orientation object from the crystallographic direction and plane parallel to lab x and z, respectively.
- Parameters:
- uvwnumpy.ndarray, shape (…,3)
Lattice direction aligned with lab frame x-direction.
- hklnumpy.ndarray, shape (…,3)
Lattice plane normal aligned with lab frame z-direction.
- family{‘triclinic’, ‘monoclinic’, ‘orthorhombic’, ‘tetragonal’, ‘hexagonal’, ‘cubic’}, optional.
Name of the crystal family. Will be inferred if ‘lattice’ is given.
- lattice{‘aP’, ‘mP’, ‘mS’, ‘oP’, ‘oS’, ‘oI’, ‘oF’, ‘tP’, ‘tI’, ‘hP’, ‘cP’, ‘cI’, ‘cF’}, optional.
Name of the Bravais lattice in Pearson notation.
- 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.
- Returns:
- newdamask.Orientation
- 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 bool, shape (self.shape)
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 bool, shape (self.shape)
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
- disorientation(other, return_operators=False)[source]#
Calculate disorientation between self 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 of int, shape (…,2), conditional
Index of symmetrically equivalent orientation that rotated vector to the SST.
Notes
Requires same crystal family for both orientations.
Examples
Disorientation between two specific orientations of hexagonal symmetry:
>>> import damask >>> a = damask.Orientation.from_Euler_angles(phi=[123,32,21],degrees=True,family='hexagonal') >>> b = damask.Orientation.from_Euler_angles(phi=[104,11,87],degrees=True,family='hexagonal') >>> a.disorientation(b) Crystal family: hexagonal Quaternion [0.976 0.189 0.018 0.103]
Plot a sample from the Mackenzie distribution.
>>> import matplotlib.pyplot as plt >>> import damask >>> N = 10000 >>> a = damask.Orientation.from_random(shape=N,family='cubic') >>> b = damask.Orientation.from_random(shape=N,family='cubic') >>> n,omega = a.disorientation(b).as_axis_angle(degrees=True,pair=True) >>> plt.hist(omega,25) [...] >>> plt.show()
- average(weights=None, return_cloud=False)[source]#
Return orientation average over last dimension.
- Parameters:
- weightsnumpy.ndarray, shape (self.shape), optional
Relative weights of orientations. Defaults to equal weights.
- return_cloudbool, optional
Return the specific (symmetrically equivalent) orientations that were averaged. Defaults to False.
- Returns:
- averageOrientation
Weighted average of original Orientation field.
- cloudOrientations, conditional
Symmetrically equivalent version of each orientation that were actually 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 lab frame vector to ensure it falls into (improper or proper) standard stereographic triangle of crystal symmetry.
- Parameters:
- vectornumpy.ndarray, shape (…,3)
Lab frame vector to align with an SST crystal frame direction. Shape of vector blends with shape of own rotation array. For example, a rotation array of shape (3,2) and a vector array of shape (2,4) 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, shape (…,3)
Rotated vector falling into SST.
- operatornumpy.ndarray of int, shape (…), conditional
Index of symmetrically equivalent orientation that rotated vector to SST.
- in_SST(vector, proper=False)[source]#
Check whether given crystal frame vector falls into standard stereographic triangle of own symmetry.
- Parameters:
- vectornumpy.ndarray, shape (…,3)
Vector to check.
- properbool, optional
Consider only vectors with z >= 0, hence combine two neighboring SSTs. Defaults to False.
- Returns:
- innumpy.ndarray, shape (…)
Whether vector falls into SST.
- IPF_color(vector, in_SST=True, proper=False)[source]#
Map lab frame vector to RGB color within standard stereographic triangle of own symmetry.
- Parameters:
- vectornumpy.ndarray, shape (…,3)
Lab frame vector to colorize. Shape of vector blends with shape of own rotation array. For example, a rotation array of shape (3,2) and a vector array of shape (2,4) result in (3,2,4) outputs.
- 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, shape (…,3)
RGB array of IPF colors.
Examples
Inverse pole figure color of the e_3 lab direction for a crystal in “Cube” orientation with cubic symmetry:
>>> import damask >>> o = damask.Orientation(family='cubic') >>> o.IPF_color([0,0,1]) array([1., 0., 0.])
Sample standard triangle for hexagonal symmetry:
>>> import damask >>> from matplotlib import pyplot as plt >>> lab = [0,0,1] >>> o = damask.Orientation.from_random(shape=500000,family='hexagonal') >>> coord = damask.util.project_equal_area(o.to_SST(lab)) >>> color = o.IPF_color(lab) >>> plt.scatter(coord[:,0],coord[:,1],color=color,s=.06) [...] >>> plt.axis('scaled') [...] >>> plt.show()
- to_lattice(*, direction=None, plane=None)[source]#
Calculate lattice vector corresponding to lab frame direction or plane normal.
- Parameters:
- direction|planenumpy.ndarray, shape (…,3)
Real space vector along direction or reciprocal space vector along plane normal. Shape of vector blends with shape of own rotation array. For example, a rotation array of shape (3,2) and a vector array of shape (2,4) result in (3,2,4) outputs.
- Returns:
- Millernumpy.ndarray, shape (…,3)
Lattice vector of direction or plane. Use util.scale_to_coprime to convert to (integer) Miller indices.
Examples
>>> import numpy as np >>> import damask >>> cubic = damask.Orientation.from_axis_angle(n_omega=[1,0,0,90],degrees=True,lattice='cI') >>> cubic.to_lattice(direction=[1, 0, 0]) array([1., 0., 0.]) >>> cubic.to_lattice(direction=[0, 1, 0]) array([0., 0., -1.]) >>> cubic.to_lattice(direction=[0, 0, 1]) array([-0., 1., 0.]) >>> tetragonal = damask.Orientation(lattice='tI',c=0.5) >>> damask.util.scale_to_coprime(tetragonal.to_lattice(direction=[1,1,1])) array([1, 1, 2]) >>> damask.util.scale_to_coprime(tetragonal.to_lattice(plane=[1,1,1])) array([2, 2, 1])
- to_frame(*, uvw=None, hkl=None, uvtw=None, hkil=None, with_symmetry=False, normalize=True)[source]#
Calculate lab frame vector along lattice direction [uvw]/[uvtw] or plane normal (hkl)/(hkil).
- Parameters:
- uvw|hkl|uvtw|hkilnumpy.ndarray, shape (…,3) or shape (…,4)
Miller(–Bravais) indices of crystallographic direction or plane normal. Shape of vector blends with shape of own rotation array. For example, a rotation array of shape (3,2) and a vector array of shape (2,4) result in (3,2,4) outputs.
- with_symmetrybool, optional
Calculate all N symmetrically equivalent vectors. Defaults to False.
- normalizebool, optional
Normalize output vector. Defaults to True.
- Returns:
- vectornumpy.ndarray, shape (…,3) or (N,…,3)
Lab frame vector (or N vectors if with_symmetry) along [uvw]/[uvtw] direction or (hkl)/(hkil) plane normal.
- Schmid(*, N_slip=None, N_twin=None)[source]#
Calculate Schmid matrix P = d ⨂ n in the lab frame for selected deformation systems.
- Parameters:
- N_slip|N_twin‘*’ or sequence of int
Number of deformation systems per family of the deformation system. Use ‘*’ to select all.
- Returns:
- Pnumpy.ndarray, shape (N,…,3,3)
Schmid matrix for each of the N deformation systems.
Examples
Schmid matrix (in lab frame) of first octahedral slip system of a face-centered cubic crystal in “Goss” orientation.
>>> import numpy as np >>> import damask >>> O = damask.Orientation.from_Euler_angles(phi=[0,45,0],degrees=True,lattice='cF') >>> O.Schmid(N_slip=[12])[0] array([[ 0.000, 0.000, 0.000], [ 0.577, -0.000, 0.816], [ 0.000, 0.000, 0.000]])
All orientations related to self by given relationship model.
- Parameters:
- modelstr
Orientation relationship model selected from self.orientation_relationships.
- targetCrystal, optional
Crystal to transform to. Providing this parameter allows specification of non-standard lattice parameters. Default is inferred from selected model and uses standard lattice parameters.
- Returns:
- relOrientation, shape (:,self.shape)
Orientations related to self according to the selected model for the orientation relationship.
Examples
Face-centered cubic orientations following from a body-centered cubic crystal in “Cube” orientation according to the Bain orientation relationship (cF -> cI).
>>> import numpy as np >>> import damask >>> damask.Orientation(lattice='cF').related('Bain') Crystal family: cubic Bravais lattice: cI a=1 m, b=1 m, c=1 m α=90°, β=90°, γ=90° Quaternions of shape (3,) [[ 6.53281482e-01 2.70598050e-01 6.53281482e-01 2.70598050e-01] [ 2.70598050e-01 -2.70598050e-01 -6.53281482e-01 -6.53281482e-01] [ 9.23879533e-01 -5.55111512e-17 -2.77555756e-17 -3.82683432e-01]]
- append(other)#
Extend array along first dimension with other array(s).
- Parameters:
- other(list of) damask.Rotation
- apply(other)#
Return self@other.
Rotate vector, second-order tensor, or fourth-order tensor. other is interpreted as an array of tensor quantities with the highest-possible order considering the shape of self. Compatible innermost dimensions will blend. For instance, shapes of (2,) and (3,3) for self and other prompt interpretation of other as a second-rank tensor and result in (2,) rotated tensors, whereas shapes of (2,1) and (3,3) for self and other result in (2,3) rotated vectors.
- Parameters:
- othernumpy.ndarray, shape (…,3), (…,3,3), or (…,3,3,3,3)
Vector or tensor on which to apply the rotation.
- Returns:
- rotatednumpy.ndarray, shape (…,3), (…,3,3), or (…,3,3,3,3)
Rotated vector or tensor, i.e. transformed to frame defined by rotation.
Examples
Application of twelve (random) rotations to a set of five vectors.
>>> import numpy as np >>> import damask >>> r = damask.Rotation.from_random(shape=(12)) >>> o = np.ones((5,3)) >>> (r@o).shape # (12) @ (5, 3) (12, 5, 3)
Application of a (random) rotation to all twelve second-rank tensors.
>>> import numpy as np >>> import damask >>> r = damask.Rotation.from_random() >>> o = np.ones((12,3,3)) >>> (r@o).shape # (1) @ (12, 3,3) (12, 3, 3)
Application of twelve (random) rotations to the corresponding twelve second-rank tensors.
>>> import numpy as np >>> import damask >>> r = damask.Rotation.from_random(shape=(12)) >>> o = np.ones((12,3,3)) >>> (r@o).shape # (12) @ (3,3) (12, 3, 3)
Application of each of three (random) rotations to all three vectors.
>>> import numpy as np >>> import damask >>> r = damask.Rotation.from_random(shape=(3)) >>> o = np.ones((3,3)) >>> (r[...,np.newaxis]@o[np.newaxis,...]).shape # (3,1) @ (1,3, 3) (3, 3, 3)
Application of twelve (random) rotations to all twelve second-rank tensors.
>>> import numpy as np >>> import damask >>> r = damask.Rotation.from_random(shape=(12)) >>> o = np.ones((12,3,3)) >>> (r@o[np.newaxis,...]).shape # (12) @ (1,12, 3,3) (12, 12, 3, 3)
- as_Euler_angles(degrees=False)#
Represent as Bunge Euler angles.
- Parameters:
- degreesbool, optional
Return angles in degrees. Defaults to False.
- Returns:
- phinumpy.ndarray, 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 >>> damask.Rotation([1,0,0,0]).as_Euler_angles() array([0., 0., 0.])
- as_Rodrigues_vector(compact=False)#
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, 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 >>> damask.Rotation([1,0,0,0]).as_Rodrigues_vector(compact=True) array([ 0., 0., 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:
- n_omeganumpy.ndarray, 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 >>> damask.Rotation([1,0,0,0]).as_axis_angle(pair=True) (array([0., 0., 1.]), array(0.))
- as_cubochoric()#
Represent as cubochoric vector.
- Returns:
- xnumpy.ndarray, 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 >>> damask.Rotation([1,0,0,0]).as_cubochoric() array([0., 0., 0.])
- as_homochoric()#
Represent as homochoric vector.
- Returns:
- hnumpy.ndarray, shape (…,3)
Homochoric vector (h_1, h_2, h_3) with ǀhǀ < (3/4*π)^(1/3).
Examples
Cube orientation as homochoric vector.
>>> import damask >>> damask.Rotation([1,0,0,0]).as_homochoric() array([0., 0., 0.])
- as_matrix()#
Represent as rotation matrix.
- Returns:
- Rnumpy.ndarray, shape (…,3,3)
Rotation matrix R with det(R) = 1, R.T ∙ R = I.
Examples
Cube orientation as rotation matrix.
>>> import damask >>> damask.Rotation([1,0,0,0]).as_matrix() array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]])
- as_quaternion()#
Represent as unit quaternion.
- Returns:
- qnumpy.ndarray, shape (…,4)
Unit quaternion (q_0, q_1, q_2, q_3) in positive real hemisphere, i.e. ǀqǀ = 1, q_0 ≥ 0.
- property basis_real#
Return 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#
Return reciprocal (dual) crystal basis.
- broadcast_to(shape, mode='right')#
Broadcast array.
- Parameters:
- shape(sequence of) int
Shape of broadcasted array, needs to be compatible with the original shape.
- 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.
- Parameters:
- order{‘C’, ‘F’, ‘A’}, optional
‘C’ flattens in row-major (C-style) order. ‘F’ flattens in column-major (Fortran-style) order. ‘A’ flattens in column-major order if object is Fortran contiguous in memory, row-major order otherwise. Defaults to ‘C’.
- Returns:
- flatteneddamask.Rotation
Rotation flattened to single dimension.
- property immutable#
Return immutable lattice parameters.
- kinematics(mode)#
Return crystal kinematics systems.
- Parameters:
- mode{‘slip’,’twin’}
Deformation mode.
- Returns:
- direction_planedictionary
Directions and planes of deformation mode families.
- property lattice_points#
Return lattice points.
- misorientation(other)#
Calculate misorientation to other Rotation.
- Parameters:
- otherdamask.Rotation
Rotation to which the misorientation is computed. Compatible innermost dimensions will blend.
- Returns:
- gdamask.Rotation
Misorientation.
- property orientation_relationships#
Return labels of orientation relationships.
- property parameters#
Return lattice parameters a, b, c, alpha, beta, gamma.
- property ratio#
Return axes ratios of own lattice.
- relation_operations(model, target=None)#
Crystallographic orientation relationships for phase transformations.
- Parameters:
- modelstr
Name of orientation relationship.
- targetCrystal, optional
Crystal to transform to. Providing this parameter allows specification of non-standard lattice parameters. Default is inferred from selected model and uses standard lattice parameters.
- Returns:
- operations(string, damask.Rotation)
Resulting lattice and 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
- reshape(newshape, order='C')#
Reshape array.
- Parameters:
- newshape(sequence of) int
The new shape should be compatible with the original shape. If an integer, then the result will be a 1-D array of that length.
- order{‘C’, ‘F’, ‘A’}, optional
‘C’ flattens in row-major (C-style) order. ‘F’ flattens in column-major (Fortran-style) order. ‘A’ flattens in column-major order if object is Fortran contiguous in memory, row-major order otherwise. Defaults to ‘C’.
- Returns:
- reshapeddamask.Rotation
Rotation of given shape.
- property standard_triangle#
Corners of the standard triangle.
Notes
Not yet defined for monoclinic.
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 ... }
- property symmetry_operations#
Return symmetry operations.
References
U.F. Kocks et al., Texture and Anisotropy: Preferred Orientations in Polycrystals and their Effect on Materials Properties. Cambridge University Press 1998. Table II
tensor#
Tensor mathematics.
All routines operate on numpy.ndarrays of shape (…,3,3).
- damask.tensor.deviatoric(T)[source]#
Calculate deviatoric part of a tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the deviatoric part is computed.
- Returns:
- T’numpy.ndarray, shape (…,3,3)
Deviatoric part of T.
- damask.tensor.eigenvalues(T_sym)[source]#
Eigenvalues, i.e. principal components, of a symmetric tensor.
- Parameters:
- T_symnumpy.ndarray, shape (…,3,3)
Symmetric tensor of which the eigenvalues are computed.
- Returns:
- lambdanumpy.ndarray, shape (…,3)
Eigenvalues of T_sym sorted in ascending order, each repeated according to its multiplicity.
- damask.tensor.eigenvectors(T_sym, RHS=False)[source]#
Eigenvectors of a symmetric tensor.
- Parameters:
- T_symnumpy.ndarray, shape (…,3,3)
Symmetric tensor of which the eigenvectors are computed.
- RHS: bool, optional
Enforce right-handed coordinate system. Defaults to False.
- Returns:
- xnumpy.ndarray, shape (…,3,3)
Eigenvectors of T_sym sorted in ascending order of their associated eigenvalues.
- damask.tensor.spherical(T, tensor=True)[source]#
Calculate spherical part of a tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the spherical part is computed.
- tensorbool, optional
Map spherical part onto identity tensor. Defaults to True.
- Returns:
- pnumpy.ndarray, shape (…,3,3)
unless tensor == False: shape (…,) Spherical part of tensor T. p is an isotropic tensor.
- damask.tensor.symmetric(T)[source]#
Symmetrize tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the symmetrized values are computed.
- Returns:
- T_symnumpy.ndarray, shape (…,3,3)
Symmetrized tensor T.
- damask.tensor.transpose(T)[source]#
Transpose tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the transpose is computed.
- Returns:
- T.Tnumpy.ndarray, shape (…,3,3)
Transpose of tensor T.
mechanics#
Finite-strain continuum mechanics.
All routines operate on numpy.ndarrays of shape (…,3,3).
- damask.mechanics.deformation_Cauchy_Green_left(F)[source]#
Calculate left Cauchy-Green deformation tensor (Finger deformation tensor).
- Parameters:
- Fnumpy.ndarray, shape (…,3,3)
Deformation gradient.
- Returns:
- Bnumpy.ndarray, shape (…,3,3)
Left Cauchy-Green deformation tensor.
Notes
\[\vb{B} = \vb{F} \vb{F}^\text{T}\]
- damask.mechanics.deformation_Cauchy_Green_right(F)[source]#
Calculate right Cauchy-Green deformation tensor.
- Parameters:
- Fnumpy.ndarray, shape (…,3,3)
Deformation gradient.
- Returns:
- Cnumpy.ndarray, shape (…,3,3)
Right Cauchy-Green deformation tensor.
Notes
\[\vb{C} = \vb{F}^\text{T} \vb{F}\]
- damask.mechanics.equivalent_strain_Mises(epsilon)[source]#
Calculate the Mises equivalent of a strain tensor.
- Parameters:
- epsilonnumpy.ndarray, shape (…,3,3)
Symmetric strain tensor of which the von Mises equivalent is computed.
- Returns:
- epsilon_vMnumpy.ndarray, shape (…)
Von Mises equivalent strain of epsilon.
Notes
The von Mises equivalent of a strain tensor is defined as:
\[\epsilon_\text{vM} = \sqrt{\frac{2}{3}\,\epsilon^\prime_{ij} \epsilon^\prime_{ij}}\]where \(\vb*{\epsilon}^\prime\) is the deviatoric part of the strain tensor.
- damask.mechanics.equivalent_stress_Mises(sigma)[source]#
Calculate the Mises equivalent of a stress tensor.
- Parameters:
- sigmanumpy.ndarray, shape (…,3,3)
Symmetric stress tensor of which the von Mises equivalent is computed.
- Returns:
- sigma_vMnumpy.ndarray, shape (…)
Von Mises equivalent stress of sigma.
Notes
The von Mises equivalent of a stress tensor is defined as:
\[\sigma_\text{vM} = \sqrt{\frac{3}{2}\,\sigma^\prime_{ij} \sigma^\prime_{ij}}\]where \(\vb*{\sigma}^\prime\) is the deviatoric part of the stress tensor.
- damask.mechanics.maximum_shear(T_sym)[source]#
Calculate the maximum shear component of a symmetric tensor.
- Parameters:
- T_symnumpy.ndarray, shape (…,3,3)
Symmetric tensor of which the maximum shear is computed.
- Returns:
- gamma_maxnumpy.ndarray, shape (…)
Maximum shear of T_sym.
- damask.mechanics.rotation(T)[source]#
Calculate the rotational part of a tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the rotational part is computed.
- Returns:
- Rdamask.Rotation, shape (…)
Rotational part of the vector.
Notes
The rotational part is calculated from the polar decomposition:
\[\vb{R} = \vb{T} \vb{U}^{-1} = \vb{V}^{-1} \vb{T}\]where \(\vb{V}\) and \(\vb{U}\) are the left and right stretch tensor, respectively.
- damask.mechanics.strain(F, t, m)[source]#
Calculate strain tensor (Seth–Hill family).
- Parameters:
- Fnumpy.ndarray, shape (…,3,3)
Deformation gradient.
- t{‘V’, ‘U’}
Type of the polar decomposition, ‘V’ for left stretch tensor or ‘U’ for right stretch tensor.
- mfloat
Order of the strain.
- Returns:
- epsilonnumpy.ndarray, shape (…,3,3)
Strain of F.
Notes
The strain is defined as:
\[\begin{split}\vb*{\epsilon}_V^{(m)} = \frac{1}{2m} (\vb{V}^{2m} - \vb{I}) \\\\ \vb*{\epsilon}_U^{(m)} = \frac{1}{2m} (\vb{U}^{2m} - \vb{I})\end{split}\]References
- damask.mechanics.stress_Cauchy(P, F)[source]#
Calculate the Cauchy stress (true stress).
Resulting tensor is symmetrized as the Cauchy stress needs to be symmetric.
- Parameters:
- Pnumpy.ndarray, shape (…,3,3)
First Piola-Kirchhoff stress.
- Fnumpy.ndarray, shape (…,3,3)
Deformation gradient.
- Returns:
- sigmanumpy.ndarray, shape (…,3,3)
Cauchy stress.
- damask.mechanics.stress_second_Piola_Kirchhoff(P, F)[source]#
Calculate the second Piola-Kirchhoff stress.
Resulting tensor is symmetrized as the second Piola-Kirchhoff stress needs to be symmetric.
- Parameters:
- Pnumpy.ndarray, shape (…,3,3)
First Piola-Kirchhoff stress.
- Fnumpy.ndarray, shape (…,3,3)
Deformation gradient.
- Returns:
- Snumpy.ndarray, shape (…,3,3)
Second Piola-Kirchhoff stress.
- damask.mechanics.stretch_left(T)[source]#
Calculate left stretch of a tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the left stretch is computed.
- Returns:
- Vnumpy.ndarray, shape (…,3,3)
Left stretch tensor from Polar decomposition of T.
Notes
The left stretch tensor is calculated from the polar decomposition:
\[\vb{V} = \vb{T} \vb{R}^\text{T}\]where \(\vb{R}\) is a rotation.
- damask.mechanics.stretch_right(T)[source]#
Calculate right stretch of a tensor.
- Parameters:
- Tnumpy.ndarray, shape (…,3,3)
Tensor of which the right stretch is computed.
- Returns:
- Unumpy.ndarray, shape (…,3,3)
Left stretch tensor from Polar decomposition of T.
Notes
The right stretch tensor is calculated from the polar decomposition:
\[\vb{U} = \vb{R}^\text{T} \vb{T}\]where \(\vb{R}\) is a rotation.