diff --git a/CHANGELOG.md b/CHANGELOG.md index 6af6d50b02..49a542e933 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -81,6 +81,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 - Fixed DRC parsing for quoted categories in klayout plugin. - Removed mode solver warnings about evaluating permittivity of a `Medium2D`. - Maximum number of grid points in an EME simulation is now based solely on transverse grid points. Maximum number of EME cells is unchanged. +- Fixed handling of polygons with holes (interiors) in `subdivide()` function. The function now properly converts polygons with interiors into `PolySlab` geometries using subtraction operations with `GeometryGroup`. ### Removed - Removed deprecated `use_complex_fields` parameter from `TwoPhotonAbsorption` and `KerrNonlinearity`. Parameters `beta` and `n2` are now real-valued only, as is `n0` if specified. diff --git a/tests/test_components/test_geometry.py b/tests/test_components/test_geometry.py index 33a396ff2b..d00a11515b 100644 --- a/tests/test_components/test_geometry.py +++ b/tests/test_components/test_geometry.py @@ -1158,6 +1158,71 @@ def test_subdivide(): subdivisions = subdivide(geom=overlapping_boxes, structures=[background_structure, box_sliver]) +def test_subdivide_geometry_group_with_polygon_holes(): + """Test that unionized geometry containing a hole works correctly.""" + mm = 1000.0 + + # Create four boxes arranged to form a cross pattern with a square hole in the middle + box_a = td.Box.from_bounds((-10 * mm, -10 * mm, 0 * mm), (-5 * mm, 5 * mm, 0 * mm)) + box_b = td.Box.from_bounds((-10 * mm, 5 * mm, 0 * mm), (5 * mm, 10 * mm, 0 * mm)) + box_c = td.Box.from_bounds((5 * mm, -5 * mm, 0 * mm), (10 * mm, 10 * mm, 0 * mm)) + box_d = td.Box.from_bounds((-5 * mm, -10 * mm, 0 * mm), (10 * mm, -5 * mm, 0 * mm)) + + geom_group = td.GeometryGroup(geometries=(box_a, box_b, box_c, box_d)) + geom_structures_group = [td.Structure(geometry=geom_group, medium=td.PEC2D)] + + feed_pin_bottom = -2 * mm + feed_pin_top = 0 * mm + feed_pin_center = 0.5 * (feed_pin_top + feed_pin_bottom) + feed_pin_length = feed_pin_top - feed_pin_bottom + feed_center_x = -7.5 * mm + feed_center_y = -7.5 * mm + rfeed = 1.0 * mm + + feed_pin = td.Structure( + geometry=td.Cylinder( + center=(feed_center_x, feed_center_y, feed_pin_center), + radius=rfeed, + length=feed_pin_length, + axis=2, + ), + medium=td.PECMedium(), + ) + + structures_list = [feed_pin, *geom_structures_group] + + freq = 1500 * 1e6 + dl = (td.C_0 / freq) / 300.0 + + mesh_overrides = [ + td.MeshOverrideStructure( + geometry=td.Box( + center=(0, 0, 0 * mm), + size=(20 * mm, 20 * mm, 6 * mm), + ), + dl=[dl, dl, dl], + ) + ] + + sim = td.Simulation( + size=[100 * mm, 100 * mm, 30 * mm], + grid_spec=td.GridSpec.auto( + min_steps_per_wvl=20, + wavelength=td.C_0 / freq, + override_structures=mesh_overrides, + ), + structures=structures_list, + run_time=1e-13, + ) + + contains_difference_operation = False + for structure in sim._finalized.structures: + geo = structure.geometry + if isinstance(geo, td.ClipOperation) and geo.operation == "difference": + contains_difference_operation = True + assert contains_difference_operation + + @pytest.mark.parametrize("snap_location", [SnapLocation.Boundary, SnapLocation.Center]) @pytest.mark.parametrize( "snap_behavior", diff --git a/tidy3d/components/geometry/base.py b/tidy3d/components/geometry/base.py index cb394ee2fa..441a6bb7ee 100644 --- a/tidy3d/components/geometry/base.py +++ b/tidy3d/components/geometry/base.py @@ -5,7 +5,7 @@ import functools import pathlib from abc import ABC, abstractmethod -from collections.abc import Iterable +from collections.abc import Iterable, Sequence from os import PathLike from typing import TYPE_CHECKING, Any, Callable, Optional, Union @@ -721,11 +721,14 @@ def evaluate_inf_shape(shape: Shapely) -> Shapely: if not any(np.isinf(b) for b in shape.bounds): return shape + def _processed_coords(coords: Sequence[tuple[Any, ...]]) -> list[tuple[float, ...]]: + evaluated = Geometry._evaluate_inf(np.array(coords)) + return [tuple(point) for point in evaluated.tolist()] + if shape.geom_type == "Polygon": - return shapely.Polygon( - Geometry._evaluate_inf(np.array(shape.exterior.coords)), - [Geometry._evaluate_inf(np.array(g.coords)) for g in shape.interiors], - ) + shell = _processed_coords(shape.exterior.coords) + holes = [_processed_coords(g.coords) for g in shape.interiors] + return shapely.Polygon(shell, holes) if shape.geom_type in {"Point", "LineString", "LinearRing"}: return shape.__class__(Geometry._evaluate_inf(np.array(shape.coords))) if shape.geom_type in { diff --git a/tidy3d/components/geometry/utils_2d.py b/tidy3d/components/geometry/utils_2d.py index 08e4785fdc..28c12c1380 100644 --- a/tidy3d/components/geometry/utils_2d.py +++ b/tidy3d/components/geometry/utils_2d.py @@ -7,7 +7,7 @@ import numpy as np import shapely -from tidy3d.components.geometry.base import Box, ClipOperation, Geometry +from tidy3d.components.geometry.base import Box, ClipOperation, Geometry, GeometryGroup from tidy3d.components.geometry.polyslab import _MIN_POLYGON_AREA, PolySlab from tidy3d.components.grid.grid import Grid from tidy3d.components.scene import Scene @@ -122,10 +122,28 @@ def subdivide( """ - def shapely_to_polyslab(polygon: shapely.Polygon, axis: Axis, center: float) -> PolySlab: - xx, yy = polygon.exterior.coords.xy - vertices = list(zip(xx, yy)) - return PolySlab(slab_bounds=(center, center), vertices=vertices, axis=axis) + def shapely_to_polyslab(polygon: shapely.Polygon, axis: Axis, center: float) -> Geometry: + def ring_vertices(ring: shapely.LinearRing) -> list[tuple[float, float]]: + xx, yy = ring.coords.xy + return list(zip(xx, yy)) + + polyslab = PolySlab( + slab_bounds=(center, center), + vertices=ring_vertices(polygon.exterior), + axis=axis, + ) + if len(polygon.interiors) == 0: + return polyslab + + interiors = [ + PolySlab( + slab_bounds=(center, center), + vertices=ring_vertices(interior), + axis=axis, + ) + for interior in polygon.interiors + ] + return polyslab - GeometryGroup(geometries=interiors) def to_multipolygon(shapely_geometry: Shapely) -> shapely.MultiPolygon: return shapely.MultiPolygon(ClipOperation.to_polygon_list(shapely_geometry))