Miscellaneous#
Auxiliary functionality, mainly used internally by the pre- and post-processing tools.
Table#
- class damask.Table(shapes={}, data=None, comments=None)[source]#
Manipulate multi-dimensional spreadsheet-like data.
- Attributes:
- labels
Methods
allclose
(other[, rtol, atol, equal_nan])Test whether all values are approximately equal to corresponding ones of other Table.
append
(other)Append other table vertically (similar to numpy.vstack).
copy
()Return deepcopy(self).
delete
(label)Delete column data.
get
(label)Get column data.
isclose
(other[, rtol, atol, equal_nan])Report where values are approximately equal to corresponding ones of other Table.
join
(other)Append other table horizontally (similar to numpy.hstack).
load
(fname)Load from ASCII table file.
load_ang
(fname[, shapes])Load from ANG file.
rename
(old, new[, info])Rename column data.
save
(fname[, with_labels])Save as plain text file.
set
(label, data[, info])Add new or replace existing column data.
sort_by
(labels[, ascending])Sort table by data of given columns.
- copy()#
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 Table.
- Parameters:
- otherdamask.Table
Table 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
Mask indicating where corresponding table values 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 Table.
- Parameters:
- otherdamask.Table
Table 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 corresponding values are close between both tables.
- static load(fname)[source]#
Load from ASCII table file.
Initial comments are marked by ‘#’. The first non-comment line contains the column labels.
Vector data column labels are indicated by ‘1_v, 2_v, …, n_v’.
Tensor data column labels are indicated by ‘3x3:1_T, 3x3:2_T, …, 3x3:9_T’.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to read.
- Returns:
- loadeddamask.Table
Table data from file.
- static load_ang(fname, shapes={'CI': 1, 'ID': 1, 'IQ': 1, 'eu': 3, 'fit': 1, 'intensity': 1, 'pos': 2})[source]#
Load from ANG file.
Regular ANG files feature the following columns:
Euler angles (Bunge notation) in radians, 3 floats, label ‘eu’.
Spatial position in meters, 2 floats, label ‘pos’.
Image quality, 1 float, label ‘IQ’.
Confidence index, 1 float, label ‘CI’.
Phase ID, 1 int, label ‘ID’.
SEM signal, 1 float, label ‘intensity’.
Fit, 1 float, label ‘fit’.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to read.
- shapesdict with str:int pairs, optional
Column labels and their width. Defaults to standard TSL ANG format.
- Returns:
- loadeddamask.Table
Table data from file.
- get(label)[source]#
Get column data.
- Parameters:
- labelstr
Column label.
- Returns:
- datanumpy.ndarray
Array of column data.
- set(label, data, info=None)[source]#
Add new or replace existing column data.
- Parameters:
- labelstr
Column label.
- datanumpy.ndarray
Column data. First dimension needs to match number of rows.
- infostr, optional
Human-readable information about the data.
- Returns:
- updateddamask.Table
Updated table.
- delete(label)[source]#
Delete column data.
- Parameters:
- labelstr
Column label.
- Returns:
- updateddamask.Table
Updated table.
- rename(old, new, info=None)[source]#
Rename column data.
- Parameters:
- label_old(iterable of) str
Old column labels.
- label_new(iterable of) str
New column labels.
- Returns:
- updateddamask.Table
Updated table.
- sort_by(labels, ascending=True)[source]#
Sort table by data of given columns.
- Parameters:
- labelstr or list
Column labels for sorting.
- ascendingbool or list, optional
Set sort order. Defaults to True.
- Returns:
- updateddamask.Table
Updated table.
- append(other)[source]#
Append other table vertically (similar to numpy.vstack).
Requires matching labels/shapes and order.
- Parameters:
- otherdamask.Table
Table to append.
- Returns:
- updateddamask.Table
Updated table.
Crystal#
- class damask.Crystal(*, family=None, lattice=None, a=None, b=None, c=None, alpha=None, beta=None, gamma=None, degrees=False)[source]#
Representation of a crystal as (general) crystal family or (more specific) as a scaled Bravais lattice.
Examples
Cubic crystal family:
>>> import damask >>> (cubic := damask.Crystal(family='cubic')) Crystal family: cubic
Body-centered cubic Bravais lattice with parameters of iron:
>>> import damask >>> (Fe := damask.Crystal(lattice='cI', a=287e-12)) Crystal family: cubic Bravais lattice: cI a=2.87e-10 m, b=2.87e-10 m, c=2.87e-10 m α=90°, β=90°, γ=90°
- Attributes:
basis_real
Return orthogonal real space crystal basis.
basis_reciprocal
Return reciprocal (dual) crystal basis.
immutable
Return immutable lattice parameters.
lattice_points
Return lattice points.
orientation_relationships
Return labels of orientation relationships.
parameters
Return lattice parameters a, b, c, alpha, beta, gamma.
ratio
Return axes ratios of own lattice.
standard_triangle
Corners of the standard triangle.
symmetry_operations
Return symmetry operations.
Methods
kinematics
(mode)Return crystal kinematics systems.
relation_operations
(model[, target])Crystallographic orientation relationships for phase transformations.
to_frame
(*[, uvw, hkl, uvtw, hkil])Calculate crystal frame vector corresponding to lattice direction [uvw]/[uvtw] or plane normal (hkl)/(hkil).
to_lattice
(*[, direction, plane])Calculate lattice vector corresponding to crystal frame direction or plane normal.
- property parameters#
Return lattice parameters a, b, c, alpha, beta, gamma.
- property immutable#
Return immutable lattice parameters.
- property orientation_relationships#
Return labels of orientation relationships.
- 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
- property ratio#
Return axes ratios of own lattice.
- 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.
- property lattice_points#
Return lattice points.
- to_lattice(*, direction=None, plane=None)[source]#
Calculate lattice vector corresponding to crystal frame direction or plane normal.
- Parameters:
- direction|planenumpy.ndarray, shape (…,3)
Real space vector along direction or reciprocal space vector along plane normal.
- Returns:
- Millernumpy.ndarray, 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, uvtw=None, hkil=None)[source]#
Calculate crystal frame vector corresponding to 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.
- Returns:
- vectornumpy.ndarray, shape (…,3)
Crystal frame vector in real space along [uvw]/[uvtw] direction or in reciprocal space along (hkl)/(hkil) plane normal.
Examples
Crystal frame vector (real space) of Magnesium corresponding to [1,1,0] direction:
>>> import damask >>> Mg = damask.Crystal(lattice='hP', a=321e-12, c=521e-12) >>> Mg.to_frame(uvw=[1, 1, 0]) array([1.60500000e-10, 2.77994155e-10, 0.00000000e+00])
Crystal frame vector (reciprocal space) of Titanium along (1,0,0) plane normal:
>>> import damask >>> Ti = damask.Crystal(lattice='hP', a=295e-12, c=469e-12) >>> Ti.to_frame(hkl=(1, 0, 0)) array([ 3.38983051e+09, 1.95711956e+09, -4.15134508e-07])
- kinematics(mode)[source]#
Return crystal kinematics systems.
- Parameters:
- mode{‘slip’,’twin’}
Deformation mode.
- Returns:
- direction_planedictionary
Directions and planes of deformation mode families.
- relation_operations(model, target=None)[source]#
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
VTK#
- class damask.VTK(vtk_data)[source]#
Spatial visualization (and potentially manipulation).
High-level interface to VTK.
- Attributes:
Methods
as_ASCII
()ASCII representation of the VTK data.
delete
(label)Delete either cell or point data.
from_image_data
(cells, size[, origin])Create VTK of type vtkImageData.
from_poly_data
(points)Create VTK of type polyData.
from_rectilinear_grid
(grid)Create VTK of type vtkRectilinearGrid.
from_unstructured_grid
(nodes, connectivity, ...)Create VTK of type vtkUnstructuredGrid.
get
(label)Get either cell or point data.
load
(fname[, dataset_type])Load from VTK file.
save
(fname[, parallel, compress])Save as VTK file.
set
([label, data, info, table])Add new or replace existing point or cell data.
show
([label, colormap])Render.
copy
- property comments#
Return the comments.
- property N_points#
Number of points in vtkdata.
- property N_cells#
Number of cells in vtkdata.
- property labels#
Labels of datasets.
- static from_image_data(cells, size, origin=array([0., 0., 0.]))[source]#
Create VTK of type vtkImageData.
This is the common type for grid solver results.
- Parameters:
- cellssequence of int, len (3)
Number of cells along each dimension.
- sizesequence of float, len (3)
Edge length along each dimension.
- originsequence of float, len (3), optional
Coordinates of grid origin.
- Returns:
- newdamask.VTK
VTK-based geometry without nodal or cell data.
- static from_unstructured_grid(nodes, connectivity, cell_type)[source]#
Create VTK of type vtkUnstructuredGrid.
This is the common type for mesh solver results.
- Parameters:
- nodesnumpy.ndarray, shape (:,3)
Spatial position of the nodes.
- connectivitynumpy.ndarray of np.dtype = np.int64
Cell connectivity (0-based), first dimension determines #Cells, second dimension determines #Nodes/Cell.
- cell_typestr
Name of the vtkCell subclass. Tested for TRIANGLE, QUAD, TETRA, and HEXAHEDRON.
- Returns:
- newdamask.VTK
VTK-based geometry without nodal or cell data.
- static from_poly_data(points)[source]#
Create VTK of type polyData.
This is the common type for point-wise data.
- Parameters:
- pointsnumpy.ndarray, shape (:,3)
Spatial position of the points.
- Returns:
- newdamask.VTK
VTK-based geometry without nodal or cell data.
- static from_rectilinear_grid(grid)[source]#
Create VTK of type vtkRectilinearGrid.
- Parameters:
- gridsequence of sequences of floats, len (3)
Grid coordinates along x, y, and z directions.
- Returns:
- newdamask.VTK
VTK-based geometry without nodal or cell data.
- static load(fname, dataset_type=None)[source]#
Load from VTK file.
- Parameters:
- fnamestr or pathlib.Path
Filename to read. Valid extensions are .vti, .vtu, .vtp, .vtr, and .vtk.
- dataset_type{‘ImageData’, ‘UnstructuredGrid’, ‘PolyData’, ‘RectilinearGrid’}, optional
Name of the vtkDataSet subclass when opening a .vtk file.
- Returns:
- loadeddamask.VTK
VTK-based geometry from file.
- save(fname, parallel=True, compress=True)[source]#
Save as VTK file.
- Parameters:
- fnamestr or pathlib.Path
Filename to write.
- parallelbool, optional
Write data in parallel background process. Defaults to True.
- compressbool, optional
Compress with zlib algorithm. Defaults to True.
- set(label=None, data=None, info=None, *, table=None)[source]#
Add new or replace existing point or cell data.
Data can either be a numpy.array, which requires a corresponding label, or a damask.Table.
- Parameters:
- labelstr, optional
Label of data array.
- datanumpy.ndarray or numpy.ma.MaskedArray, optional
Data to add or replace. First array dimension needs to match either number of cells or number of points.
- infostr, optional
Human-readable information about the data.
- table: damask.Table, optional
Data to add or replace. Each table label is individually considered. Number of rows needs to match either number of cells or number of points.
- Returns:
- updateddamask.VTK
Updated VTK-based geometry.
Notes
If the number of cells equals the number of points, the data is added to both.
- get(label)[source]#
Get either cell or point data.
Cell data takes precedence over point data, i.e. this function assumes that labels are unique among cell and point data.
- Parameters:
- labelstr
Data label.
- Returns:
- datanumpy.ndarray
Data stored under the given label.
- delete(label)[source]#
Delete either cell or point data.
Cell data takes precedence over point data, i.e. this function assumes that labels are unique among cell and point data.
- Parameters:
- labelstr
Data label.
- Returns:
- updateddamask.VTK
Updated VTK-based geometry.
- show(label=None, colormap='cividis')[source]#
Render.
- Parameters:
- labelstr, optional
Label of the dataset to show.
- colormapdamask.Colormap or str, optional
Colormap for visualization of dataset. Defaults to ‘cividis’.
Notes
The first component is shown when visualizing vector datasets (this includes tensor datasets as they are flattened).
util#
Miscellaneous helper functionality.
- damask.util.srepr(msg, glue='\n', quote=False)[source]#
Join (quoted) items with glue string.
- Parameters:
- msg(sequence of) object with __repr__
Items to join.
- gluestr, optional
Glue used for joining operation. Defaults to ‘n’.
- quotebool, optional
Quote items. Defaults to False.
- Returns:
- joinedstr
String representation of the joined and quoted items.
- damask.util.emph(msg)[source]#
Format with emphasis.
- Parameters:
- msg(sequence of) object with __repr__
Message to format.
- Returns:
- formattedstr
Formatted string representation of the joined items.
- damask.util.deemph(msg)[source]#
Format with deemphasis.
- Parameters:
- msg(sequence of) object with __repr__
Message to format.
- Returns:
- formattedstr
Formatted string representation of the joined items.
- damask.util.warn(msg)[source]#
Format for warning.
- Parameters:
- msg(sequence of) object with __repr__
Message to format.
- Returns:
- formattedstr
Formatted string representation of the joined items.
- damask.util.strikeout(msg)[source]#
Format as strikeout.
- Parameters:
- msg(iterable of) object with __repr__
Message to format.
- Returns:
- formattedstr
Formatted string representation of the joined items.
- damask.util.run(cmd, wd='./', env=None, timeout=None)[source]#
Run a command.
- Parameters:
- cmdstr
Command to be executed.
- wdstr, optional
Working directory of process. Defaults to ‘./’.
- envdict, optional
Environment for execution.
- timeoutinteger, optional
Timeout in seconds.
- Returns:
- stdout, stderr(str, str)
Output of the executed command.
- damask.util.open_text(fname, mode='r')[source]#
Open a text file with Unix line endings
If a path or string is given, a context manager ensures that the file handle is closed. If a file handle is given, it remains unmodified.
- Parameters:
- fnamefile, str, or pathlib.Path
Name or handle of file.
- mode: {‘r’,’w’}, optional
Access mode: ‘r’ead or ‘w’rite, defaults to ‘r’.
- Returns:
- ffile handle
- damask.util.execution_stamp(class_name, function_name=None)[source]#
Timestamp the execution of a (function within a) class.
- damask.util.show_progress(iterable, N_iter=None, prefix='', bar_length=50)[source]#
Decorate a loop with a progress bar.
Use similar like enumerate.
- Parameters:
- iterableiterable
Iterable to be decorated.
- N_iterint, optional
Total number of iterations. Required if iterable is not a sequence.
- prefixstr, optional
Prefix string. Defaults to ‘’.
- bar_lengthint, optional
Length of progress bar in characters. Defaults to 50.
- damask.util.scale_to_coprime(v, N_significant=9)[source]#
Scale vector to co-prime (relatively prime) integers.
- Parameters:
- vsequence of float, len (:)
Vector to scale.
- N_significant: int, optional
Number of significant digits to consider. Defaults to 9.
- Returns:
- mnumpy.ndarray, shape (:)
Vector scaled to co-prime numbers.
- damask.util.project_equal_angle(vector, direction='z', normalize=True, keepdims=False)[source]#
Apply equal-angle projection to vector.
- Parameters:
- vectornumpy.ndarray, shape (…,3)
Vector coordinates to be projected.
- direction{‘x’, ‘y’, ‘z’}
Projection direction. Defaults to ‘z’.
- normalizebool
Ensure unit length of input vector. Defaults to True.
- keepdimsbool
Maintain three-dimensional output coordinates. Defaults to False.
- Returns:
- coordinatesnumpy.ndarray, shape (…,2 | 3)
Projected coordinates.
Notes
Two-dimensional output uses right-handed frame spanned by the next and next-next axis relative to the projection direction, e.g. x-y when projecting along z and z-x when projecting along y.
Examples
>>> import damask >>> import numpy as np >>> project_equal_angle(np.ones(3)) array([0.3660, 0.3660]) >>> project_equal_angle(np.ones(3),direction='x',normalize=False,keepdims=True) array([0. , 0.5, 0.5]) >>> project_equal_angle([0,1,1],direction='y',normalize=True,keepdims=False) array([0.4142, 0. ])
- damask.util.project_equal_area(vector, direction='z', normalize=True, keepdims=False)[source]#
Apply equal-area projection to vector.
- Parameters:
- vectornumpy.ndarray, shape (…,3)
Vector coordinates to be projected.
- direction{‘x’, ‘y’, ‘z’}
Projection direction. Defaults to ‘z’.
- normalizebool
Ensure unit length of input vector. Defaults to True.
- keepdimsbool
Maintain three-dimensional output coordinates. Defaults to False.
- Returns:
- coordinatesnumpy.ndarray, shape (…,2 | 3)
Projected coordinates.
Notes
Two-dimensional output uses right-handed frame spanned by the next and next-next axis relative to the projection direction, e.g. x-y when projecting along z and z-x when projecting along y.
Examples
>>> import damask >>> import numpy as np >>> project_equal_area(np.ones(3)) array([0.4597, 0.4597]) >>> project_equal_area(np.ones(3),direction='x',normalize=False,keepdims=True) array([0. , 0.7071, 0.7071]) >>> project_equal_area([0,1,1],direction='y',normalize=True,keepdims=False) array([0.5412, 0. ])
- damask.util.hybrid_IA(dist, N, rng_seed=None)[source]#
Hybrid integer approximation.
- Parameters:
- distnumpy.ndarray
Distribution to be approximated.
- Nint
Number of samples to draw.
- rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS.
- Returns:
- histnumpy.ndarray, shape (N)
Integer approximation of the distribution.
- damask.util.shapeshifter(fro, to, mode='left', keep_ones=False)[source]#
Return dimensions that reshape ‘fro’ to become broadcastable to ‘to’.
- Parameters:
- frotuple
Original shape of array.
- totuple
Target shape of array after broadcasting. len(to) cannot be less than len(fro).
- mode{‘left’, ‘right’}, optional
Indicates whether new axes are preferably added to either left or right of the original shape. Defaults to ‘left’.
- keep_onesbool, optional
Treat ‘1’ in fro as literal value instead of dimensional placeholder. Defaults to False.
- Returns:
- new_dimstuple
Dimensions for reshape.
Examples
>>> import numpy as np >>> from damask import util >>> a = np.ones((3,4,2)) >>> b = np.ones(4) >>> b_extended = b.reshape(util.shapeshifter(b.shape,a.shape)) >>> (a * np.broadcast_to(b_extended,a.shape)).shape (3, 4, 2)
- damask.util.shapeblender(a, b, keep_ones=False)[source]#
Return a shape that overlaps the rightmost entries of ‘a’ with the leftmost of ‘b’.
- Parameters:
- atuple
Shape of first array.
- btuple
Shape of second array.
- keep_onesbool, optional
Treat innermost ‘1’s as literal value instead of dimensional placeholder. Defaults to False.
Examples
>>> shapeblender((3,2),(3,2)) (3, 2) >>> shapeblender((4,3),(3,2)) (4, 3, 2) >>> shapeblender((4,4),(3,2)) (4, 4, 3, 2) >>> shapeblender((1,2),(1,2,3)) (1, 2, 3) >>> shapeblender((),(2,2,1)) (2, 2, 1) >>> shapeblender((1,),(2,2,1)) (2, 2, 1) >>> shapeblender((1,),(2,2,1),True) (1, 2, 2, 1)
- damask.util.extend_docstring(docstring=None, **kwargs)[source]#
Decorator: Extend the function’s docstring.
- Parameters:
- docstringstr or callable, optional
Docstring to extend. Defaults to that of decorated function.
- adopted_*str or callable, optional
Additional information to insert into/append to respective section.
Notes
Return type will become own type if docstring is callable.
- damask.util.pass_on(keyword, target, wrapped=None)[source]#
Decorator: Combine signatures of ‘wrapped’ and ‘target’ functions and pass on output of ‘target’ as ‘keyword’ argument.
- Parameters:
- keywordstr
Keyword added to **kwargs of the decorated function passing on the result of ‘target’.
- targetcallable
The output of this function is passed to the decorated function as ‘keyword’ argument.
- wrapped: callable, optional
Signature of ‘wrapped’ function combined with that of ‘target’ yields the overall signature of decorated function.
Notes
The keywords used by ‘target’ will be prioritized if they overlap with those of the decorated function. Functions ‘target’ and ‘wrapped’ are assumed to only have keyword arguments.
- damask.util.DREAM3D_base_group(fname)[source]#
Determine the base group of a DREAM.3D file.
The base group is defined as the group (folder) that contains a ‘SPACING’ dataset in a ‘_SIMPL_GEOMETRY’ group.
- Parameters:
- fnamestr, pathlib.Path, or _h5py.File
Filename of the DREAM.3D (HDF5) file.
- Returns:
- pathstr
Path to the base group.
- damask.util.DREAM3D_cell_data_group(fname)[source]#
Determine the cell data group of a DREAM.3D file.
The cell data group is defined as the group (folder) that contains a dataset in the base group whose length matches the total number of points as specified in ‘_SIMPL_GEOMETRY/DIMENSIONS’.
- Parameters:
- fnamestr, pathlib.Path, or h5py.File
Filename of the DREAM.3D (HDF5) file.
- Returns:
- pathstr
Path to the cell data group.
- damask.util.Bravais_to_Miller(*, uvtw=None, hkil=None)[source]#
Transform 4 Miller–Bravais indices to 3 Miller indices of crystal direction [uvw] or plane normal (hkl).
- Parameters:
- uvtw|hkilnumpy.ndarray, shape (…,4)
Miller–Bravais indices of crystallographic direction [uvtw] or plane normal (hkil).
- Returns:
- uvw|hklnumpy.ndarray, shape (…,3)
Miller indices of [uvw] direction or (hkl) plane normal.
- damask.util.Miller_to_Bravais(*, uvw=None, hkl=None)[source]#
Transform 3 Miller indices to 4 Miller–Bravais indices of crystal direction [uvtw] or plane normal (hkil).
- Parameters:
- uvw|hklnumpy.ndarray, shape (…,3)
Miller indices of crystallographic direction [uvw] or plane normal (hkl).
- Returns:
- uvtw|hkilnumpy.ndarray, shape (…,4)
Miller–Bravais indices of [uvtw] direction or (hkil) plane normal.
- damask.util.dict_prune(d)[source]#
Recursively remove empty dictionaries.
- Parameters:
- ddict
Dictionary to prune.
- Returns:
- pruneddict
Pruned dictionary.
- damask.util.dict_flatten(d)[source]#
Recursively remove keys of single-entry dictionaries.
- Parameters:
- ddict
Dictionary to flatten.
- Returns:
- flatteneddict
Flattened dictionary.
- class damask.util.ProgressBar(total, prefix, bar_length)[source]#
Report progress of an interation as a status bar.
Works for 0-based loops, ETA is estimated by linear extrapolation.
Methods
update
grid_filters#
Filters for operations on regular grids.
The grids are defined as (x,y,z,…) where x is fastest and z is slowest. This convention is consistent with the layout in grid vti files.
When converting to/from a plain list (e.g. storage in ASCII table), the following operations are required for tensorial data:
D3 = D1.reshape(cells+(-1,),order=’F’).reshape(cells+(3,3))
D1 = D3.reshape(cells+(-1,)).reshape(-1,9,order=’F’)
- damask.grid_filters.curl(size, f)[source]#
Calculate curl of a vector or tensor field in Fourier space.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- fnumpy.ndarray, shape (:,:,:,3) or (:,:,:,3,3)
Periodic field of which the curl is calculated.
- Returns:
- ∇ × fnumpy.ndarray, shape (:,:,:,3) or (:,:,:,3,3)
Curl of f.
- damask.grid_filters.divergence(size, f)[source]#
Calculate divergence of a vector or tensor field in Fourier space.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- fnumpy.ndarray, shape (:,:,:,3) or (:,:,:,3,3)
Periodic field of which the divergence is calculated.
- Returns:
- ∇ · fnumpy.ndarray, shape (:,:,:,1) or (:,:,:,3)
Divergence of f.
- damask.grid_filters.gradient(size, f)[source]#
Calculate gradient of a scalar or vector field in Fourier space.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- fnumpy.ndarray, shape (:,:,:,1) or (:,:,:,3)
Periodic field of which the gradient is calculated.
- Returns:
- ∇ fnumpy.ndarray, shape (:,:,:,3) or (:,:,:,3,3)
Gradient of f.
- damask.grid_filters.coordinates0_point(cells, size, origin=array([0., 0., 0.]))[source]#
Cell center positions (undeformed).
- Parameters:
- cellssequence of int, len (3)
Number of cells.
- sizesequence of float, len (3)
Physical size of the periodic field.
- originsequence of float, len(3), optional
Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].
- Returns:
- x_p_0numpy.ndarray, shape (:,:,:,3)
Undeformed cell center coordinates.
- damask.grid_filters.displacement_fluct_point(size, F)[source]#
Cell center displacement field from fluctuation part of the deformation gradient field.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- Returns:
- u_p_fluctnumpy.ndarray, shape (:,:,:,3)
Fluctuating part of the cell center displacements.
- damask.grid_filters.displacement_avg_point(size, F)[source]#
Cell center displacement field from average part of the deformation gradient field.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- Returns:
- u_p_avgnumpy.ndarray, shape (:,:,:,3)
Average part of the cell center displacements.
- damask.grid_filters.displacement_point(size, F)[source]#
Cell center displacement field from deformation gradient field.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- Returns:
- u_pnumpy.ndarray, shape (:,:,:,3)
Cell center displacements.
- damask.grid_filters.coordinates_point(size, F, origin=array([0., 0., 0.]))[source]#
Cell center positions.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- originsequence of float, len(3), optional
Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].
- Returns:
- x_pnumpy.ndarray, shape (:,:,:,3)
Cell center coordinates.
- damask.grid_filters.cellsSizeOrigin_coordinates0_point(coordinates0, ordered=True)[source]#
Return grid ‘DNA’, i.e. cells, size, and origin from 1D array of point positions.
- Parameters:
- coordinates0numpy.ndarray, shape (:,3)
Undeformed cell center coordinates.
- orderedbool, optional
Expect coordinates0 data to be ordered (x fast, z slow). Defaults to True.
- Returns:
- cells, size, originThree numpy.ndarray, each of shape (3)
Information to reconstruct grid.
- damask.grid_filters.coordinates0_node(cells, size, origin=array([0., 0., 0.]))[source]#
Nodal positions (undeformed).
- Parameters:
- cellssequence of int, len (3)
Number of cells.
- sizesequence of float, len (3)
Physical size of the periodic field.
- originsequence of float, len(3), optional
Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].
- Returns:
- x_n_0numpy.ndarray, shape (:,:,:,3)
Undeformed nodal coordinates.
- damask.grid_filters.displacement_fluct_node(size, F)[source]#
Nodal displacement field from fluctuation part of the deformation gradient field.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- Returns:
- u_n_fluctnumpy.ndarray, shape (:,:,:,3)
Fluctuating part of the nodal displacements.
- damask.grid_filters.displacement_avg_node(size, F)[source]#
Nodal displacement field from average part of the deformation gradient field.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- Returns:
- u_n_avgnumpy.ndarray, shape (:,:,:,3)
Average part of the nodal displacements.
- damask.grid_filters.displacement_node(size, F)[source]#
Nodal displacement field from deformation gradient field.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- Returns:
- u_nnumpy.ndarray, shape (:,:,:,3)
Nodal displacements.
- damask.grid_filters.coordinates_node(size, F, origin=array([0., 0., 0.]))[source]#
Nodal positions.
- Parameters:
- sizesequence of float, len (3)
Physical size of the periodic field.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field.
- originsequence of float, len(3), optional
Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].
- Returns:
- x_nnumpy.ndarray, shape (:,:,:,3)
Nodal coordinates.
- damask.grid_filters.cellsSizeOrigin_coordinates0_node(coordinates0, ordered=True)[source]#
Return grid ‘DNA’, i.e. cells, size, and origin from 1D array of nodal positions.
- Parameters:
- coordinates0numpy.ndarray, shape (:,3)
Undeformed nodal coordinates.
- orderedbool, optional
Expect coordinates0 data to be ordered (x fast, z slow). Defaults to True.
- Returns:
- cells, size, originThree numpy.ndarray, each of shape (3)
Information to reconstruct grid.
- damask.grid_filters.point_to_node(cell_data)[source]#
Interpolate periodic point data to nodal data.
- Parameters:
- cell_datanumpy.ndarray, shape (:,:,:,…)
Data defined on the cell centers of a periodic grid.
- Returns:
- node_datanumpy.ndarray, shape (:,:,:,…)
Data defined on the nodes of a periodic grid.
- damask.grid_filters.node_to_point(node_data)[source]#
Interpolate periodic nodal data to point data.
- Parameters:
- node_datanumpy.ndarray, shape (:,:,:,…)
Data defined on the nodes of a periodic grid.
- Returns:
- cell_datanumpy.ndarray, shape (:,:,:,…)
Data defined on the cell centers of a periodic grid.
- damask.grid_filters.coordinates0_valid(coordinates0)[source]#
Check whether coordinates form a regular grid.
- Parameters:
- coordinates0numpy.ndarray, shape (:,3)
Array of undeformed cell coordinates.
- Returns:
- validbool
Whether the coordinates form a regular grid.
- damask.grid_filters.unravel_index(idx)[source]#
Convert flat indices to coordinate indices.
- Parameters:
- idxnumpy.ndarray, shape (:,:,:)
Grid of flat indices.
- Returns:
- unravellednumpy.ndarray, shape (:,:,:,3)
Grid of coordinate indices.
Examples
Unravel a linearly increasing sequence of material indices on a 3 × 2 × 1 grid.
>>> import numpy as np >>> import damask >>> seq = np.arange(6).reshape((3,2,1),order='F') >>> (coord_idx := damask.grid_filters.unravel_index(seq)) array([[[[0, 0, 0]], [[0, 1, 0]]], [[[1, 0, 0]], [[1, 1, 0]]], [[[2, 0, 0]], [[2, 1, 0]]]]) >>> coord_idx[1,1,0] array([1, 1, 0])
- damask.grid_filters.ravel_index(idx)[source]#
Convert coordinate indices to flat indices.
- Parameters:
- idxnumpy.ndarray, shape (:,:,:,3)
Grid of coordinate indices.
- Returns:
- ravellednumpy.ndarray, shape (:,:,:)
Grid of flat indices.
Examples
Ravel a reversed sequence of coordinate indices on a 2 × 2 × 1 grid.
>>> import numpy as np >>> import damask >>> (rev := np.array([[1,1,0],[0,1,0],[1,0,0],[0,0,0]]).reshape((2,2,1,3))) array([[[[1, 1, 0]], [[0, 1, 0]]], [[[1, 0, 0]], [[0, 0, 0]]]]) >>> (flat_idx := damask.grid_filters.ravel_index(rev)) array([[[3], [2]], [[1], [0]]])
- damask.grid_filters.regrid(size, F, cells)[source]#
Map a deformed grid A back to a rectilinear grid B.
The size of grid B is chosen as the average deformed size of grid A.
- Parameters:
- sizesequence of float, len (3)
Physical size of grid A.
- Fnumpy.ndarray, shape (:,:,:,3,3)
Deformation gradient field on grid A.
- cellssequence of int, len (3)
Cell count along x,y,z of grid B.
- Returns:
- idxnumpy.ndarray of int, shape (cells)
Flat index of closest point on deformed grid A for each point on grid B.
solver.Marc#
- class damask.solver.Marc(version='2023.4', marc_root='/opt/msc', damask_root='/tmp/12923_3b8e5bf3600ceb007b03a3bbd2b2e4e0d7e07133/DAMASK')[source]#
Wrapper to run DAMASK with MSC Marc.
- Attributes:
- library_path
- tools_path
Methods
submit_job
(model, job[, compile, ...])Assemble command line arguments and call Marc executable.
- submit_job(model, job, compile=False, optimization='', env=None)[source]#
Assemble command line arguments and call Marc executable.
- Parameters:
- modelstr
Name of model.
- jobstr
Name of job.
- compilebool, optional
Compile DAMASK_Marc user subroutine (and save for future use). Defaults to False.
- optimization{‘’, ‘l’, ‘h’}, optional
Optimization level ‘’: -O0, ‘l’: -O1, or ‘h’: -O3. Defaults to ‘’.
- envdict, optional
Environment for execution.