tensormesh.functional

ops

sym(a: Tensor) Tensor[source]
\[\text{sym}(A)_{\cdots ij} = \frac{1}{2} (A_{\cdots i} + A_{\cdots j})\]

Examples

>>> x = torch.tensor([1., 2.])
>>> sym(x)
tensor([[1.0000, 1.5000],
        [1.5000, 2.0000]])
>>> x = torch.tensor([1., 2., 3.])
>>> sym(x)
tensor([[1.0000, 1.5000, 2.0000],
        [1.5000, 2.0000, 2.5000],
        [2.0000, 2.5000, 3.0000]])
Parameters:

a (Tensor) – \([..., D]\), where \(D\) is the dimension of the matrix

Returns:

\([..., D, D]\), where \(D\) is the dimension of the matrix

Return type:

Tensor

skew(x: Tensor, sign: bool = True, at_least2d: bool = False) Tensor[source]

Compute the skew-symmetric matrix from a vector.

For 2D:

\[\begin{split}\text{skew}\left(\begin{bmatrix} v_1 \\ v_2 \end{bmatrix}\right) = \begin{cases} \begin{bmatrix} -v_2 \\ v_1 \end{bmatrix} & \text{if sign=True} \\[1em] \begin{bmatrix} v_2 \\ v_1 \end{bmatrix} & \text{if sign=False} \end{cases}\end{split}\]

For 3D:

\[\begin{split}\text{skew}(v) = \begin{cases} \begin{bmatrix} 0 & -v_3 & v_2 \\ v_3 & 0 & -v_1 \\ -v_2 & v_1 & 0 \end{bmatrix} & \text{if sign=True} \\[1em] \begin{bmatrix} 0 & v_3 & v_2 \\ v_3 & 0 & v_1 \\ v_2 & v_1 & 0 \end{bmatrix} & \text{if sign=False} \end{cases}\end{split}\]

Examples

>>> x = torch.tensor([1., 2.])
>>> skew(x)
tensor([-2.,  1.])

>>> skew(x, sign=False)
tensor([2., 1.])

>>> skew(x, at_least2d=True)
tensor([[-2.,  1.]])

>>> x = torch.tensor([1., 2., 3.])
>>> skew(x)
tensor([[ 0., -3.,  2.],
       [ 3.,  0., -1.],
       [-2.,  1.,  0.]])

>>> skew(x, sign=False)
tensor([[0., 3., 2.],
       [3., 0., 1.],
       [2., 1., 0.]])
Parameters:
  • x (Tensor) – 1D Tensor of shape [2] or [3], representing a vector in \(\mathbb{R}^2\) or \(\mathbb{R}^3\)

  • sign (bool, optional) – If True, use negative signs in skew matrix. Default is True.

  • at_least2d (bool, optional) – If True, ensure output is at least 2D for 2D case. Default is False.

Returns:

For 2D case:
  • 1D Tensor of shape [2] if at_least2d=False

  • 2D Tensor of shape [1,2] if at_least2d=True

For 3D case:
  • 2D Tensor of shape [3,3]

The skew-symmetric matrix representation of the input vector.

Return type:

Tensor

sqrt(x: Tensor) Tensor[source]

Square root function that returns 0 for negative inputs.

This function computes the square root of the input tensor, but clamps negative values to 0 first. This avoids NaN values that would occur from taking the square root of negative numbers.

\[\begin{split}\sqrt{x} = \begin{cases} \sqrt{x} & \text{if } x \geq 0 \\ 0 & \text{if } x < 0 \end{cases}\end{split}\]

Examples

>>> x = torch.tensor([-1.0, 0.0, 4.0])
>>> sqrt(x)
tensor([0.0000, 0.0000, 2.0000])
Parameters:

x (Tensor) – \([...]\)

Returns:

\([...]\)

Return type:

Tensor

divide(x: Tensor, y: Tensor) Tensor[source]

Safe division function that returns 0 for division by zero.

This function performs element-wise division of x by y, but returns 0 wherever y is 0. This avoids NaN/Inf values that would occur from dividing by zero.

\[\begin{split}\frac{x}{y} = \begin{cases} \frac{x}{y} & \text{if } y \neq 0 \\ 0 & \text{if } y = 0 \end{cases}\end{split}\]

Examples

>>> x = torch.tensor([1.0, 2.0, 3.0])
>>> y = torch.tensor([2.0, 0.0, 4.0])
>>> divide(x, y)
tensor([0.5000, 0.0000, 0.7500])
Parameters:
Returns:

\([...]\)

Return type:

Tensor

elasticity

strain(gradu: Tensor) Tensor[source]
\[\varepsilon_{ij} = \frac{1}{2}(\nabla u_{ij} + \nabla u_{ji}), \quad \varepsilon,\nabla u \in \mathbb{R}^{d \times d}\]

where:

  • \(\nabla u \in \mathbb{R}^{d \times d}\) is the displacement gradient tensor

  • \(\varepsilon \in \mathbb{R}^{d \times d}\) is the strain tensor

Parameters:

gradu (Tensor) – 1D Tensor of shape [d] gradient of displacement field u with respect to spatial coordinates

Returns:

2D Tensor of shape [d, d]

Return type:

Tensor

isotropic_stress(strain: Tensor, E: float | Tensor = 70.0, nu: float | Tensor = 0.3) Tensor[source]
\[\sigma_{ij} = \frac{E}{1+\nu} \left(\varepsilon_{ij} + \frac{\nu}{1-2\nu} \text{tr}(\varepsilon) \delta_{ij}\right), \quad \sigma,\varepsilon \in \mathbb{R}^{d \times d}\]

where:

  • \(\sigma \in \mathbb{R}^{d \times d}\) is the stress tensor

  • \(\varepsilon \in \mathbb{R}^{d \times d}\) is the strain tensor

  • \(E \in \mathbb{R}\) is Young’s modulus

  • \(\nu \in \mathbb{R}\) is Poisson’s ratio

  • \(\delta_{ij} \in \mathbb{R}^{d \times d}\) is the Kronecker delta tensor

  • \(\text{tr}: \mathbb{R}^{d \times d} \rightarrow \mathbb{R}\) denotes the trace operator

Parameters:
  • strain (Tensor) – 2D Tensor of shape [d, d] strain tensor

  • E (Union[float,Tensor]) – if torch.Tensor, 0D Tensor of shape [] Young’s modulus

  • nu (Union[float,Tensor]) – if torch.Tensor, 0D Tensor of shape [] Poisson’s ratio

Returns:

2D Tensor of shape [d, d]

Return type:

Tensor

deviatoric_stress(stress: Tensor) Tensor[source]
\[s_{ij} = \sigma_{ij} - \frac{1}{d} \text{tr}(\sigma) \delta_{ij}, \quad s,\sigma \in \mathbb{R}^{d \times d}\]

where:

  • \(s_{ij} \in \mathbb{R}^{d \times d}\) is the deviatoric stress tensor

  • \(\sigma_{ij} \in \mathbb{R}^{d \times d}\) is the stress tensor

  • \(d \in \mathbb{N}\) is the dimension

  • \(\delta_{ij} \in \mathbb{R}^{d \times d}\) is the Kronecker delta

  • \(\text{tr}: \mathbb{R}^{d \times d} \rightarrow \mathbb{R}\) denotes the trace operator

Parameters:

stress (Tensor) – 2D Tensor of shape [d, d] stress tensor

Returns:

2D Tensor of shape [d, d]

Return type:

Tensor

deviatoric_stress_norm(stress: Tensor) Tensor[source]

Euclidean norm of the deviatoric stress tensor.

\[\|s\| = \sqrt{\frac{3}{2} s:s}\]

where \(s\) is the deviatoric stress tensor, \(\|s\|\) is its norm, and : denotes the double-dot product.

Parameters:

stress (Tensor) – 2D Tensor of shape [d, d] stress tensor

Returns:

0D Tensor of shape []

Return type:

Tensor

voigt_shape_grad(gradu: Tensor) Tensor[source]

Convert displacement gradient to Voigt notation for strain calculation.

For 2D:

\[\begin{split}B \in \mathbb{R}^{3 \times (2\times 3)} = \begin{bmatrix} \frac{\partial u_1}{\partial x_1} & 0 & \mid & \frac{\partial u_2}{\partial x_1} & 0 & \mid & \frac{\partial u_3}{\partial x_1} & 0 \\ 0 & \frac{\partial u_1}{\partial x_2} & \mid & 0 & \frac{\partial u_2}{\partial x_2} & \mid & 0 & \frac{\partial u_3}{\partial x_2} \\ \frac{\partial u_1}{\partial x_2} & \frac{\partial u_1}{\partial x_1} & \mid & \frac{\partial u_2}{\partial x_2} & \frac{\partial u_2}{\partial x_1} & \mid & \frac{\partial u_3}{\partial x_2} & \frac{\partial u_3}{\partial x_1} \end{bmatrix}\end{split}\]

representing [\(\varepsilon_{xx}, \varepsilon_{yy}, \gamma_{xy}\)]

where:

  • \(u_i \in \mathbb{R}\) is the displacement component in direction i

  • \(x_i \in \mathbb{R}\) is the spatial coordinate in direction i

  • \(\frac{\partial u_i}{\partial x_j} \in \mathbb{R}\) is the partial derivative of displacement i with respect to coordinate j

For 3D:

\[\begin{split}B \in \mathbb{R}^{6 \times (3\times 3)} = \begin{bmatrix} \frac{\partial u_1}{\partial x_1} & 0 & 0 & \mid & \frac{\partial u_2}{\partial x_1} & 0 & 0 & \mid & \frac{\partial u_3}{\partial x_1} & 0 & 0 \\ 0 & \frac{\partial u_1}{\partial x_2} & 0 & \mid & 0 & \frac{\partial u_2}{\partial x_2} & 0 & \mid & 0 & \frac{\partial u_3}{\partial x_2} & 0 \\ 0 & 0 & \frac{\partial u_1}{\partial x_3} & \mid & 0 & 0 & \frac{\partial u_2}{\partial x_3} & \mid & 0 & 0 & \frac{\partial u_3}{\partial x_3} \\ 0 & \frac{\partial u_1}{\partial x_3} & \frac{\partial u_1}{\partial x_2} & \mid & 0 & \frac{\partial u_2}{\partial x_3} & \frac{\partial u_2}{\partial x_2} & \mid & 0 & \frac{\partial u_3}{\partial x_3} & \frac{\partial u_3}{\partial x_2} \\ \frac{\partial u_1}{\partial x_3} & 0 & \frac{\partial u_1}{\partial x_1} & \mid & \frac{\partial u_2}{\partial x_3} & 0 & \frac{\partial u_2}{\partial x_1} & \mid & \frac{\partial u_3}{\partial x_3} & 0 & \frac{\partial u_3}{\partial x_1} \\ \frac{\partial u_1}{\partial x_2} & \frac{\partial u_1}{\partial x_1} & 0 & \mid & \frac{\partial u_2}{\partial x_2} & \frac{\partial u_2}{\partial x_1} & 0 & \mid & \frac{\partial u_3}{\partial x_2} & \frac{\partial u_3}{\partial x_1} & 0 \end{bmatrix}\end{split}\]

representing [\(\varepsilon_{xx}, \varepsilon_{yy}, \varepsilon_{zz}, \gamma_{yz}, \gamma_{xz}, \gamma_{xy}\)]

where \(\varepsilon\) denotes normal strain and \(\gamma\) denotes shear strain components

where:

  • \(u_i \in \mathbb{R}\) is the displacement component in direction i

  • \(x_i \in \mathbb{R}\) is the spatial coordinate in direction i

  • \(\frac{\partial u_i}{\partial x_j} \in \mathbb{R}\) is the partial derivative of displacement i with respect to coordinate j

Parameters:

gradu (Tensor) – 1D Tensor of shape [d], where dim must be 2 or 3 gradient of displacement

Returns:

2D Tensor with shape [3, 2] if d=2 or [6, 3] if d=3,

Return type:

Tensor

voigt_stiffness(E: float | Tensor, nu: float | Tensor, dim: int = 2) Tensor[source]

For 2D:

\[\begin{split}\mathbb{C} \in \mathbb{R}^{3 \times 3} = \begin{bmatrix} \lambda + 2\mu & \lambda & 0 \\ \lambda & \lambda + 2\mu & 0 \\ 0 & 0 & \mu \end{bmatrix}\end{split}\]

representing [\(\varepsilon_{xx}, \varepsilon_{yy}, \gamma_{xy}\)]

For 3D:

\[\begin{split}\mathbb{C} \in \mathbb{R}^{6 \times 6} = \begin{bmatrix} \lambda + 2\mu & \lambda & \lambda & 0 & 0 & 0 \\ \lambda & \lambda + 2\mu & \lambda & 0 & 0 & 0 \\ \lambda & \lambda & \lambda + 2\mu & 0 & 0 & 0 \\ 0 & 0 & 0 & \mu & 0 & 0 \\ 0 & 0 & 0 & 0 & \mu & 0 \\ 0 & 0 & 0 & 0 & 0 & \mu \end{bmatrix}\end{split}\]

representing [\(\varepsilon_{xx}, \varepsilon_{yy}, \varepsilon_{zz}, \gamma_{yz}, \gamma_{xz}, \gamma_{xy}\)]

where:

  • \(\lambda \in \mathbb{R} = \frac{E\nu}{(1+\nu)(1-2\nu)}\)

  • \(\mu \in \mathbb{R} = \frac{E}{2(1+\nu)}\)

Parameters:
  • E (Union[float,Tensor]) – if torch.Tensor, 0D Tensor Young’s modulus

  • nu (float) – if torch.Tensor, 0D Tensor Poisson’s ratio

  • dim (int) – Spatial dimension of the problem (2 or 3). For 2D problems use dim=2, for 3D problems use dim=3. Default is 2.

Returns:

2D Tensor of shape \([d(d+1)/2, d(d+1)/2]\), where \(d\) is the spatial dimension

Return type:

Tensor

voigt_shape_val(u: Tensor, dim: int) Tensor[source]

Convert shape functions to Voigt notation matrix for strain-displacement relations.

For 2D:

\[\begin{split}N_1, N_2, N_3 \in \mathbb{R} \quad \text{(shape values)} \\ N \in \mathbb{R}^{2 \times (2\times 3)}= \begin{bmatrix} N_1 & 0 &\mid& N_2 & 0 &\mid& N_3 & 0 \\ 0 & N_1 &\mid& 0 & N_2 &\mid& 0 & N_3 \end{bmatrix}\end{split}\]

For 3D:

\[\begin{split}N_1, N_2, N_3 \in \mathbb{R} \quad \text{(shape values)} \\ N \in \mathbb{R}^{3 \times (3\times 3)}= \begin{bmatrix} N_1 & 0 & 0 &\mid& N_2 & 0 & 0 &\mid& N_3 & 0 & 0 \\ 0 & N_1 & 0 &\mid& 0 & N_2 & 0 &\mid& 0 & N_3 & 0 \\ 0 & 0 & N_1 &\mid& 0 & 0 & N_2 &\mid& 0 & 0 & N_3 \end{bmatrix}\end{split}\]
Parameters:
  • u (Tensor) – 0D Tensor of shape [], shape value

  • dim (int) – Spatial dimension. Must be 2 or 3.

Returns:

2D Tensor of shape [d, d], where \(d\in \{2,3\}\)

Return type:

Tensor

plastic

update_plastic_stress(gradu: Tensor, strain: Tensor, stress: Tensor, E: float | Tensor = 70.0, yield_stress: float | Tensor = 250.0, strain_fn: Callable[[Tensor], Tensor] = strain, stress_fn: Callable[[Tensor, Tensor | float], Tensor] = isotropic_stress) Tensor[source]

Update stress tensor using plastic constitutive model.

The plastic model follows von Mises yield criterion with perfect plasticity:

\[ \begin{align}\begin{aligned}\sigma_{\text{trial}} = \sigma + \mathbb{C}:\Delta\varepsilon\\f(\sigma_{\text{trial}}) = \|\text{dev}(\sigma_{\text{trial}})\| - \sigma_y\\\Delta\gamma = \frac{\langle f(\sigma_{\text{trial}}) \rangle}{\|\text{dev}(\sigma_{\text{trial}})\|}\\\sigma = \sigma_{\text{trial}} - \Delta\gamma\, \text{dev}(\sigma_{\text{trial}})\end{aligned}\end{align} \]

where:

  • \(\sigma\) is the stress tensor in \(\mathbb{R}^{D \times D}\)

  • \(\mathbb{C}\) is the elasticity tensor in \(\mathbb{R}^{D \times D \times D \times D}\)

  • \(\varepsilon\) is the strain tensor in \(\mathbb{R}^{D \times D}\)

  • \(\sigma_y\) is the yield stress scalar in \(\mathbb{R}\)

  • \(\text{dev}\) denotes the deviatoric part operator \(\mathbb{R}^{D \times D} \rightarrow \mathbb{R}^{D \times D}\)

  • \(\|\cdot\|\) is the von Mises norm operator \(\mathbb{R}^{D \times D} \rightarrow \mathbb{R}\)

  • \(\langle \cdot \rangle\) denotes the positive part operator \(\mathbb{R} \rightarrow \mathbb{R}\)

The model uses a trial elastic predictor followed by plastic correction if yielding occurs. If the trial stress exceeds the yield surface, it is projected back onto the yield surface.

Parameters:
  • gradu (Tensor) – 1D Tensor of shape [d], where d is the spatial dimension. Gradient of displacement field with respect to spatial coordinates.

  • strain (Tensor) – 2D Tensor of shape [d, d], where d is the spatial dimension. Current strain tensor at the start of the timestep.

  • stress (Tensor) – 2D Tensor of shape [d, d], where d is the spatial dimension. Current stress tensor at the start of the timestep.

  • E (Union[float, Tensor], default=70.0) – Young’s modulus. If tensor, must be 0D scalar tensor. Controls the elastic stiffness of the material.

  • yield_stress (Union[float, Tensor], default=250.0) – Yield stress threshold. If tensor, must be 0D scalar tensor. Material yields plastically when von Mises stress exceeds this value.

  • strain_fn (Callable[[Tensor], Tensor], default=strain) –

    Function to compute strain tensor from displacement gradient. Default uses small strain assumption:

    \[\varepsilon_{ij} = \frac{1}{2}(\nabla u_{ij} + \nabla u_{ji}), \quad \varepsilon,\nabla u \in \mathbb{R}^{d \times d}\]

  • stress_fn (Callable[[Tensor, Union[float,Tensor]], Tensor], default=isotropic_stress) –

    Function to compute stress tensor from strain tensor and Young’s modulus. Default uses isotropic linear elasticity:

    \[\sigma_{ij} = \lambda \text{tr}(\varepsilon)\delta_{ij} + 2\mu\varepsilon_{ij}, \quad \sigma,\varepsilon \in \mathbb{R}^{d \times d}\]

    where \(\lambda = \frac{E\nu}{(1+\nu)(1-2\nu)}\), \(\mu = \frac{E}{2(1+\nu)}\), \(E,\nu \in \mathbb{R}\), and \(\delta_{ij}\) is the Kronecker delta

Returns:

2D Tensor of shape [d, d], where d is the spatial dimension. Updated stress tensor after plastic correction.

Return type:

Tensor