diff --git a/CHANGELOG.md b/CHANGELOG.md index d25fbad153cb..17b07122794c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 * Added `compas_rhino.install_with_pip` with corresponding command line utility `install_in_rhino`. * Added support for `.stp` file extension in addition to `.step` for `RhinoBrep.from_step()` and `RhinoBrep.to_step()` methods. +* Added `volume()` method to `compas.datastructures.Mesh` for computing the volume of closed meshes using signed volume of triangles. ### Changed diff --git a/src/compas/datastructures/mesh/mesh.py b/src/compas/datastructures/mesh/mesh.py index 246023a48315..81e371780d70 100644 --- a/src/compas/datastructures/mesh/mesh.py +++ b/src/compas/datastructures/mesh/mesh.py @@ -42,6 +42,7 @@ from compas.geometry import distance_line_line from compas.geometry import distance_point_plane from compas.geometry import distance_point_point +from compas.geometry import dot_vectors from compas.geometry import length_vector from compas.geometry import midpoint_line from compas.geometry import normal_polygon @@ -3882,6 +3883,75 @@ def area(self): """ return sum(self.face_area(fkey) for fkey in self.faces()) + def volume(self, copy=True, unify_cycles=True): + """Calculate the volume of the mesh. + + Parameters + ---------- + copy : bool, optional + If True, a copy of the mesh is made before computation to avoid modifying the original. + Default is True. + unify_cycles : bool, optional + If True, face cycles are unified to ensure consistent orientation. + Default is True. + + Returns + ------- + float | None + The volume of the mesh if the mesh is closed, None otherwise. + + Notes + ----- + The volume is computed using the signed volume of tetrahedra formed by each + triangulated face and the origin. This method works for both convex and + non-convex meshes, as long as they are closed and properly oriented. + + The volume is only meaningful for closed meshes. For open meshes, this method + returns None. + + When faces are non-convex, the triangulation might not be correct, since it uses + the centroid of the face. For accurate results with non-convex faces, consider + using a mesh with triangulated faces. + + By default, the mesh is copied internally and face cycles are unified to ensure + correct orientation before computing the volume. These operations can be disabled + by setting ``copy=False`` and ``unify_cycles=False`` for performance in cases where + the mesh is already correctly oriented or when the original mesh can be modified. + + Examples + -------- + >>> from compas.datastructures import Mesh + >>> mesh = Mesh.from_polyhedron(6) # Create a cube + >>> volume = mesh.volume() + >>> volume is not None + True + + """ + if not self.is_closed(): + return None + + # Make a copy to avoid modifying the original mesh + mesh_to_use = self.copy() if copy else self + + # Unify cycles to ensure consistent face orientation + if unify_cycles: + mesh_to_use.unify_cycles() + + # Use built-in triangulation to get triangulated faces + vertices, faces = mesh_to_use.to_vertices_and_faces(triangulated=True) + + volume = 0.0 + for face in faces: + # Each face is now a triangle (3 vertices) + a, b, c = [vertices[i] for i in face] + # Signed volume of tetrahedron formed by triangle and origin + # V = (1/6) * (a dot (b cross c)) where a, b, c are the vertices + bc = cross_vectors(b, c) + vol = dot_vectors(a, bc) / 6.0 + volume += vol + + return abs(volume) + def centroid(self): """Calculate the mesh centroid. diff --git a/tests/compas/datastructures/test_mesh.py b/tests/compas/datastructures/test_mesh.py index f88d92bb30f8..319756f284b8 100644 --- a/tests/compas/datastructures/test_mesh.py +++ b/tests/compas/datastructures/test_mesh.py @@ -1069,6 +1069,54 @@ def test_normal(): ) +def test_volume(): + import math + + # Test with a box (3x4x5) + box = Box.from_width_height_depth(3, 4, 5) + mesh = Mesh.from_shape(box) + volume = mesh.volume() + expected_volume = 3 * 4 * 5 # 60 + assert TOL.is_close(volume, expected_volume) + + # Test with a smaller box (2x2x2) + box2 = Box.from_width_height_depth(2, 2, 2) + mesh2 = Mesh.from_shape(box2) + volume2 = mesh2.volume() + expected_volume2 = 2 * 2 * 2 # 8 + assert TOL.is_close(volume2, expected_volume2) + + # Test with a tetrahedron from polyhedron + # Regular tetrahedron with edge length ~1.633 has volume = edge^3 / (6*sqrt(2)) + tet = Mesh.from_polyhedron(4) + volume = tet.volume() + # Expected volume for the platonic tetrahedron from polyhedron(4) + expected_tet_volume = 0.5132002392796675 + assert TOL.is_close(volume, expected_tet_volume) + + # Test with a sphere approximation + sphere_mesh = Mesh.from_shape(Sphere(radius=1.0), u=32, v=32) + volume = sphere_mesh.volume() + assert volume is not None + expected_sphere_volume = (4.0 / 3.0) * math.pi + # Allow for ~1% error due to discretization + assert TOL.is_close(volume, expected_sphere_volume, rtol=0.02) + + # Test with an open mesh (should return None) + mesh = Mesh.from_obj(compas.get("faces.obj")) + volume = mesh.volume() + assert volume is None + + # Test optional parameters + box3 = Box.from_width_height_depth(2, 3, 4) + mesh3 = Mesh.from_shape(box3) + + # Test with copy=False and unify_cycles=False (should still work for well-oriented mesh) + volume3 = mesh3.volume(copy=False, unify_cycles=False) + expected_volume3 = 2 * 3 * 4 # 24 + assert TOL.is_close(volume3, expected_volume3) + + # -------------------------------------------------------------------------- # vertex geometry # --------------------------------------------------------------------------