Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions CHANGELOG.md
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
`.to_dual_graph()` methods. These allow Mesh conversion to 0D point clouds, 1D
edge graphs, and 1D dual graphs, respectively, when connectivity information
is not needed.
- Adds von Mises stress (`physicsnemo.metrics.cae.von_mises`).
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
- Adds von Mises stress (`physicsnemo.metrics.cae.von_mises`).
- Adds von Mises stress (`physicsnemo.metrics.cae.von_mises`) as a CAE metric.


### Changed

Expand Down
121 changes: 121 additions & 0 deletions physicsnemo/metrics/cae/von_mises.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,121 @@
# SPDX-FileCopyrightText: Copyright (c) 2023 - 2026 NVIDIA CORPORATION & AFFILIATES.
# SPDX-FileCopyrightText: All rights reserved.
# SPDX-License-Identifier: Apache-2.0
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

"""Von Mises stress metric for structural mechanics applications.

Von Mises stress is a scalar equivalent stress derived from the
stress tensor and used to predict yielding in ductile materials
under multiaxial loading, based on the distortional energy (J₂)
criterion.
"""

from typing import Union

import numpy as np
import torch


def von_mises_from_stress_tensor(
stress: Union[torch.Tensor, np.ndarray],
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

In new code, we should always prefer modern | union types over Union. (Both here and throughout)

) -> Union[torch.Tensor, np.ndarray]:
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

These should all use jaxtyping to be compliant with repo rules. (Both here and throughout, on all type signatures for tensors)

"""Compute von Mises (equivalent) stress from a symmetric stress tensor.

For a 3D symmetric stress tensor σ with components σxx, σyy, σzz, τxy, τyz, τxz,
the von Mises stress is:

σ_vm = sqrt(0.5 * ((σxx-σyy)² + (σyy-σzz)² + (σzz-σxx)²) + 3*(τxy² + τyz² + τxz²))

This is equivalent to sqrt(3*J2) where J2 is the second invariant of the
deviatoric stress tensor. Shear components are symmetrized (e.g. τxy = 0.5*(σxy+σyx))
for robustness to slight numerical asymmetry.

Parameters
----------
stress : Union[Tensor, np.ndarray]
Symmetric stress tensor. Supported shapes:
- (3, 3): single 3x3 tensor
- (N, 3, 3): batch of N 3x3 tensors
Components are ordered as [σxx, σxy, σxz; σxy, σyy, σyz; σxz, σyz, σzz].

Returns
-------
Union[Tensor, np.ndarray]
Von Mises stress. Scalar for (3,3) input, shape (N,) for (N,3,3) input.
Same type (torch.Tensor or np.ndarray) as input.

Raises
------
ValueError
If stress has invalid shape (not 3x3 or Nx3x3).

Examples
--------
>>> import torch
>>> from physicsnemo.metrics.cae.von_mises import von_mises_from_stress_tensor
>>> # Uniaxial tension: σxx = 100 MPa, all others zero
>>> stress = torch.tensor([[100., 0., 0.], [0., 0., 0.], [0., 0., 0.]])
>>> von_mises_from_stress_tensor(stress)
tensor(100.)
>>> # Hydrostatic stress (no distortion): σxx = σyy = σzz = 50
>>> stress = torch.eye(3) * 50
>>> von_mises_from_stress_tensor(stress)
tensor(0.)
"""
if isinstance(stress, torch.Tensor):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It might be worth using the Array API standard to deduplicate the parallel numpy and torch code here? This is the official effort coordinated across NumPy/Torch/JAX/others to allow for multiple dispatch on arrays.

https://data-apis.org/array-api-compat/

I've used this in other projects - it's excellent and makes for very clean reusable implementations.

_validate_stress_shape(stress.shape)
sxx = stress[..., 0, 0]
syy = stress[..., 1, 1]
szz = stress[..., 2, 2]
# Symmetrize shear components for robustness to slight asymmetry
sxy = 0.5 * (stress[..., 0, 1] + stress[..., 1, 0])
syz = 0.5 * (stress[..., 1, 2] + stress[..., 2, 1])
sxz = 0.5 * (stress[..., 0, 2] + stress[..., 2, 0])
vm_sq = 0.5 * ((sxx - syy) ** 2 + (syy - szz) ** 2 + (szz - sxx) ** 2) + 3.0 * (
sxy**2 + syz**2 + sxz**2
)
return torch.sqrt(torch.clamp(vm_sq, min=0.0))

_validate_stress_shape(np.shape(stress))
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

For readability, and for better error handling, I'd recommend always prefixing switch-statement clauses with an:

elif isinstance(stress, np.ndarray):

And then having an:

else:
    raise TypeError("...")

In other words, switch statements should always declare their switches.

stress = np.asarray(stress, dtype=np.float64)
sxx = stress[..., 0, 0]
syy = stress[..., 1, 1]
szz = stress[..., 2, 2]
sxy = 0.5 * (stress[..., 0, 1] + stress[..., 1, 0])
syz = 0.5 * (stress[..., 1, 2] + stress[..., 2, 1])
sxz = 0.5 * (stress[..., 0, 2] + stress[..., 2, 0])
vm_sq = 0.5 * ((sxx - syy) ** 2 + (syy - szz) ** 2 + (szz - sxx) ** 2) + 3.0 * (
sxy**2 + syz**2 + sxz**2
)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This can be vectorized, if desired, using einsum. This may yield performance benefits if this is used downstream in performance-critical applications - up to you.

It's also possible to write this code in a way that's dimensionally-generic (e.g., works on 2x2 stress tensors for 2D problems), which could be nice down the road for reusability?

return np.sqrt(np.maximum(vm_sq, 0.0))


def _validate_stress_shape(shape):
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Missing type signature here

"""Raise ValueError if stress shape is invalid."""
ndim = len(shape)
if ndim == 2:
if shape != (3, 3):
raise ValueError(
f"Expected stress shape (3, 3) for single tensor, got {shape}"
)
elif ndim == 3:
if shape[1:] != (3, 3):
raise ValueError(
f"Expected stress shape (N, 3, 3) for batched tensor, got {shape}"
)
else:
raise ValueError(
f"Expected stress with 2 or 3 dimensions, got {ndim} (shape {shape})"
)
135 changes: 135 additions & 0 deletions test/metrics/test_metrics_von_mises.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,135 @@
# SPDX-FileCopyrightText: Copyright (c) 2023 - 2026 NVIDIA CORPORATION & AFFILIATES.
# SPDX-FileCopyrightText: All rights reserved.
# SPDX-License-Identifier: Apache-2.0
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.

"""Tests for von Mises stress metrics."""

import numpy as np
import pytest
import torch

from physicsnemo.metrics.cae.von_mises import von_mises_from_stress_tensor


def test_von_mises_uniaxial_tension_torch():
"""Uniaxial tension: σxx = 100, von Mises = 100."""
stress = torch.tensor([[100.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]])
vm = von_mises_from_stress_tensor(stress)
assert vm.shape == ()
assert torch.allclose(vm, torch.tensor(100.0))


def test_von_mises_uniaxial_tension_numpy():
"""Uniaxial tension with numpy input."""
stress = np.array([[100.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]])
vm = von_mises_from_stress_tensor(stress)
assert vm.shape == ()
assert np.allclose(vm, 100.0)


def test_von_mises_hydrostatic_torch():
"""Hydrostatic stress: σxx = σyy = σzz = p, von Mises = 0 (no distortion)."""
p = 50.0
stress = torch.eye(3) * p
vm = von_mises_from_stress_tensor(stress)
assert torch.allclose(vm, torch.tensor(0.0))


def test_von_mises_pure_shear_torch():
"""Pure shear: τxy = τ, von Mises = sqrt(3)*τ."""
tau = 50.0
stress = torch.tensor(
[[0.0, tau, 0.0], [tau, 0.0, 0.0], [0.0, 0.0, 0.0]],
dtype=torch.float64,
)
vm = von_mises_from_stress_tensor(stress)
expected = np.sqrt(3) * tau
assert torch.allclose(vm, torch.tensor(expected))


def test_von_mises_batched_torch(device):
"""Batched stress tensors."""
# Two samples: uniaxial 100 and uniaxial 200
stress = torch.zeros(2, 3, 3, device=device)
stress[0, 0, 0] = 100.0
stress[1, 0, 0] = 200.0
vm = von_mises_from_stress_tensor(stress)
assert vm.shape == (2,)
assert torch.allclose(vm[0], torch.tensor(100.0, device=device))
assert torch.allclose(vm[1], torch.tensor(200.0, device=device))


def test_von_mises_batched_numpy():
"""Batched numpy input."""
stress = np.zeros((3, 3, 3))
stress[0, 0, 0] = 10.0
stress[1, 1, 1] = 20.0
stress[2, 2, 2] = 30.0
vm = von_mises_from_stress_tensor(stress)
assert vm.shape == (3,)
# For uniaxial: vm = principal stress
assert np.allclose(vm[0], 10.0)
assert np.allclose(vm[1], 20.0)
assert np.allclose(vm[2], 30.0)


def test_von_mises_invalid_shape_2d():
"""Invalid 2D shape raises ValueError."""
stress = torch.ones(2, 4)
with pytest.raises(ValueError, match="Expected stress shape"):
von_mises_from_stress_tensor(stress)


def test_von_mises_invalid_shape_3d():
"""Invalid 3D shape raises ValueError."""
stress = torch.ones(5, 2, 2)
with pytest.raises(ValueError, match="Expected stress shape"):
von_mises_from_stress_tensor(stress)


def test_von_mises_invalid_ndim():
"""Invalid number of dimensions raises ValueError."""
stress = torch.ones(2, 3, 3, 3)
with pytest.raises(ValueError, match="Expected stress with 2 or 3 dimensions"):
von_mises_from_stress_tensor(stress)


def test_von_mises_autograd():
"""Gradients flow through for torch tensors."""
stress = torch.tensor(
[[100.0, 10.0, 0.0], [10.0, 20.0, 0.0], [0.0, 0.0, 0.0]],
requires_grad=True,
)
vm = von_mises_from_stress_tensor(stress)
vm.backward()
assert stress.grad is not None
assert stress.grad.shape == stress.shape


def test_von_mises_asymmetric_symmetrized():
"""Asymmetric shear (σxy ≠ σyx) is symmetrized; result matches symmetric case."""
# Symmetric: σxy = σyx = 50
stress_sym = torch.tensor(
[[0.0, 50.0, 0.0], [50.0, 0.0, 0.0], [0.0, 0.0, 0.0]],
dtype=torch.float64,
)
vm_sym = von_mises_from_stress_tensor(stress_sym)
# Asymmetric: σxy=60, σyx=40 -> symmetrized τxy = 50, same as symmetric
stress_asym = torch.tensor(
[[0.0, 60.0, 0.0], [40.0, 0.0, 0.0], [0.0, 0.0, 0.0]],
dtype=torch.float64,
)
vm_asym = von_mises_from_stress_tensor(stress_asym)
assert torch.allclose(vm_sym, vm_asym)
Loading