Miscellaneous

Auxiliary functionality, mainly used internally by the pre- and post-processing tools.

Table

class damask.Table(data, shapes, comments=None)[source]

Manipulate multi-dimensional spreadsheet-like data.

Attributes
labels

Methods

add(label, data[, info])

Add column data.

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()

Create deep copy.

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)

Load from ang file.

rename(old, new[, info])

Rename column data.

save(fname)

Save as plain text file.

set(label, data[, info])

Set column data.

sort_by(labels[, ascending])

Sort table by values of given labels.

copy()

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 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 containing 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 for reading.

Returns
loadeddamask.Table

Table data from file.

static load_ang(fname)[source]

Load from ang file.

A valid TSL ang file has to have 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 for reading.

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]

Set column data.

Parameters
labelstr

Column label.

datanumpy.ndarray

New data.

infostr, optional

Human-readable information about the new data.

Returns
updateddamask.Table

Updated table.

add(label, data, info=None)[source]

Add column data.

Parameters
labelstr

Column label.

datanumpy.ndarray

Modified data.

infostr, optional

Human-readable information about the modified 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_oldstr or iterable of str

Old column label(s).

label_newstr or iterable of str

New column label(s).

Returns
updateddamask.Table

Updated table.

sort_by(labels, ascending=True)[source]

Sort table by values of given labels.

Parameters
labelstr or list

Column labels for sorting.

ascendingbool or list, optional

Set sort order.

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.

join(other)[source]

Append other table horizontally (similar to numpy.hstack).

Requires matching number of rows and no common labels.

Parameters
otherdamask.Table

Table to join.

Returns
updateddamask.Table

Updated table.

save(fname)[source]

Save as plain text file.

Parameters
fnamefile, str, or pathlib.Path

Filename or file for writing.


VTK

class damask.VTK(vtk_data)[source]

Spatial visualization (and potentially manipulation).

High-level interface to VTK.

Methods

add(data[, label])

Add data to either cells or points.

add_comments(comments)

Add comments.

from_image_data(cells, size[, origin])

Create VTK of type vtk.vtkImageData.

from_poly_data(points)

Create VTK of type vtk.polyData.

from_rectilinear_grid(grid, size[, origin])

Create VTK of type vtk.vtkRectilinearGrid.

from_unstructured_grid(nodes, connectivity, …)

Create VTK of type vtk.vtkUnstructuredGrid.

get(label)

Get either cell or point data.

get_comments()

Return the comments.

load(fname[, dataset_type])

Load from VTK file.

save(fname[, parallel, compress])

Save as VTK file.

set_comments(comments)

Set comments.

show()

Render.

static from_image_data(cells, size, origin=array([0., 0., 0.]))[source]

Create VTK of type vtk.vtkImageData.

This is the common type for grid solver results.

Parameters
cellsiterable of int, len (3)

Number of cells along each dimension.

sizeiterable of float, len (3)

Physical length along each dimension.

originiterable of float, len (3), optional

Coordinates of grid origin.

Returns
newdamask.VTK

VTK-based geometry without nodal or cell data.

static from_rectilinear_grid(grid, size, origin=array([0., 0., 0.]))[source]

Create VTK of type vtk.vtkRectilinearGrid.

Parameters
griditerable of int, len (3)

Number of cells along each dimension.

sizeiterable of float, len (3)

Physical length along each dimension.

originiterable 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 vtk.vtkUnstructuredGrid.

This is the common type for mesh solver results.

Parameters
nodesnumpy.ndarray of shape (:,3)

Spatial position of the nodes.

connectivitynumpy.ndarray of np.dtype = int

Cell connectivity (0-based), first dimension determines #Cells, second dimension determines #Nodes/Cell.

cell_typestr

Name of the vtk.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 vtk.polyData.

This is the common type for point-wise data.

Parameters
pointsnumpy.ndarray of shape (:,3)

Spatial position of the points.

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 for reading. Valid extensions are .vti, .vtr, .vtu, .vtp, and .vtk.

dataset_type{‘vtkImageData’, ‘’vtkRectilinearGrid’, ‘vtkUnstructuredGrid’, ‘vtkPolyData’}, optional

Name of the vtk.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 for writing.

parallelboolean, optional

Write data in parallel background process. Defaults to True.

compressbool, optional

Compress with zlib algorithm. Defaults to True.

add(data, label=None)[source]

Add data to either cells or points.

Parameters
datanumpy.ndarray or numpy.ma.MaskedArray

Data to add. First dimension needs to match either number of cells or number of points.

labelstr

Data label.

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.

get_comments()[source]

Return the comments.

set_comments(comments)[source]

Set comments.

Parameters
commentsstr or list of str

Comments.

add_comments(comments)[source]

Add comments.

Parameters
commentsstr or list of str

Comments to add.

show()[source]

Render.

See http://compilatrix.com/article/vtk-1 for further ideas.


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 vtr 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
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

fnumpy.ndarray of shape (:,:,:,3) or (:,:,:,3,3)

Periodic field of which the curl is calculated.

Returns
∇ × fnumpy.ndarray

Curl of f.

damask.grid_filters.divergence(size, f)[source]

Calculate divergence of a vector or tensor field in Fourier space.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

fnumpy.ndarray of shape (:,:,:,3) or (:,:,:,3,3)

Periodic field of which the divergence is calculated.

Returns
∇ · fnumpy.ndarray

Divergence of f.

damask.grid_filters.gradient(size, f)[source]

Calculate gradient of a scalar or vector fieldin Fourier space.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

fnumpy.ndarray of shape (:,:,:,1) or (:,:,:,3)

Periodic field of which the gradient is calculated.

Returns
∇ fnumpy.ndarray

Divergence of f.

damask.grid_filters.coordinates0_point(cells, size, origin=array([0., 0., 0.]))[source]

Cell center positions (undeformed).

Parameters
cellsnumpy.ndarray of shape (3)

Number of cells.

sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

originnumpy.ndarray, optional

Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].

Returns
x_p_0numpy.ndarray

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
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

Returns
u_p_fluctnumpy.ndarray

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
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

Returns
u_p_avgnumpy.ndarray

Average part of the cell center displacements.

damask.grid_filters.displacement_point(size, F)[source]

Cell center displacement field from deformation gradient field.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

Returns
u_pnumpy.ndarray

Cell center displacements.

damask.grid_filters.coordinates_point(size, F, origin=array([0., 0., 0.]))[source]

Cell center positions.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

originnumpy.ndarray of shape (3), optional

Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].

Returns
x_pnumpy.ndarray

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 of shape (:,3)

Undeformed cell 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
cellsnumpy.ndarray of shape (3)

Number of cells.

sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

originnumpy.ndarray of shape (3), optional

Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].

Returns
x_n_0numpy.ndarray

Undeformed nodal coordinates.

damask.grid_filters.displacement_fluct_node(size, F)[source]

Nodal displacement field from fluctuation part of the deformation gradient field.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

Returns
u_n_fluctnumpy.ndarray

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
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

Returns
u_n_avgnumpy.ndarray

Average part of the nodal displacements.

damask.grid_filters.displacement_node(size, F)[source]

Nodal displacement field from deformation gradient field.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

Returns
u_pnumpy.ndarray

Nodal displacements.

damask.grid_filters.coordinates_node(size, F, origin=array([0., 0., 0.]))[source]

Nodal positions.

Parameters
sizenumpy.ndarray of shape (3)

Physical size of the periodic field.

Fnumpy.ndarray

Deformation gradient field.

originnumpy.ndarray of shape (3), optional

Physical origin of the periodic field. Defaults to [0.0,0.0,0.0].

Returns
x_nnumpy.ndarray

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 of 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 of shape (:,:,:,…)

Data defined on the cell centers of a periodic grid.

Returns
node_datanumpy.ndarray of 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 of shape (:,:,:,…)

Data defined on the nodes of a periodic grid.

Returns
cell_datanumpy.ndarray of 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

Array of undeformed cell coordinates.

Returns
validbool

Whether the coordinates form a regular grid.

damask.grid_filters.regrid(size, F, cells)[source]

Return mapping from coordinates in deformed configuration to a regular grid.

Parameters
sizenumpy.ndarray of shape (3)

Physical size.

Fnumpy.ndarray of shape (:,:,:,3,3)

Deformation gradient field.

cellsnumpy.ndarray of shape (3)

Cell count along x,y,z of remapping grid.


tensor

Tensor mathematics.

All routines operate on numpy.ndarrays of shape (…,3,3).

damask.tensor.deviatoric(T)[source]

Calculate deviatoric part of a tensor.

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

Tensor of which the deviatoric part is computed.

Returns
T’numpy.ndarray of shape (…,3,3)

Deviatoric part of T.

damask.tensor.eigenvalues(T_sym)[source]

Eigenvalues, i.e. principal components, of a symmetric tensor.

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

Symmetric tensor of which the eigenvalues are computed.

Returns
lambdanumpy.ndarray of shape (…,3)

Eigenvalues of T_sym sorted in ascending order, each repeated according to its multiplicity.

damask.tensor.eigenvectors(T_sym, RHS=False)[source]

Eigenvectors of a symmetric tensor.

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

Symmetric tensor of which the eigenvectors are computed.

RHS: bool, optional

Enforce right-handed coordinate system. Defaults to False.

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

Eigenvectors of T_sym sorted in ascending order of their associated eigenvalues.

damask.tensor.spherical(T, tensor=True)[source]

Calculate spherical part of a tensor.

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

Tensor of which the spherical part is computed.

tensorbool, optional

Map spherical part onto identity tensor. Defaults to True.

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

unless tensor == False: shape (…,) Spherical part of tensor T. p is an isotropic tensor.

damask.tensor.symmetric(T)[source]

Symmetrize tensor.

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

Tensor of which the symmetrized values are computed.

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

Symmetrized tensor T.

damask.tensor.transpose(T)[source]

Transpose tensor.

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

Tensor of which the transpose is computed.

Returns
T.Tnumpy.ndarray of shape (…,3,3)

Transpose of tensor T.


mechanics

Finite-strain continuum mechanics.

All routines operate on numpy.ndarrays of shape (…,3,3).

damask.mechanics.deformation_Cauchy_Green_left(F)[source]

Calculate left Cauchy-Green deformation tensor (Finger deformation tensor).

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

Deformation gradient.

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

Left Cauchy-Green deformation tensor.

damask.mechanics.deformation_Cauchy_Green_right(F)[source]

Calculate right Cauchy-Green deformation tensor.

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

Deformation gradient.

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

Right Cauchy-Green deformation tensor.

damask.mechanics.equivalent_strain_Mises(epsilon)[source]

Calculate the Mises equivalent of a strain tensor.

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

Symmetric strain tensor of which the von Mises equivalent is computed.

Returns
epsilon_vMnumpy.ndarray of shape (…)

Von Mises equivalent strain of epsilon.

damask.mechanics.equivalent_stress_Mises(sigma)[source]

Calculate the Mises equivalent of a stress tensor.

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

Symmetric stress tensor of which the von Mises equivalent is computed.

Returns
sigma_vMnumpy.ndarray of shape (…)

Von Mises equivalent stress of sigma.

damask.mechanics.maximum_shear(T_sym)[source]

Calculate the maximum shear component of a symmetric tensor.

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

Symmetric tensor of which the maximum shear is computed.

Returns
gamma_maxnumpy.ndarray of shape (…)

Maximum shear of T_sym.

damask.mechanics.rotation(T)[source]

Calculate the rotational part of a tensor.

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

Tensor of which the rotational part is computed.

Returns
Rdamask.Rotation of shape (…)

Rotational part of the vector.

damask.mechanics.strain(F, t, m)[source]

Calculate strain tensor (Seth–Hill family).

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

Deformation gradient.

t{‘V’, ‘U’}

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

mfloat

Order of the strain.

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

Strain of F.

References

https://en.wikipedia.org/wiki/Finite_strain_theory https://de.wikipedia.org/wiki/Verzerrungstensor

damask.mechanics.stress_Cauchy(P, F)[source]

Calculate the Cauchy stress (true stress).

Resulting tensor is symmetrized as the Cauchy stress needs to be symmetric.

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

First Piola-Kirchhoff stress.

Fnumpy.ndarray of shape (…,3,3)

Deformation gradient.

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

Cauchy stress.

damask.mechanics.stress_second_Piola_Kirchhoff(P, F)[source]

Calculate the second Piola-Kirchhoff stress.

Resulting tensor is symmetrized as the second Piola-Kirchhoff stress needs to be symmetric.

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

First Piola-Kirchhoff stress.

Fnumpy.ndarray of shape (…,3,3)

Deformation gradient.

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

Second Piola-Kirchhoff stress.

damask.mechanics.stretch_left(T)[source]

Calculate left stretch of a tensor.

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

Tensor of which the left stretch is computed.

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

Left stretch tensor from Polar decomposition of T.

damask.mechanics.stretch_right(T)[source]

Calculate right stretch of a tensor.

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

Tensor of which the right stretch is computed.

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

Left stretch tensor from Polar decomposition of T.


solver.Marc

class damask.solver.Marc(msc_version=2020, msc_root='/opt/msc', damask_root='/home/m/DAMASK')[source]

Wrapper to run DAMASK with MSCMarc.

Attributes
library_path
tools_path

Methods

submit_job