✨ FEAT: Add reflection transformation with verification tests

CHANGELOG.md

@@ -7,6 +7,9 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 ## [Unreleased]
 
+### Added
+- The method `Geometry.reflected` can be used to create a reflected copy of any geometry off a plane. As for other transformations, for efficiency, `reflected` `PolySlab` directly returns an updated `PolySlab` object rather than a `Transformed` object, except when the normal of the plane of reflection has a non-zero component along the slab axis, in which case `Transformed` is still returned.
+
 ## [2.8.3] - 2025-04-24
 
 ### Added

tests/test_components/test_geometry.py

@@ -439,6 +439,7 @@ def test_planar_transform(axis):
     geo = (
         td.Box(size=(3 * axis, 2 * abs(axis - 1), 4 * (2 - axis)))
         .rotated(2.0, axis)
+        .reflected((axis, 2 * (axis - 1), 3 * (axis - 2)))
         .translated(-1, 2, 3)
         .scaled(1.4, -1.2, 1.3)
     )
@@ -479,6 +480,30 @@ def test_transforms():
     assert len(geo.intersections_plane(x=0)) == 1
     assert len(geo.intersections_plane(z=0)) == 3
 
+    # Test reflection of a Box across the XY plane and verify point inclusion.
+    xyz = (np.array([1, 1, 1, 3]), np.array([1, 1, 1, 3]), np.array([1, -1, -1.5, 3]))
+    geo = td.Box(center=(1, 1, 1), size=(2, 2, 2))
+    assert (geo.inside(*xyz) == (True, False, False, False)).all()
+    geo = geo.reflected((0, 0, 1))
+    assert (geo.inside(*xyz) == (False, True, True, False)).all()
+
+    # Test Sphere multiple reflections not influencing point inclusion.
+    xyz = (np.array([1, 2, 2, 1]), np.array([2, 1, 2, 3]), np.array([2, -2, -0.5, -2]))
+    geo = td.Sphere(radius=3.5)
+    assert (geo.inside(*xyz) == (True, True, True, False)).all()
+    geo = geo.reflected((2, 3, 1)).reflected((1, 2, 3))
+    assert (geo.inside(*xyz) == (True, True, True, False)).all()
+
+    # Test PolySlab reflection across non-axis plane and verify point inclusion.
+    xyz = (np.array([0, 1.5, -1.5, -1.5]), np.array([0, 1.5, -1.5, -2.5]), np.array([0, 0, 0, 0]))
+    geo = td.PolySlab(
+        vertices=[(1, 0), (3, 2), (2, 2), (0, 0), (2, -2), (3, -2)],
+        slab_bounds=(-1, 1),
+    )
+    assert (geo.inside(*xyz) == (True, True, False, False)).all()
+    geo = geo.reflected((1, 1, 0))
+    assert (geo.inside(*xyz) == (True, False, True, True)).all()
+
 
 def test_polyslab_transforms():
     # More tests on PolySlab tranforms matching direct Transformed
@@ -501,6 +526,10 @@ def test_polyslab_transforms():
     geo = geo.rotated(0.3, (0, -0.2, 0))
     assert geo.type != geo_trans.type
     assert np.allclose(geo.inside(*xyz), geo_trans.inside(*xyz))
+    geo_trans = td.Transformed(geometry=geo, transform=td.Transformed.reflection((1, 0, 2)))
+    geo = geo.reflected((1, 0, 2))
+    assert geo.type != geo_trans.type
+    assert np.allclose(geo.inside(*xyz), geo_trans.inside(*xyz))
 
 
 def test_general_rotation():
@@ -514,6 +543,21 @@ def test_general_rotation():
     assert np.allclose(td.Transformed.rotation(0.3, 2), td.Transformed.rotation(-0.3, [0, 0, -4]))
 
 
+def test_general_reflection():
+    # Magnitude of normal direction does not affect the transformation.
+    assert np.allclose(td.Transformed.reflection((1, 1, 0)), td.Transformed.reflection((5, 5, 0)))
+    assert np.allclose(td.Transformed.reflection((1, 0, 1)), td.Transformed.reflection((5, 0, 5)))
+    assert np.allclose(td.Transformed.reflection((0, 1, 1)), td.Transformed.reflection((0, 5, 5)))
+    # Negative normal direction means the same transformation.
+    assert np.allclose(td.Transformed.reflection((1, 1, 0)), td.Transformed.reflection((-1, -1, 0)))
+    assert np.allclose(td.Transformed.reflection((1, 0, 1)), td.Transformed.reflection((-1, 0, -1)))
+    assert np.allclose(td.Transformed.reflection((0, 1, 1)), td.Transformed.reflection((0, -1, -1)))
+    # Magnitude and sign of normal direction does not affect the transformation.
+    assert np.allclose(td.Transformed.reflection((1, 1, 0)), td.Transformed.reflection((-5, -5, 0)))
+    assert np.allclose(td.Transformed.reflection((1, 0, 1)), td.Transformed.reflection((-5, 0, -5)))
+    assert np.allclose(td.Transformed.reflection((0, 1, 1)), td.Transformed.reflection((0, -5, -5)))
+
+
 def test_flattening():
     flat = list(
         flatten_groups(

tests/test_components/test_structure.py

@@ -138,6 +138,31 @@ def test_invalid_polyslab(axis):
     geo7 = td.GeometryGroup(geometries=[(ps - box).rotated(np.pi / 4, j)]).rotated(-np.pi / 4, j)
     _ = td.Structure(geometry=geo7, medium=medium)
 
+    # normal with zero component along the slab axis
+    n1 = (axis, 2 * (axis - 1), 3 * (axis - 2))
+    # normal with non-zero component along the slab axis
+    n2 = (2 * (axis - 1), 3 * (axis - 2), axis)
+
+    geo8 = ps.reflected(n1)
+    _ = td.Structure(geometry=geo8, medium=medium)
+
+    geo9 = ps.reflected(n2).scaled(2, 2, 2).translated(-1, 0.5, 2).reflected(n2)
+    _ = td.Structure(geometry=geo9, medium=medium)
+
+    geo10 = ps.reflected(n1)
+    _ = td.Structure(geometry=geo10, medium=medium)
+
+    geo11 = ps.reflected(n2)
+    with pytest.raises(pd.ValidationError):
+        _ = td.Structure(geometry=geo11, medium=medium)
+
+    geo12 = td.GeometryGroup(geometries=[ps]).reflected(n2)
+    with pytest.raises(pd.ValidationError):
+        _ = td.Structure(geometry=geo12, medium=medium)
+
+    geo13 = td.GeometryGroup(geometries=[(ps - box).reflected(n2)]).reflected(n2)
+    _ = td.Structure(geometry=geo13, medium=medium)
+
 
 def test_validation_of_structures_with_2d_materials():
     med2d = td.Medium2D(ss=td.PEC, tt=td.PEC)

tidy3d/components/geometry/base.py

@@ -30,7 +30,7 @@
 from ..autograd import AutogradFieldMap, TracedCoordinate, TracedSize, get_static
 from ..autograd.derivative_utils import DerivativeInfo, integrate_within_bounds
 from ..base import Tidy3dBaseModel, cached_property
-from ..transformation import RotationAroundAxis
+from ..transformation import ReflectionFromPlane, RotationAroundAxis
 from ..types import (
     ArrayFloat2D,
     ArrayFloat3D,
@@ -975,6 +975,22 @@ def rotated(self, angle: float, axis: Union[Axis, Coordinate]) -> Geometry:
         """
         return Transformed(geometry=self, transform=Transformed.rotation(angle, axis))
 
+    def reflected(self, normal: Coordinate) -> Geometry:
+        """Return a reflected copy of this geometry.
+
+        Parameters
+        ----------
+        normal : Tuple[float, float, float]
+            The 3D normal vector of the plane of reflection. The plane is assumed
+                to pass through the origin (0,0,0).
+
+        Returns
+        -------
+        :class:`Geometry`
+            Reflected copy of this geometry.
+        """
+        return Transformed(geometry=self, transform=Transformed.reflection(normal))
+
     """ Field and coordinate transformations """
 
     @staticmethod
@@ -2877,6 +2893,25 @@ def rotation(angle: float, axis: Union[Axis, Coordinate]) -> MatrixReal4x4:
         transform[:3, :3] = RotationAroundAxis(angle=angle, axis=axis).matrix
         return transform
 
+    @staticmethod
+    def reflection(normal: Coordinate) -> MatrixReal4x4:
+        """Return a reflection matrix.
+
+        Parameters
+        ----------
+        normal : Tuple[float, float, float]
+            Normal of the plane of reflection.
+
+        Returns
+        -------
+        numpy.ndarray
+            Transform matrix with shape (4, 4).
+        """
+
+        transform = np.eye(4)
+        transform[:3, :3] = ReflectionFromPlane(normal=normal).matrix
+        return transform
+
     @staticmethod
     def preserves_axis(transform: MatrixReal4x4, axis: Axis) -> bool:
         """Indicate if the transform preserves the orientation of a given axis.

tidy3d/components/geometry/polyslab.py

@@ -19,7 +19,7 @@
 from ..autograd.derivative_utils import DerivativeInfo, DerivativeSurfaceMesh
 from ..autograd.types import TracedFloat
 from ..base import cached_property, skip_if_fields_missing
-from ..transformation import RotationAroundAxis
+from ..transformation import ReflectionFromPlane, RotationAroundAxis
 from ..types import (
     ArrayFloat1D,
     ArrayFloat2D,
@@ -1726,6 +1726,33 @@ def rotated(self, angle: float, axis: Union[Axis, Coordinate]) -> PolySlab:
 
         return super().rotated(angle=angle, axis=axis)
 
+    def reflected(self, normal: Coordinate) -> PolySlab:
+        """Return a reflected copy of this geometry.
+
+            Parameters
+            ----------
+            normal : Tuple[float, float, float]
+                The 3D normal vector of the plane of reflection. The plane is assumed
+                    to pass through the origin (0,0,0).
+
+            Returns
+            -------
+        -------
+        :class:`PolySlab`
+            Reflected copy of this ``PolySlab``.
+        """
+        if math.isclose(normal[self.axis], 0):
+            _, plane_axs = self.pop_axis((0, 1, 2), self.axis)
+            verts_3d = np.zeros((3, self.vertices.shape[0]))
+            verts_3d[plane_axs[0], :] = self.vertices[:, 0]
+            verts_3d[plane_axs[1], :] = self.vertices[:, 1]
+            reflection = ReflectionFromPlane(normal=normal)
+            reflected_vertices = reflection.reflect_vector(verts_3d)
+            reflected_vertices = reflected_vertices[plane_axs, :].T
+            return self.updated_copy(vertices=reflected_vertices)
+
+        return super().reflected(normal=normal)
+
 
 class ComplexPolySlabBase(PolySlab):
     """Interface for dividing a complex polyslab where self-intersecting polygon can

tidy3d/components/transformation.py

@@ -127,4 +127,78 @@ def matrix(self) -> TensorReal:
         return R
 
 
+class AbstractReflection(ABC, Tidy3dBaseModel):
+    """Abstract reflection of vectors and tensors."""
+
+    @cached_property
+    @abstractmethod
+    def matrix(self) -> TensorReal:
+        """Reflection matrix."""
+
+    def reflect_vector(self, vector: ArrayFloat2D) -> ArrayFloat2D:
+        """Reflect a vector/point or a list of vectors/points.
+
+        Parameters
+        ----------
+        vector : ArrayLike[float]
+            Array of shape ``(3, ...)``.
+
+        Returns
+        -------
+        Coordinate
+            Reflected vector.
+        """
+
+        if len(vector.shape) == 1:
+            return self.matrix @ vector
+
+        return np.tensordot(self.matrix, vector, axes=1)
+
+    def reflect_tensor(self, tensor: TensorReal) -> TensorReal:
+        """Reflect a tensor.
+
+        Parameters
+        ----------
+        tensor : ArrayLike[float]
+            Array of shape ``(3, 3)``.
+
+        Returns
+        -------
+        TensorReal
+            Reflected tensor.
+        """
+
+        return np.matmul(self.matrix, np.matmul(tensor, self.matrix.T))
+
+
+class ReflectionFromPlane(AbstractReflection):
+    """Reflection of vectors and tensors around a given vector."""
+
+    normal: Coordinate = pd.Field(
+        (1, 0, 0),
+        title="Normal of the reflecting plane",
+        description="A vector that specifies the normal of the plane of reflection",
+    )
+
+    @pd.validator("normal")
+    def _guarantee_nonzero_normal(cls, val):
+        norm = np.linalg.norm(val)
+        if np.isclose(norm, 0):
+            raise ValidationError(
+                "The norm of vector 'normal' cannot be zero. Please provide a proper normal vector."
+            )
+        return val
+
+    @cached_property
+    def matrix(self) -> TensorReal:
+        """Reflection matrix."""
+
+        norm = np.linalg.norm(self.normal)
+        n = self.normal / norm
+        R = np.eye(3) - 2 * np.outer(n, n)
+
+        return R
+
+
 RotationType = Union[RotationAroundAxis]
+ReflectionType = Union[ReflectionFromPlane]