Pre-Processing#
The following classes and modules support the setup of DAMASK simulations.
YAML#
- class damask.YAML(config=None, **kwargs)[source]#
YAML-based configuration.
- Attributes:
is_complete
Check for completeness.
is_valid
Check for valid file layout.
Methods
clear
()copy
()Return deepcopy(self).
delete
(keys)Remove configuration keys.
fromkeys
(iterable[, value])Create a new dictionary with keys from iterable and values set to value.
get
(key[, default])Return the value for key if key is in the dictionary, else default.
items
()keys
()load
(fname)Load from YAML file.
pop
(key[, default])If key is not found, default is returned if given, otherwise KeyError is raised
popitem
(/)Remove and return a (key, value) pair as a 2-tuple.
save
(fname, **kwargs)Save to YAML file.
setdefault
(key[, default])Insert key with a value of default if key is not in the dictionary.
update
([E, ]**F)If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
values
()- copy()#
Return deepcopy(self).
Create deep copy.
- delete(keys)[source]#
Remove configuration keys.
- Parameters:
- keysiterable or scalar
Label of the key(s) to remove.
- Returns:
- updateddamask.YAML
Updated configuration.
- classmethod load(fname)[source]#
Load from YAML file.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to read.
- Returns:
- loadeddamask.YAML
YAML from file.
- save(fname, **kwargs)[source]#
Save to YAML file.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to write.
- **kwargsdict
Keyword arguments parsed to yaml.dump.
- abstract property is_complete#
Check for completeness.
- abstract property is_valid#
Check for valid file layout.
- clear() None. Remove all items from D. #
- fromkeys(iterable, value=None, /)#
Create a new dictionary with keys from iterable and values set to value.
- get(key, default=None, /)#
Return the value for key if key is in the dictionary, else default.
- items() a set-like object providing a view on D's items #
- keys() a set-like object providing a view on D's keys #
- pop(key, default=<unrepresentable>, /)#
If key is not found, default is returned if given, otherwise KeyError is raised
- popitem(/)#
Remove and return a (key, value) pair as a 2-tuple.
Pairs are returned in LIFO (last-in, first-out) order. Raises KeyError if the dict is empty.
- setdefault(key, default=None, /)#
Insert key with a value of default if key is not in the dictionary.
Return the value for key if key is in the dictionary, else default.
- update([E, ]**F) None. Update D from dict/iterable E and F. #
If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
- values() an object providing a view on D's values #
ConfigMaterial#
- class damask.ConfigMaterial(config=None, *, homogenization=None, phase=None, material=None)[source]#
Material configuration.
Manipulate material configurations for storage in YAML format. A complete material configuration file has the entries ‘material’, ‘phase’, and ‘homogenization’.
- Attributes:
is_complete
Check for completeness.
is_valid
Check for valid content.
Methods
clear
()copy
()Return deepcopy(self).
delete
(keys)Remove configuration keys.
from_table
(table, *[, homogenization, ...])Generate from an ASCII table.
fromkeys
(iterable[, value])Create a new dictionary with keys from iterable and values set to value.
get
(key[, default])Return the value for key if key is in the dictionary, else default.
items
()keys
()load
(fname)Load from YAML file.
load_DREAM3D
(fname[, grain_data, cell_data, ...])Load DREAM.3D (HDF5) file.
material_add
(*[, homogenization, phase, v, ...])Add material entries.
material_rename_homogenization
(mapping[, ID])Change homogenization name in material.
material_rename_phase
(mapping[, ID, constituent])Change phase name in material.
pop
(key[, default])If key is not found, default is returned if given, otherwise KeyError is raised
popitem
(/)Remove and return a (key, value) pair as a 2-tuple.
save
(fname, **kwargs)Save to YAML file.
setdefault
(key[, default])Insert key with a value of default if key is not in the dictionary.
update
([E, ]**F)If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
values
()- static load_DREAM3D(fname, grain_data=None, cell_data=None, cell_ensemble_data='CellEnsembleData', phases='Phases', Euler_angles='EulerAngles', phase_names='PhaseName', base_group=None)[source]#
Load DREAM.3D (HDF5) file.
Data in DREAM.3D files can be stored per cell (‘CellData’) and/or per grain (‘Grain Data’). Per default, i.e. if ‘grain_data’ is None, cell-wise data is assumed.
- Parameters:
- fnamestr or pathlib.Path
Filename of the DREAM.3D (HDF5) file.
- grain_datastr
Name of the group (folder) containing grain-wise data. Defaults to None, in which case cell-wise data is used.
- cell_datastr
Name of the group (folder) containing cell-wise data. Defaults to None in wich case it is automatically detected.
- cell_ensemble_datastr
Name of the group (folder) containing data of cell ensembles. This group is used to inquire the name of the phases. Phases will get numeric IDs if this group is not found. Defaults to ‘CellEnsembleData’.
- phasesstr
Name of the dataset containing the phase ID (cell-wise or grain-wise). Defaults to ‘Phases’.
- Euler_anglesstr
Name of the dataset containing the crystallographic orientation as Euler angles in radians (cell-wise or grain-wise). Defaults to ‘EulerAngles’.
- phase_namesstr
Name of the dataset containing the phase names. Phases will get numeric IDs if this dataset is not found. Defaults to ‘PhaseName’.
- base_groupstr
Path to the group (folder) that contains geometry (_SIMPL_GEOMETRY), and grain- or cell-wise data. Defaults to None, in which case it is set as the path that contains _SIMPL_GEOMETRY/SPACING.
- Returns:
- loadeddamask.ConfigMaterial
Material configuration from file.
Notes
A grain-wise material configuration is based on segmented data from the DREAM.3D file. This data is typically available when the microstructure was synthetically created. In cell-wise representations, cells having the same orientation and phase are grouped. Since synthetically created microstructures have typically no in-grain scatter, cell-wise grids can appear to be segmented.
damask.GeomGrid.load_DREAM3D creates the corresponding grid-based geometry definition. Since the numbering of materials in cell-wise and grain-wise grids is different, it is imperative to use the same mode for both load_DREAM3D functions. That means, if the “grain_data” argument is used for this function, the correct grid configuration is only obtained if the “feature_IDs” argument is used when calling damask.GeomGrid.load_DREAM3D.
Homogenization and phase entries are emtpy and need to be defined separately.
- static from_table(table, *, homogenization=None, phase=None, v=None, O=None, V_e=None)[source]#
Generate from an ASCII table.
- Parameters:
- tabledamask.Table
Table that contains material information.
- homogenization: (array-like) of str, optional
Homogenization label.
- phase: (array-like) of str, optional
Phase label (per constituent).
- v: (array-like) of float or str, optional
Constituent volume fraction (per constituent). Defaults to 1/N_constituent.
- O: (array-like) of damask.Rotation or np.array/list of shape(4) or str, optional
Orientation as unit quaternion (per constituent).
- V_e: (array-like) of np.array/list of shape(3,3) or str, optional
Left elastic stretch (per constituent).
- Returns:
- newdamask.ConfigMaterial
Material configuration from values in table.
Notes
If the value of an argument is a string that is a column label, data from the table is used to fill the corresponding entry in the material configuration. Otherwise, the value is used directly.
First index of array-like values that are defined per constituent runs over materials, whereas second index runs over constituents.
Examples
>>> import damask >>> from damask import ConfigMaterial as cm >>> t = damask.Table.load('small.txt') >>> t 3:pos pos pos 4:qu qu qu qu phase homog 0 0 0 0 1.0 0.0 0.0 0.0 Aluminum SX 1 1 0 0 0.0 1.0 0.0 0.0 Steel SX 2 1 1 0 0.0 1.0 0.0 0.0 Steel SX
>>> cm.from_table(t,O='qu',phase='phase',homogenization='homog') homogenization: {SX: null} phase: {Aluminum: null, Steel: null} material: - constituents: - phase: Aluminum O: [1.0, 0.0, 0.0, 0.0] v: 1.0 homogenization: SX - constituents: - phase: Steel O: [0.0, 1.0, 0.0, 0.0] v: 1.0 homogenization: SX
>>> cm.from_table(t,O='qu',phase='phase',homogenization='single_crystal') homogenization: {single_crystal: null} phase: {Aluminum: null, Steel: null} material: - constituents: - phase: Aluminum O: [1.0, 0.0, 0.0, 0.0] v: 1.0 homogenization: single_crystal - constituents: - phase: Steel O: [0.0, 1.0, 0.0, 0.0] v: 1.0 homogenization: single_crystal
- property is_complete#
Check for completeness.
Only the general file layout is considered. This check does not consider whether specific parameters for a particular phase/homogenization model are missing.
- Returns:
- completebool
Whether the material.yaml definition is complete.
- property is_valid#
Check for valid content.
Only the generic file content is considered. This check does not consider whether parameters for a particular phase/homogenization mode are out of bounds.
- Returns:
- validbool
Whether the material.yaml definition is valid.
- material_rename_phase(mapping, ID=None, constituent=None)[source]#
Change phase name in material.
- Parameters:
- mapping: dict
Mapping from old name to new name
- ID: list of ints, optional
Limit renaming to selected material IDs.
- constituent: list of ints, optional
Limit renaming to selected constituents.
- Returns:
- updateddamask.ConfigMaterial
Updated material configuration.
- material_rename_homogenization(mapping, ID=None)[source]#
Change homogenization name in material.
- Parameters:
- mapping: dict
Mapping from old name to new name
- ID: list of ints, optional
Limit renaming to selected homogenization IDs.
- Returns:
- updateddamask.ConfigMaterial
Updated material configuration.
- material_add(*, homogenization=None, phase=None, v=None, O=None, V_e=None)[source]#
Add material entries.
- Parameters:
- homogenization: (array-like) of str, optional
Homogenization label.
- phase: (array-like) of str, optional
Phase label (per constituent).
- v: (array-like) of float, optional
Constituent volume fraction (per constituent). Defaults to 1/N_constituent.
- O: (array-like) of damask.Rotation or np.array/list of shape(4), optional
Orientation as unit quaternion (per constituent).
- V_e: (array-like) of np.array/list of shape(3,3), optional
Left elastic stretch (per constituent).
- Returns:
- updateddamask.ConfigMaterial
Updated material configuration.
Notes
First index of array-like values that are defined per constituent runs over materials, whereas second index runs over constituents.
Examples
Create two grains of ferrite and one grain of martensite, each with random orientation:
>>> import damask >>> m = damask.ConfigMaterial() >>> m = m.material_add(phase = ['Ferrite','Martensite','Ferrite'], ... O = damask.Rotation.from_random(3,rng_seed=20191102), ... homogenization = 'SX') >>> m homogenization: {SX: null} phase: {Ferrite: null, Martensite: null} material: - constituents: - phase: Ferrite O: [0.0047, -0.9582, 0.1084, 0.2645] v: 1.0 homogenization: SX - constituents: - phase: Martensite O: [0.9147, -0.1907, 0.2901, -0.2068] v: 1.0 homogenization: SX - constituents: - phase: Ferrite O: [0.1068, -0.4427, 0.1369, 0.8797] v: 1.0 homogenization: SX
Create five materials that each approximate a duplex stainless steel microstructure with three austenite and one relatively bigger ferrite grain of random orientation each:
>>> import numpy as np >>> import damask >>> m = damask.ConfigMaterial() >>> N_materials = 5 >>> m = m.material_add(phase = np.array([['Austenite']*3+['Ferrite']]), ... O = damask.Rotation.from_random((N_materials,4),rng_seed=20191102), ... v = np.array([[0.2]*3+[0.4]]), ... homogenization = 'Taylor') >>> m homogenization: {Taylor: null} phase: {Austenite: null, Ferrite: null} material: - constituents: - phase: Austenite v: 0.2 O: [0.004702411137213036, -0.9582446864633862, 0.1084379916089085, 0.2645490694937509] - phase: Austenite v: 0.2 O: [0.9147097460704486, -0.19068436891182194, 0.29014401444532145, -0.20678975501215882] - phase: Austenite v: 0.2 O: [0.10677819003833185, -0.4427133706883004, 0.13690394495734726, 0.879693468999888] - phase: Ferrite v: 0.4 O: [0.8664338002923555, 0.04448357787828491, -0.4945927532088464, 0.05188149461403649] homogenization: Taylor - constituents: - phase: Austenite v: 0.2 O: [0.5621873738314133, 0.0028841916095125584, -0.817023371343172, -0.1281009321680984] - phase: Austenite v: 0.2 O: [0.1566777437467901, -0.8117282158019414, 0.5096142534839398, 0.23841707348975383] - phase: Austenite v: 0.2 O: [0.3559036203819333, 0.1946923701552408, 0.058744995087853975, -0.9121274689178566] - phase: Ferrite v: 0.4 O: [0.467387781713959, -0.35644325887489176, 0.8031986430613528, 0.09679258489963502] homogenization: Taylor - constituents: - phase: Austenite v: 0.2 O: [0.4399087544327661, 0.12802483830067418, -0.8257167208737983, 0.32906203886337354] - phase: Austenite v: 0.2 O: [0.12410381094181624, -0.5125024631828828, -0.8493860709598213, 0.021972068647108236] - phase: Austenite v: 0.2 O: [0.03909373022192218, 0.4596226773046959, 0.42809626138739537, 0.7771436583738773] - phase: Ferrite v: 0.4 O: [0.737821660605232, 0.38809925187040367, -0.012129167758963711, 0.5521331824196455] homogenization: Taylor - constituents: - phase: Austenite v: 0.2 O: [0.4924738838478857, -0.0534798919571679, -0.6570981342247908, 0.5681825559468784] - phase: Austenite v: 0.2 O: [0.13073521303792138, 0.2534173177988532, -0.9582490914178947, -0.021133998872519554] - phase: Austenite v: 0.2 O: [0.1633346595899539, 0.6775968809652247, -0.07127256805012916, -0.71351557581203] - phase: Ferrite v: 0.4 O: [0.7658044627436773, -0.5327872540278646, 0.1102330397070761, 0.34282640467772235] homogenization: Taylor - constituents: - phase: Austenite v: 0.2 O: [0.25814496892598815, -0.6159961898524933, -0.5080223627084379, 0.543896265930874] - phase: Austenite v: 0.2 O: [0.8497433829153472, 0.4264182767672584, 0.05570674517418605, -0.3049596612218108] - phase: Austenite v: 0.2 O: [0.5146112784760113, 0.529467219604771, 0.661078636611197, 0.13347183839881469] - phase: Ferrite v: 0.4 O: [0.18430893147208752, 0.012407731059331692, -0.5551804816056372, -0.8109567798802285] homogenization: Taylor
- clear() None. Remove all items from D. #
- copy()#
Return deepcopy(self).
Create deep copy.
- delete(keys)#
Remove configuration keys.
- Parameters:
- keysiterable or scalar
Label of the key(s) to remove.
- Returns:
- updateddamask.YAML
Updated configuration.
- fromkeys(iterable, value=None, /)#
Create a new dictionary with keys from iterable and values set to value.
- get(key, default=None, /)#
Return the value for key if key is in the dictionary, else default.
- items() a set-like object providing a view on D's items #
- keys() a set-like object providing a view on D's keys #
- classmethod load(fname)#
Load from YAML file.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to read.
- Returns:
- loadeddamask.YAML
YAML from file.
- pop(key, default=<unrepresentable>, /)#
If key is not found, default is returned if given, otherwise KeyError is raised
- popitem(/)#
Remove and return a (key, value) pair as a 2-tuple.
Pairs are returned in LIFO (last-in, first-out) order. Raises KeyError if the dict is empty.
- save(fname, **kwargs)#
Save to YAML file.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to write.
- **kwargsdict
Keyword arguments parsed to yaml.dump.
- setdefault(key, default=None, /)#
Insert key with a value of default if key is not in the dictionary.
Return the value for key if key is in the dictionary, else default.
- update([E, ]**F) None. Update D from dict/iterable E and F. #
If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
- values() an object providing a view on D's values #
GeomGrid#
- class damask.GeomGrid(material, size, origin=array([0., 0., 0.]), initial_conditions=None, comments=None)[source]#
Geometry definition for grid solvers.
Create and manipulate geometry definitions for storage as VTK ImageData files (‘.vti’ extension). A grid has a physical size, a coordinate origin, and contains the material ID (indexing an entry in ‘material.yaml’) as well as initial condition fields.
- Attributes:
N_materials
Number of (unique) material indices within grid.
cells
Cell counts along x,y,z direction.
initial_conditions
Fields of initial conditions.
material
Material indices.
origin
Vector to grid origin in meter.
size
Edge lengths of grid in meter.
Methods
add_primitive
(dimension, center, exponent[, ...])Insert a primitive geometric object at a given position.
assemble
(idx)Assemble new grid from index map.
canvas
([cells, offset, fill])Crop or enlarge/pad grid.
clean
([distance, selection, ...])Smooth grid by selecting most frequent material ID within given stencil at each location.
copy
()Return deepcopy(self).
flip
(directions)Flip grid along given directions.
from_Laguerre_tessellation
(cells, size, ...)Create grid from Laguerre tessellation.
from_Voronoi_tessellation
(cells, size, seeds)Create grid from Voronoi tessellation.
from_minimal_surface
(cells, size, surface[, ...])Create grid from definition of triply-periodic minimal surface.
from_table
(table, coordinates, labels)Create grid from ASCII table.
get_grain_boundaries
([periodic, directions])Create VTK unstructured grid containing grain boundaries.
load
(fname)Load from VTK ImageData file with material IDs stored as 'material'.
load_ASCII
(fname)Load from geom file.
load_DREAM3D
(fname[, feature_IDs, ...])Load DREAM.3D (HDF5) file.
load_Neper
(fname)Load from Neper VTK file.
load_SPPARKS
(fname)Load from SPPARKS VTK dump.
mirror
(directions[, reflect])Mirror grid along given directions.
renumber
()Renumber sorted material indices as 0,...,N-1.
rotate
(R[, fill])Rotate grid (possibly extending its bounding box).
save
(fname[, compress])Save as VTK ImageData file.
save_ASCII
(fname)Save as geom file.
scale
(cells)Scale grid to new cell counts.
show
([colormap])Show on screen.
sort
()Sort material indices such that min(material ID) is located at (0,0,0).
substitute
(from_material, to_material)Substitute material indices.
vicinity_offset
([distance, offset, ...])Offset material ID of points in the vicinity of selected (or just other) material IDs.
- copy()#
Return deepcopy(self).
Create deep copy.
- property material#
Material indices.
- property size#
Edge lengths of grid in meter.
- property origin#
Vector to grid origin in meter.
- property initial_conditions#
Fields of initial conditions.
- property cells#
Cell counts along x,y,z direction.
- property N_materials#
Number of (unique) material indices within grid.
- static load(fname)[source]#
Load from VTK ImageData file with material IDs stored as ‘material’.
- Parameters:
- fnamestr or pathlib.Path
GeomGrid file to read. Valid extension is .vti, which will be appended if not given.
- Returns:
- loadeddamask.GeomGrid
Grid-based geometry from file.
- static load_SPPARKS(fname)[source]#
Load from SPPARKS VTK dump.
- Parameters:
- fnamestr or pathlib.Path
SPPARKS VTK dump file to read. Valid extension is .vti, which will be appended if not given.
- Returns:
- loadeddamask.GeomGrid
Grid-based geometry from file.
Notes
A SPPARKS VTI dump is equivalent to a DAMASK VTI file, but stores the materialID information as ‘Spin’ rather than ‘material’.
- static load_ASCII(fname)[source]#
Load from geom file.
Storing geometry files in ASCII format is deprecated. This function will be removed in a future version of DAMASK.
- Parameters:
- fnamestr, pathlib.Path, or file handle
Geometry file to read.
- Returns:
- loadeddamask.GeomGrid
Grid-based geometry from file.
- static load_Neper(fname)[source]#
Load from Neper VTK file.
- Parameters:
- fnamestr or pathlib.Path
Geometry file to read.
- Returns:
- loadeddamask.GeomGrid
Grid-based geometry from file.
Notes
Material indices in Neper usually start at 1 unless a buffer material with index 0 is added.
Examples
Read a periodic polycrystal generated with Neper.
>>> import damask >>> N_grains = 20 >>> cells = (32,32,32) >>> damask.util.run(f'neper -T -n {N_grains} -tesrsize {cells[0]}:{cells[1]}:{cells[2]} -periodicity all -format vtk') running 'neper -T -n 20 -tesrsize 32:32:32 -periodicity all -format vtk' ... >>> damask.GeomGrid.load_Neper(f'n{N_grains}-id1.vtk').renumber() cells: 32 × 32 × 32 size: 1.0 × 1.0 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 20
- static load_DREAM3D(fname, feature_IDs=None, cell_data=None, phases='Phases', Euler_angles='EulerAngles', base_group=None)[source]#
Load DREAM.3D (HDF5) file.
Data in DREAM.3D files can be stored per cell (‘CellData’) and/or per grain (‘Grain Data’). Per default, i.e. if ‘feature_IDs’ is None, cell-wise data is assumed.
- Parameters:
- fnamestr or pathlib.Path
Filename of the DREAM.3D (HDF5) file.
- feature_IDsstr, optional
Name of the dataset containing the mapping between cells and grain-wise data. Defaults to ‘None’, in which case cell-wise data is used.
- cell_datastr, optional
Name of the group (folder) containing cell-wise data. Defaults to None in wich case it is automatically detected.
- phasesstr, optional
Name of the dataset containing the phase ID. It is not used for grain-wise data, i.e. when feature_IDs is not None. Defaults to ‘Phases’.
- Euler_anglesstr, optional
Name of the dataset containing the crystallographic orientation as Euler angles in radians It is not used for grain-wise data, i.e. when feature_IDs is not None. Defaults to ‘EulerAngles’.
- base_groupstr, optional
Path to the group (folder) that contains geometry (_SIMPL_GEOMETRY), and grain- or cell-wise data. Defaults to None, in which case it is set as the path that contains _SIMPL_GEOMETRY/SPACING.
- Returns:
- loadeddamask.GeomGrid
Grid-based geometry from file.
Notes
A grain-wise geometry definition is based on segmented data from the DREAM.3D file. This data is typically available when the microstructure was synthetically created. In cell-wise representations, cells having the same orientation and phase are grouped. Since synthetically created microstructures have typically no in-grain scatter, cell-wise grids can appear to be segmented.
damask.ConfigMaterial.load_DREAM3D creates the corresponding material definition. Since the numbering of materials in cell-wise and grain-wise grids is different, it is imperative to use the same mode for both load_DREAM3D functions. That means, if the “feature_IDs” argument is used for this function, the correct material configuration is only obtained if the “grain_data” argument is used when calling damask.ConfigMaterial.load_DREAM3D.
- static from_table(table, coordinates, labels)[source]#
Create grid from ASCII table.
- Parameters:
- tabledamask.Table
Table that contains material information.
- coordinatesstr
Label of the vector column containing the spatial coordinates. Need to be ordered (1./x fast, 3./z slow).
- labels(sequence of) str
Label(s) of the columns containing the material definition. Each unique combination of values results in one material ID.
- Returns:
- newdamask.GeomGrid
Grid-based geometry from values in table.
- static from_Laguerre_tessellation(cells, size, seeds, weights, material=None, periodic=True)[source]#
Create grid from Laguerre tessellation.
- Parameters:
- cellssequence of int, len (3)
Cell counts along x,y,z direction.
- sizesequence of float, len (3)
Edge lengths of the grid in meter.
- seedsnumpy.ndarray of float, shape (:,3)
Position of the seed points in meter. All points need to lay within the box [(0,0,0),size].
- weightssequence of float, len (seeds.shape[0])
Weights of the seeds. Setting all weights to 1.0 gives a standard Voronoi tessellation.
- materialsequence of int, len (seeds.shape[0]), optional
Material ID of the seeds. Defaults to None, in which case materials are consecutively numbered.
- periodicbool, optional
Assume grid to be periodic. Defaults to True.
- Returns:
- newdamask.GeomGrid
Grid-based geometry from tessellation.
- static from_Voronoi_tessellation(cells, size, seeds, material=None, periodic=True)[source]#
Create grid from Voronoi tessellation.
- Parameters:
- cellssequence of int, len (3)
Cell counts along x,y,z direction.
- sizesequence of float, len (3)
Edge lengths of the grid in meter.
- seedsnumpy.ndarray of float, shape (:,3)
Position of the seed points in meter. All points need to lay within the box [(0,0,0),size].
- materialsequence of int, len (seeds.shape[0]), optional
Material ID of the seeds. Defaults to None, in which case materials are consecutively numbered.
- periodicbool, optional
Assume grid to be periodic. Defaults to True.
- Returns:
- newdamask.GeomGrid
Grid-based geometry from tessellation.
- static from_minimal_surface(cells, size, surface, threshold=0.0, periods=1, materials=(0, 1))[source]#
Create grid from definition of triply-periodic minimal surface.
- Parameters:
- cellssequence of int, len (3)
Cell counts along x,y,z direction.
- sizesequence of float, len (3)
Edge lengths of the grid in meter.
- surfacestr
Type of the minimal surface. See notes for details.
- thresholdfloat, optional.
Threshold of the minimal surface. Defaults to 0.0.
- periodsinteger, optional.
Number of periods per unit cell. Defaults to 1.
- materialssequence of int, len (2)
Material IDs. Defaults to (0,1).
- Returns:
- newdamask.GeomGrid
Grid-based geometry from definition of minimal surface.
Notes
- The following triply-periodic minimal surfaces are implemented:
Schwarz P
Double Primitive
Schwarz D
Complementary D
Double Diamond
Dprime
Gyroid
Gprime
Karcher K
Lidinoid
Neovius
Fisher-Koch S
References
S.B.G. Blanquer et al., Biofabrication 9(2):025001, 2017 https://doi.org/10.1088/1758-5090/aa6553
M. Wohlgemuth et al., Macromolecules 34(17):6083-6089, 2001 https://doi.org/10.1021/ma0019499
M.-T. Hsieh and L. Valdevit, Software Impacts 6:100026, 2020 https://doi.org/10.1016/j.simpa.2020.100026
Examples
Minimal surface of ‘Gyroid’ type.
>>> import numpy as np >>> import damask >>> damask.GeomGrid.from_minimal_surface([64]*3,np.ones(3)*1.e-4,'Gyroid') cells: 64 × 64 × 64 size: 0.0001 × 0.0001 × 0.0001 m³ origin: 0.0 0.0 0.0 m # materials: 2
Minimal surface of ‘Neovius’ type with specific material IDs.
>>> import numpy as np >>> import damask >>> damask.GeomGrid.from_minimal_surface([80]*3,np.ones(3)*5.e-4, ... 'Neovius',materials=(1,5)) cells: 80 × 80 × 80 size: 0.0005 × 0.0005 × 0.0005 m³ origin: 0.0 0.0 0.0 m # materials: 2 (min: 1, max: 5)
- save(fname, compress=True)[source]#
Save as VTK ImageData file.
- Parameters:
- fnamestr or pathlib.Path
Filename to write. Valid extension is .vti, which will be appended if not given.
- compressbool, optional
Compress with zlib algorithm. Defaults to True.
- save_ASCII(fname)[source]#
Save as geom file.
Storing geometry files in ASCII format is deprecated. This function will be removed in a future version of DAMASK.
- Parameters:
- fnamestr or file handle
Geometry file to write with extension ‘.geom’.
- compressbool, optional
Compress geometry using ‘x of y’ and ‘a to b’.
- show(colormap='cividis')[source]#
Show on screen.
- Parameters:
- colormapdamask.Colormap or str, optional
Colormap for visualization of material IDs. Defaults to ‘cividis’.
- canvas(cells=None, offset=None, fill=None)[source]#
Crop or enlarge/pad grid.
- Parameters:
- cellssequence of int, len (3), optional
Cell counts along x,y,z direction.
- offsetsequence of int, len (3), optional
Offset (measured in cells) from old to new grid. Defaults to [0,0,0].
- fillint, optional
Material ID to fill the background. Defaults to material.max()+1.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
Examples
Remove lower 1/2 of the microstructure in z-direction.
>>> import numpy as np >>> import damask >>> g = damask.GeomGrid(np.zeros([32]*3,int),np.ones(3)*1e-3) >>> g.canvas([32,32,16],[0,0,16]) cells: 32 × 32 × 16 size: 0.001 × 0.001 × 0.0005 m³ origin: 0.0 0.0 0.0005 m # materials: 1
- mirror(directions, reflect=False)[source]#
Mirror grid along given directions.
- Parameters:
- directions(sequence of) {‘x’, ‘y’, ‘z’}
Direction(s) along which the grid is mirrored.
- reflectbool, optional
Reflect (include) outermost layers. Defaults to False.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
Examples
Mirror along y-direction.
>>> import numpy as np >>> import damask >>> (g := damask.GeomGrid(np.arange(4*5*6).reshape([4,5,6]),np.ones(3))) cells: 4 × 5 × 6 size: 1.0 × 1.0 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 120 >>> g.mirror('y') cells: 4 × 8 × 6 size: 1.0 × 1.6 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 120
Reflect along x- and y-direction.
>>> g.mirror('xy',reflect=True) cells: 8 × 10 × 6 size: 2.0 × 2.0 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 120
Independence of mirroring order.
>>> g.mirror('xy') == g.mirror(['y','x']) True
- flip(directions)[source]#
Flip grid along given directions.
- Parameters:
- directions(sequence of) {‘x’, ‘y’, ‘z’}
Direction(s) along which the grid is flipped.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
Examples
Invariance of flipping order.
>>> import numpy as np >>> import damask >>> (g := damask.GeomGrid(np.arange(4*5*6).reshape([4,5,6]),np.ones(3))) cells: 4 × 5 × 6 size: 1.0 × 1.0 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 120 >>> g.flip('xyz') == g.flip(['x','z','y']) True
Invariance of flipping a (fully) mirrored grid.
>>> g.mirror('x',True) == g.mirror('x',True).flip('x') True
- rotate(R, fill=None)[source]#
Rotate grid (possibly extending its bounding box).
- Parameters:
- Rdamask.Rotation
Rotation to apply to the grid.
- fillint, optional
Material ID to fill enlarged bounding box. Defaults to material.max()+1.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
Examples
Rotation by 180° (π) is equivalent to twice flipping.
>>> import numpy as np >>> import damask >>> (g := damask.GeomGrid(np.arange(4*5*6).reshape([4,5,6]),np.ones(3))) cells: 4 × 5 × 6 size: 1.0 × 1.0 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 120 >>> g.rotate(damask.Rotation.from_axis_angle([0,0,1,180],degrees=True)) == g.flip('xy') True
- scale(cells)[source]#
Scale grid to new cell counts.
- Parameters:
- cellssequence of int, len (3)
Cell counts along x,y,z direction.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
Examples
Double grid resolution.
>>> import numpy as np >>> import damask >>> (g := damask.GeomGrid(np.zeros([32]*3,int),np.ones(3)*1e-4)) cells: 32 × 32 × 32 size: 0.0001 × 0.0001 × 0.0001 m³ origin: 0.0 0.0 0.0 m # materials: 1 >>> g.scale(g.cells*2) cells: 64 × 64 × 64 size: 0.0001 × 0.0001 × 0.0001 m³ origin: 0.0 0.0 0.0 m # materials: 1
- assemble(idx)[source]#
Assemble new grid from index map.
- Parameters:
- idxnumpy.ndarray of int, shape (:,:,:) or (:,:,:,3)
GeomGrid of flat indices or coordinate indices.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry. Cell count of resulting grid matches shape of index map.
- renumber()[source]#
Renumber sorted material indices as 0,…,N-1.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
- substitute(from_material, to_material)[source]#
Substitute material indices.
- Parameters:
- from_material(sequence of) int
Material indices to be substituted.
- to_material(sequence of) int
New material indices.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
- sort()[source]#
Sort material indices such that min(material ID) is located at (0,0,0).
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
- clean(distance=np.float64(1.7320508075688772), selection=None, invert_selection=False, periodic=True, rng_seed=None)[source]#
Smooth grid by selecting most frequent material ID within given stencil at each location.
- Parameters:
- distancefloat, optional
Voxel distance checked for presence of other materials. Defaults to sqrt(3).
- selection(sequence of) int, optional
Material IDs to consider. Defaults to all.
- invert_selectionbool, optional
Consider all material IDs except those in selection. Defaults to False.
- periodicbool, optional
Assume grid to be periodic. Defaults to True.
- 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:
- updateddamask.GeomGrid
Updated grid-based geometry.
Notes
If multiple material IDs are most frequent within a stencil, a random choice is taken.
- add_primitive(dimension, center, exponent, fill=None, R=Quaternion [1. 0. 0. 0.], inverse=False, periodic=True)[source]#
Insert a primitive geometric object at a given position.
- Parameters:
- dimensionsequence of int or float, len (3)
Dimension (diameter/side length) of the primitive. If given as integers, cell centers are addressed. If given as floats, physical coordinates are addressed.
- centersequence of int or float, len (3)
Center of the primitive. If given as integers, cell centers are addressed. If given as floats, physical coordinates are addressed.
- exponent(sequence of) float, len (3)
Exponents for the three axes. 0 gives octahedron (ǀxǀ^(2^0) + ǀyǀ^(2^0) + ǀzǀ^(2^0) < 1) 1 gives sphere (ǀxǀ^(2^1) + ǀyǀ^(2^1) + ǀzǀ^(2^1) < 1)
- fillint, optional
Fill value for primitive. Defaults to material.max()+1.
- Rdamask.Rotation, optional
Rotation of the primitive. Defaults to no rotation.
- inversebool, optional
Retain original materials within primitive and fill outside. Defaults to False.
- periodicbool, optional
Assume grid to be periodic. Defaults to True.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
Examples
Add a sphere at the center.
>>> import numpy as np >>> import damask >>> g = damask.GeomGrid(np.zeros([64]*3,int), np.ones(3)*1e-4) >>> g.add_primitive(np.ones(3)*5e-5,np.ones(3)*5e-5,1) cells: 64 × 64 × 64 size: 0.0001 × 0.0001 × 0.0001 m³ origin: 0.0 0.0 0.0 m # materials: 2
Add a cube at the origin.
>>> import numpy as np >>> import damask >>> g = damask.GeomGrid(np.zeros([64]*3,int), np.ones(3)*1e-4) >>> g.add_primitive(np.ones(3,int)*32,np.zeros(3),np.inf) cells: 64 × 64 × 64 size: 0.0001 × 0.0001 × 0.0001 m³ origin: 0.0 0.0 0.0 m # materials: 2
- vicinity_offset(distance=np.float64(1.7320508075688772), offset=None, selection=None, invert_selection=False, periodic=True)[source]#
Offset material ID of points in the vicinity of selected (or just other) material IDs.
Trigger points are variations in material ID, i.e. grain/phase boundaries or explicitly given material IDs.
- Parameters:
- distancefloat, optional
Voxel distance checked for presence of other materials. Defaults to sqrt(3).
- offsetint, optional
Offset (positive or negative) to tag material IDs. Defaults to material.max()+1.
- selection(sequence of) int, optional
Material IDs that trigger an offset. Defaults to any other than own material ID.
- invert_selectionbool, optional
Consider all material IDs except those in selection. Defaults to False.
- periodicbool, optional
Assume grid to be periodic. Defaults to True.
- Returns:
- updateddamask.GeomGrid
Updated grid-based geometry.
- get_grain_boundaries(periodic=True, directions='xyz')[source]#
Create VTK unstructured grid containing grain boundaries.
- Parameters:
- periodicbool, optional
Assume grid to be periodic. Defaults to True.
- directions(sequence of) {‘x’, ‘y’, ‘z’}, optional
Direction(s) along which the boundaries are determined. Defaults to ‘xyz’.
- Returns:
- grain_boundariesdamask.VTK
VTK-based geometry of grain boundary network.
LoadcaseGrid#
- class damask.LoadcaseGrid(config=None, *, solver=None, loadstep=None)[source]#
Load case for grid solver.
- Attributes:
is_complete
Check for completeness.
is_valid
Check for valid file layout.
Methods
clear
()copy
()Return deepcopy(self).
delete
(keys)Remove configuration keys.
fromkeys
(iterable[, value])Create a new dictionary with keys from iterable and values set to value.
get
(key[, default])Return the value for key if key is in the dictionary, else default.
items
()keys
()load
(fname)Load from YAML file.
pop
(key[, default])If key is not found, default is returned if given, otherwise KeyError is raised
popitem
(/)Remove and return a (key, value) pair as a 2-tuple.
save
(fname, **kwargs)Save to YAML file.
setdefault
(key[, default])Insert key with a value of default if key is not in the dictionary.
update
([E, ]**F)If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
values
()- save(fname, **kwargs)[source]#
Save to YAML file.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to write.
- **kwargsdict
Keyword arguments parsed to yaml.dump.
- clear() None. Remove all items from D. #
- copy()#
Return deepcopy(self).
Create deep copy.
- delete(keys)#
Remove configuration keys.
- Parameters:
- keysiterable or scalar
Label of the key(s) to remove.
- Returns:
- updateddamask.YAML
Updated configuration.
- fromkeys(iterable, value=None, /)#
Create a new dictionary with keys from iterable and values set to value.
- get(key, default=None, /)#
Return the value for key if key is in the dictionary, else default.
- abstract property is_complete#
Check for completeness.
- abstract property is_valid#
Check for valid file layout.
- items() a set-like object providing a view on D's items #
- keys() a set-like object providing a view on D's keys #
- classmethod load(fname)#
Load from YAML file.
- Parameters:
- fnamefile, str, or pathlib.Path
Filename or file to read.
- Returns:
- loadeddamask.YAML
YAML from file.
- pop(key, default=<unrepresentable>, /)#
If key is not found, default is returned if given, otherwise KeyError is raised
- popitem(/)#
Remove and return a (key, value) pair as a 2-tuple.
Pairs are returned in LIFO (last-in, first-out) order. Raises KeyError if the dict is empty.
- setdefault(key, default=None, /)#
Insert key with a value of default if key is not in the dictionary.
Return the value for key if key is in the dictionary, else default.
- update([E, ]**F) None. Update D from dict/iterable E and F. #
If E is present and has a .keys() method, then does: for k in E: D[k] = E[k] If E is present and lacks a .keys() method, then does: for k, v in E: D[k] = v In either case, this is followed by: for k in F: D[k] = F[k]
- values() an object providing a view on D's values #
seeds#
Functionality for generation of seed points for Voronoi or Laguerre tessellation.
- damask.seeds.from_random(size, N_seeds, cells=None, rng_seed=None)[source]#
Place seeds randomly in space.
- Parameters:
- sizesequence of float, len (3)
Edge lengths of the seeding domain.
- N_seedsint
Number of seeds.
- cellssequence of int, len (3), optional.
If given, ensures that each seed results in a grain when a standard Voronoi tessellation is performed using the given grid resolution (i.e. size/cells).
- rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS.
- Returns:
- coordsnumpy.ndarray, shape (N_seeds,3)
Seed coordinates in 3D space.
- damask.seeds.from_Poisson_disc(size, N_seeds, N_candidates, distance, periodic=True, rng_seed=None)[source]#
Place seeds following a Poisson disc distribution.
- Parameters:
- sizesequence of float, len (3)
Edge lengths of the seeding domain.
- N_seedsint
Number of seeds.
- N_candidatesint
Number of candidates to consider for finding best candidate.
- distancefloat
Minimum acceptable distance to other seeds.
- periodicbool, optional
Calculate minimum distance for periodically repeated grid. Defaults to True.
- rng_seed{None, int, array_like[ints], SeedSequence, BitGenerator, Generator}, optional
A seed to initialize the BitGenerator. Defaults to None. If None, then fresh, unpredictable entropy will be pulled from the OS.
- Returns:
- coordsnumpy.ndarray, shape (N_seeds,3)
Seed coordinates in 3D space.
- damask.seeds.from_grid(grid, selection=None, invert_selection=False, average=False, periodic=True)[source]#
Create seeds from grid description.
- Parameters:
- griddamask.GeomGrid
Grid from which the material IDs are used as seeds.
- selection(sequence of) int, optional
Material IDs to consider.
- invert_selectionbool, optional
Consider all material IDs except those in selection. Defaults to False.
- averagebool, optional
Seed corresponds to center of gravity of material ID cloud. Defaults to False.
- periodicbool, optional
Center of gravity accounts for periodic boundaries. Defaults to True.
- Returns:
- coords, materialsnumpy.ndarray, shape (:,3); numpy.ndarray, shape (:)
Seed coordinates in 3D space, material IDs.
Notes
The origin is not considered in order to obtain coordinates in a coordinate system located at the origin. This is expected by damask.GeomGrid.from_Voronoi_tessellation.
Examples
Recreate seeds from Voronoi tessellation.
>>> import numpy as np >>> import scipy.spatial >>> import damask >>> seeds = damask.seeds.from_random(np.ones(3),29,[128]*3,rng_seed=20191102) >>> (g := damask.GeomGrid.from_Voronoi_tessellation([128]*3,np.ones(3),seeds)) cells: 128 × 128 × 128 size: 1.0 × 1.0 × 1.0 m³ origin: 0.0 0.0 0.0 m # materials: 29 >>> COG,matID = damask.seeds.from_grid(g,average=True) >>> distance,ID = scipy.spatial.KDTree(COG,boxsize=g.size).query(seeds) >>> np.max(distance) / np.linalg.norm(g.size/g.cells) 10.1 >>> (ID == matID).all() True
Rotation#
- class damask.Rotation(rotation=array([1., 0., 0., 0.]))[source]#
Rotation with functionality for conversion between different representations.
The following conventions apply:
Coordinate frames are right-handed.
A rotation angle ω is taken to be positive for a counterclockwise rotation when viewing from the end point of the rotation axis towards the origin.
Rotations will be interpreted in the passive sense, i.e. as rotation of the coordinate frame.
P = -1 (as default).
References
D. Rowenhorst et al., Modelling and Simulation in Materials Science and Engineering 23:083501, 2015 https://doi.org/10.1088/0965-0393/23/8/083501
Examples
Rotate vector ‘a’ (defined in coordinate system ‘A’) to coordinates ‘b’ expressed in system ‘B’:
>>> import numpy as np >>> import damask >>> Q = damask.Rotation.from_random() >>> a = np.random.rand(3) >>> b = Q @ a >>> np.allclose(np.dot(Q.as_matrix(),a),b) True
Compound rotations R1 (first) and R2 (second):
>>> import numpy as np >>> import damask >>> R1 = damask.Rotation.from_random() >>> R2 = damask.Rotation.from_random() >>> R = R2 * R1 >>> np.allclose(R.as_matrix(), np.dot(R2.as_matrix(),R1.as_matrix())) True
- Attributes:
- quaternion
- shape
- size
Methods
allclose
(other[, rtol, atol, equal_nan])Test whether all values are approximately equal to corresponding ones of other Rotation.
append
(other)Extend array along first dimension with other array(s).
apply
(other)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])Average along last array dimension.
broadcast_to
(shape[, mode])Broadcast array.
copy
([rotation])Return deepcopy(self).
flatten
([order])Flatten array.
from_Euler_angles
(phi[, degrees])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, P])Initialize from Rodrigues–Frank vector (with angle separated from axis).
from_axis_angle
(n_omega[, degrees, normalize, P])Initialize from axis–angle pair.
from_basis
(basis[, orthonormal, reciprocal])Initialize from basis vector triplet.
from_cubochoric
(x[, P])Initialize from cubochoric vector.
from_fiber_component
(crystal, sample[, ...])Initialize with samples from a Gaussian distribution around a given direction.
from_homochoric
(h[, P])Initialize from homochoric vector.
from_matrix
(R[, normalize])Initialize from rotation matrix.
from_parallel
(a, b)Initialize from pairs of two orthogonal basis vectors.
from_quaternion
(q[, accept_homomorph, ...])Initialize from quaternion.
from_random
([shape, rng_seed])Initialize with samples from a uniform distribution.
from_spherical_component
(center, sigma[, ...])Initialize with samples from a Gaussian distribution around a given center.
isclose
(other[, rtol, atol, equal_nan])Report where values are approximately equal to corresponding ones of other Rotation.
misorientation
(other)Calculate misorientation to other Rotation.
reshape
(newshape[, order])Reshape array.
- 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 Rotation.
- Parameters:
- otherRotation
Rotation to compare against.
- rtolfloat, optional
Relative tolerance of equality.
- atolfloat, optional
Absolute tolerance of equality.
- equal_nanbool, optional
Consider matching NaN values as equal. Defaults to True.
- Returns:
- masknumpy.ndarray of bool, shape (self.shape)
Mask indicating where corresponding rotations are close.
- allclose(other, rtol=1e-05, atol=1e-08, equal_nan=True)[source]#
Test whether all values are approximately equal to corresponding ones of other Rotation.
- Parameters:
- otherRotation
Rotation to compare against.
- rtolfloat, optional
Relative tolerance of equality.
- atolfloat, optional
Absolute tolerance of equality.
- equal_nanbool, optional
Consider matching NaN values as equal. Defaults to True.
- Returns:
- answerbool
Whether all values are close between both rotations.
- apply(other)#
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)
- append(other)[source]#
Extend array along first dimension with other array(s).
- Parameters:
- other(list of) damask.Rotation
- flatten(order='C')[source]#
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.
- reshape(newshape, order='C')[source]#
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.
- broadcast_to(shape, mode='right')[source]#
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.
- average(weights=None)[source]#
Average along last array dimension.
- Parameters:
- weightsnumpy.ndarray, shape (self.shape), optional
Relative weight of each rotation.
- Returns:
- averagedamask.Rotation
Weighted average of original Rotation field.
References
F. Landis Markley et al., Journal of Guidance, Control, and Dynamics 30(4):1193-1197, 2007 https://doi.org/10.2514/1.28949
- misorientation(other)[source]#
Calculate misorientation to other Rotation.
- Parameters:
- otherdamask.Rotation
Rotation to which the misorientation is computed. Compatible innermost dimensions will blend.
- Returns:
- gdamask.Rotation
Misorientation.
- as_quaternion()[source]#
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.
- as_Euler_angles(degrees=False)[source]#
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_axis_angle(degrees=False, pair=False)[source]#
Represent as axis–angle pair.
- Parameters:
- degreesbool, optional
Return rotation angle in degrees. Defaults to False.
- pairbool, optional
Return tuple of axis and angle. Defaults to False.
- Returns:
- 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_matrix()[source]#
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_Rodrigues_vector(compact=False)[source]#
Represent as Rodrigues–Frank vector with separate axis and angle argument.
- Parameters:
- compactbool, optional
Return three-component Rodrigues–Frank vector, i.e. axis and angle argument are not separated.
- Returns:
- rhonumpy.ndarray, 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_homochoric()[source]#
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_cubochoric()[source]#
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.])
- static from_quaternion(q, accept_homomorph=False, normalize=False, P=-1)[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.
- Returns:
- newdamask.Rotation
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.]]
- static from_Euler_angles(phi, degrees=False)[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.
- Returns:
- newdamask.Rotation
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.]
- static from_axis_angle(n_omega, degrees=False, normalize=False, P=-1)[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.
- Returns:
- newdamask.Rotation
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. ]]
- static from_basis(basis, orthonormal=True, reciprocal=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.
- Returns:
- newdamask.Rotation
- static from_matrix(R, normalize=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.
- Returns:
- newdamask.Rotation
Examples
>>> import damask >>> damask.Rotation.from_matrix([[1,0,0],[0,0,-1],[0,1,0]]) Quaternion [ 0.707 -0.707 0. 0. ]
- static from_parallel(a, b)[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.
- Returns:
- newdamask.Rotation
Examples
>>> import damask >>> damask.Rotation.from_parallel([[2,0,0],[0,1,0]],[[1,0,0],[0,2,0]]) Quaternion [1. 0. 0. 0.]
- static from_Rodrigues_vector(rho, normalize=False, P=-1)[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.
- Returns:
- newdamask.Rotation
Examples
>>> import damask >>> damask.Rotation.from_Rodrigues_vector([0,0,1,1]) Quaternion [0.707 0. 0. 0.707]
- static from_homochoric(h, P=-1)[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.
- Returns:
- newdamask.Rotation
- static from_cubochoric(x, P=-1)[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.
- Returns:
- newdamask.Rotation
- static from_random(shape=None, rng_seed=None)[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.
- Returns:
- newdamask.Rotation
- static from_ODF(weights, phi, shape=None, degrees=False, fractions=True, rng_seed=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.
- Returns:
- newdamask.Rotation
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
- static from_spherical_component(center, sigma, shape=None, degrees=False, rng_seed=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.
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)
- static from_fiber_component(crystal, sample, sigma=0.0, shape=None, degrees=False, rng_seed=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.
- Returns:
- newdamask.Rotation
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)