diff --git a/tests/test_IO.py b/tests/test_IO.py index c75591c543..f82d78bfbc 100644 --- a/tests/test_IO.py +++ b/tests/test_IO.py @@ -23,7 +23,7 @@ def test_simulation_preserve_types(): st = GaussianPulse(freq0=1.0, fwidth=1.0) sim_all = Simulation( - size=(1.0, 1.0, 1.0), + size=(10.0, 10.0, 10.0), grid_size=(1, 1, 1), structures=[ Structure(geometry=Box(size=(1, 1, 1)), medium=Medium()), @@ -40,7 +40,13 @@ def test_simulation_preserve_types(): ], sources=[ VolumeSource(size=(0, 0, 0), source_time=st, polarization="Ex"), - PlaneWave(size=(inf, inf, 0), source_time=st, direction="+", polarization="Ex"), + PlaneWave( + size=(inf, inf, 0), + center=(0, 0, 4.5), + source_time=st, + direction="+", + polarization="Ex", + ), GaussianBeam( size=(inf, inf, 0), source_time=st, diff --git a/tests/test_components.py b/tests/test_components.py index 58f657c82b..ee0b8de152 100644 --- a/tests/test_components.py +++ b/tests/test_components.py @@ -7,6 +7,33 @@ from tidy3d.log import ValidationError, SetupError +def assert_log_level(caplog, log_level_expected): + """ensure something got logged if log_level is not None. + note: I put this here rather than utils.py because if we import from utils.py, + it will validate the sims there and those get included in log. + """ + + # get log output + logs = caplog.record_tuples + + # there's a log but the log level is not None (problem) + if logs and not log_level_expected: + raise Exception + + # we expect a log but none is given (problem) + if log_level_expected and not logs: + raise Exception + + # both expected and got log, check the log levels match + if logs and log_level_expected: + for log in logs: + log_level = log[1] + if log_level == log_level_expected: + # log level was triggered, exit + return + raise Exception + + def test_sim(): """make sure a simulation can be initialized""" @@ -116,7 +143,128 @@ def place_box(center_offset): def test_sim_grid_size(): size = (1, 1, 1) - s = Simulation(size=size, grid_size=(1.0, 1.0, 1.0)) + _ = Simulation(size=size, grid_size=(1.0, 1.0, 1.0)) + + +@pytest.mark.parametrize("fwidth,log_level", [(0.001, None), (3, 30)]) +def test_sim_frequency_range(caplog, fwidth, log_level): + # small fwidth should be inside range, large one should throw warning + + size = (1, 1, 1) + medium = Medium(frequency_range=(2, 3)) + box = Structure(geometry=Box(size=(0.1, 0.1, 0.1)), medium=medium) + src = VolumeSource( + source_time=GaussianPulse(freq0=2.4, fwidth=fwidth), + size=(0, 0, 0), + polarization="Ex", + ) + _ = Simulation(size=(1, 1, 1), grid_size=(0.1, 0.1, 0.1), structures=[box], sources=[src]) + + assert_log_level(caplog, log_level) + + +@pytest.mark.parametrize("grid_size,log_level", [(0.001, None), (3, 30)]) +def test_sim_grid_size(caplog, grid_size, log_level): + # small fwidth should be inside range, large one should throw warning + + medium = Medium(permittivity=2, frequency_range=(2e14, 3e14)) + box = Structure(geometry=Box(size=(0.1, 0.1, 0.1)), medium=medium) + src = VolumeSource( + source_time=GaussianPulse(freq0=2.5e14, fwidth=1e13), + size=(0, 0, 0), + polarization="Ex", + ) + _ = Simulation(size=(1, 1, 1), grid_size=(0.1, 0.1, grid_size), structures=[box], sources=[src]) + + assert_log_level(caplog, log_level) + + +@pytest.mark.parametrize("box_size,log_level", [(0.001, None), (9.9, 30), (20, None)]) +def test_sim_structure_gap(caplog, box_size, log_level): + """Make sure the gap between a structure and PML is not too small compared to lambda0.""" + medium = Medium(permittivity=2) + box = Structure(geometry=Box(size=(box_size, box_size, box_size)), medium=medium) + src = VolumeSource( + source_time=GaussianPulse(freq0=3e14, fwidth=1e13), + size=(0, 0, 0), + polarization="Ex", + ) + sim = Simulation( + size=(10, 10, 10), + grid_size=(0.1, 0.1, 0.1), + structures=[box], + sources=[src], + pml_layers=[PML(num_layers=5), PML(num_layers=5), PML(num_layers=5)], + ) + assert_log_level(caplog, log_level) + + +def test_sim_plane_wave_error(): + """ "Make sure we error if plane wave is not intersecting homogeneous region of simulation.""" + + medium_bg = Medium(permittivity=2) + medium_air = Medium(permittivity=1) + + box = Structure(geometry=Box(size=(0.1, 0.1, 0.1)), medium=medium_air) + + box_transparent = Structure(geometry=Box(size=(0.1, 0.1, 0.1)), medium=medium_bg) + + src = PlaneWave( + source_time=GaussianPulse(freq0=2.5e14, fwidth=1e13), + center=(0, 0, 0), + size=(inf, inf, 0), + direction="+", + polarization="Ex", + ) + + # with transparent box continue + _ = Simulation( + size=(1, 1, 1), + grid_size=(0.1, 0.1, 0.1), + medium=medium_bg, + structures=[box_transparent], + sources=[src], + ) + + # with non-transparent box, raise + with pytest.raises(SetupError): + _ = Simulation( + size=(1, 1, 1), + grid_size=(0.1, 0.1, 0.1), + medium=medium_bg, + structures=[box_transparent, box], + sources=[src], + ) + + +@pytest.mark.parametrize( + "box_size,log_level", + [((0.1, 0.1, 0.1), None), ((1, 0.1, 0.1), 30), ((0.1, 1, 0.1), 30), ((0.1, 0.1, 1), 30)], +) +def test_sim_structure_extent(caplog, box_size, log_level): + """Make sure we warn if structure extends exactly to simulation edges.""" + + box = Structure(geometry=Box(size=box_size), medium=Medium(permittivity=2)) + sim = Simulation(size=(1, 1, 1), grid_size=(0.1, 0.1, 0.1), structures=[box]) + + assert_log_level(caplog, log_level) + + +def test_num_mediums(): + """Make sure we error if too many mediums supplied.""" + + structures = [] + for i in range(200): + structures.append( + Structure(geometry=Box(size=(1, 1, 1)), medium=Medium(permittivity=i + 1)) + ) + sim = Simulation(size=(1, 1, 1), grid_size=(0.1, 0.1, 0.1), structures=structures) + + with pytest.raises(SetupError): + structures.append( + Structure(geometry=Box(size=(1, 1, 1)), medium=Medium(permittivity=i + 2)) + ) + sim = Simulation(size=(1, 1, 1), grid_size=(0.1, 0.1, 0.1), structures=structures) """ geometry """ @@ -259,7 +407,6 @@ def test_medium_dispersion_create(): def eps_compare(medium: Medium, expected: Dict, tol: float = 1e-5): for freq, val in expected.items(): - # print(f"Expected: {medium.eps_model(freq)}, got: {val}") assert np.abs(medium.eps_model(freq) - val) < tol @@ -340,7 +487,7 @@ def test_modes(): """ names """ -def test_names_default(): +def _test_names_default(): """makes sure default names are set""" sim = Simulation( diff --git a/tests/utils.py b/tests/utils.py index 8604534363..ad35a18469 100644 --- a/tests/utils.py +++ b/tests/utils.py @@ -34,7 +34,7 @@ def prepend_tmp(path): SIM_MONITORS = Simulation( - size=(2.0, 2.0, 2.0), + size=(10.0, 10.0, 10.0), grid_size=(0.1, 0.1, 0.1), monitors=[ FieldMonitor(size=(1, 1, 1), center=(0, 1, 0), freqs=[1, 2, 5, 7, 8], name="field_freq"), @@ -52,7 +52,7 @@ def prepend_tmp(path): ) SIM_FULL = Simulation( - size=(2.0, 2.0, 2.0), + size=(10.0, 10.0, 10.0), grid_size=(0.1, 0.1, 0.1), run_time=40e-11, structures=[ diff --git a/tidy3d/components/simulation.py b/tidy3d/components/simulation.py index 24db7c67cb..eaef4d2ea4 100644 --- a/tidy3d/components/simulation.py +++ b/tidy3d/components/simulation.py @@ -10,26 +10,30 @@ from mpl_toolkits.axes_grid1 import make_axes_locatable from descartes import PolygonPatch -from .validators import assert_unique_names, assert_objects_in_sim_bounds, set_names +from .validators import assert_unique_names, assert_objects_in_sim_bounds from .geometry import Box from .types import Symmetry, Ax, Shapely, FreqBound, GridSize from .grid import Coords1D, Grid, Coords from .medium import Medium, MediumType, AbstractMedium from .structure import Structure -from .source import SourceType +from .source import SourceType, PlaneWave from .monitor import MonitorType from .pml import PMLTypes, PML from .viz import StructMediumParams, StructEpsParams, PMLParams, SymParams, add_ax_if_none from ..constants import inf, C_0 -from ..log import log, Tidy3dKeyError +from ..log import log, Tidy3dKeyError, SetupError # for docstring examples from .geometry import Sphere, Cylinder, PolySlab # pylint:disable=unused-import from .source import VolumeSource, GaussianPulse # pylint:disable=unused-import from .monitor import FieldMonitor, FluxMonitor, Monitor # pylint:disable=unused-import -# technically this is creating a circular import issue because it calls tidy3d/__init__.py -# from .. import __version__ as version_number + +# minimum number of grid points allowed per central wavelength in a medium +MIN_GRIDS_PER_WVL = 6.0 + +# maximum number of mediums supported +MAX_NUM_MEDIUMS = 200 class Simulation(Box): # pylint:disable=too-many-public-methods @@ -170,9 +174,8 @@ def set_none_to_zero_layers(cls, val): _monitors_in_bounds = assert_objects_in_sim_bounds("monitors") # assign names to unnamed structures, sources, and mediums - _structure_names = set_names("structures") - _source_names = set_names("sources") - + # _structure_names = set_names("structures") + # _source_names = set_names("sources") # @pydantic.validator("structures", allow_reuse=True, always=True) # def set_medium_names(cls, val, values): # """check for intersection of each structure with simulation bounds.""" @@ -190,12 +193,209 @@ def set_none_to_zero_layers(cls, val): _unique_monitor_names = assert_unique_names("monitors") # _unique_medium_names = assert_unique_names("structures", check_mediums=True) - # TODO: - # - check sources in medium freq range - # - check PW in homogeneous medium - # - check nonuniform grid covers the whole simulation domain - # - check any structures close to PML (in lambda) without intersecting. - # - check any structures.bounds == simulation.bounds -> warn + # _few_enough_mediums = validate_num_mediums() + # _structures_not_at_edges = validate_structure_bounds_not_at_edges() + # _gap_size_ok = validate_pml_gap_size() + # _medium_freq_range_ok = validate_medium_frequency_range() + # _resolution_fine_enough = validate_resolution() + # _plane_waves_in_homo = validate_plane_wave_intersections() + + @pydantic.validator("structures", always=True) + def _validate_num_mediums(cls, val): + """Error if too many mediums present.""" + + if val is None: + return val + + mediums = {structure.medium for structure in val} + if len(mediums) > MAX_NUM_MEDIUMS: + raise SetupError( + f"Tidy3d only supports {MAX_NUM_MEDIUMS} distinct mediums." + f"{len(mediums)} were supplied." + ) + + return val + + @pydantic.validator("structures", always=True) + def _structures_not_at_edges(cls, val, values): + """Warn if any structures lie at the simulation boundaries.""" + + if val is None: + return val + + sim_box = Box(size=values.get("size"), center=values.get("center")) + sim_bound_min, sim_bound_max = sim_box.bounds + sim_bounds = list(sim_bound_min) + list(sim_bound_max) + + for structure in val: + struct_bound_min, struct_bound_max = structure.geometry.bounds + struct_bounds = list(struct_bound_min) + list(struct_bound_max) + + for sim_val, struct_val in zip(sim_bounds, struct_bounds): + + if np.isclose(sim_val, struct_val): + log.warning( + f"structure\n\n{structure}\n\nbounds extend exactly to simulation " + "edges. This can cause unexpected behavior. If intending to extend " + "the structure to infinity along one dimension, use td.inf as a " + "size variable instead to make this explicit." + ) + + return val + + @pydantic.validator("pml_layers", always=True) + def _structures_not_close_pml(cls, val, values): # pylint:disable=too-many-locals + """Warn if any structures lie at the simulation boundaries.""" + + if val is None: + return val + + sim_box = Box(size=values.get("size"), center=values.get("center")) + sim_bound_min, sim_bound_max = sim_box.bounds + + structures = values.get("structures") + sources = values.get("sources") + if (not structures) or (not sources): + return val + + def warn(structure): + """Warning message for a structure too close to PML.""" + log.warning( + f"a structure\n\n{structure}\n\nwas detected as being less " + "than half of a central wavelength from a PML on - side" + "To avoid inaccurate results, please increase gap between " + "any structures and PML or fully extend structure through the pml." + ) + + for structure in structures: + struct_bound_min, struct_bound_max = structure.geometry.bounds + + for source in sources: + fmin_src, fmax_src = source.source_time.frequency_range + f_average = (fmin_src + fmax_src) / 2.0 + lambda0 = C_0 / f_average + + for sim_val, struct_val, pml in zip(sim_bound_min, struct_bound_min, val): + if pml.num_layers > 0 and struct_val > sim_val: + if abs(sim_val - struct_val) < lambda0 / 2: + warn(structure) + + for sim_val, struct_val, pml in zip(sim_bound_max, struct_bound_max, val): + if pml.num_layers > 0 and struct_val < sim_val: + if abs(sim_val - struct_val) < lambda0 / 2: + warn(structure) + + return val + + @pydantic.validator("sources", always=True) + def _warn_sources_mediums_frequency_range(cls, val, values): + """Warn user if any sources have frequency range outside of medium frequency range.""" + + if val is None: + return val + + structures = values.get("structures") + structures = [] if not structures else structures + medium_bg = values.get("medium") + mediums = [medium_bg] + [structure.medium for structure in structures] + + for source in val: + fmin_src, fmax_src = source.source_time.frequency_range + for medium in mediums: + + # skip mediums that have no freq range (all freqs valid) + if medium.frequency_range is None: + continue + + # make sure medium frequency range includes the source frequency range + fmin_med, fmax_med = medium.frequency_range + if fmin_med > fmin_src or fmax_med < fmax_src: + log.warning( + f"A medium in the simulation:\n\n({medium})\n\nhas a frequency " + "range that does not fully cover the spetrum of a source:" + f"\n\n({source})\n\nThis can cause innacuracies in the " + "simulation results." + ) + return val + + @pydantic.validator("sources", always=True) + def _warn_grid_size_too_small(cls, val, values): # pylint:disable=too-many-locals + """Warn user if any grid size is too large compared to minimum wavelength in material.""" + + if val is None: + return val + + structures = values.get("structures") + structures = [] if not structures else structures + medium_bg = values.get("medium") + grid_size = values.get("grid_size") + mediums = [medium_bg] + [structure.medium for structure in structures] + + for source in val: + fmin_src, fmax_src = source.source_time.frequency_range + f_average = (fmin_src + fmax_src) / 2.0 + + for medium in mediums: + + eps_material = medium.eps_model(f_average) + n_material, _ = medium.eps_complex_to_nk(eps_material) + lambda_min = C_0 / f_average / n_material + + for dl in grid_size: + if isinstance(dl, float): + if dl > lambda_min / MIN_GRIDS_PER_WVL: + log.warning( + f"One of the grid sizes with a value of {dl:.4f} (um) " + "was detected as being too large when compared to the " + "central wavelength of a source within one of the " + f"simulation mediums {lambda_min:.4f} (um). " + "To avoid inaccuracies, it is reccomended the grid size is " + "reduced." + ) + # TODO: nonuniform grid + + return val + + @pydantic.validator("sources", always=True) + def _plane_wave_homogeneous(cls, val, values): + """Error if plane wave intersects""" + + if val is None: + return val + + # list of structures including background as a Box() + structure_bg = Structure( + geometry=Box( + size=values.get("size"), + center=values.get("center"), + ), + medium=values.get("medium"), + ) + + structures = values.get("structures") + structures = [] if not structures else structures + total_structures = [structure_bg] + structures + + # for each plane wave in the sources list + for source in val: + if isinstance(source, PlaneWave): + + # get all merged structures on the plane + normal_axis_index = source.size.index(0.0) + dim = "xyz"[normal_axis_index] + pos = source.center[normal_axis_index] + xyz_kwargs = {dim: pos} + structures_merged = cls._filter_structures_plane(total_structures, **xyz_kwargs) + + # make sure there is no more than one medium in the returned list + mediums = {medium for medium, _ in structures_merged} + if len(mediums) > 1: + raise SetupError( + f"{len(mediums)} different mediums detected on plane " + "intersecting a plane wave source. Plane must be homogeneous." + ) + + return val """ Accounting """ @@ -345,7 +545,7 @@ def plot_structures( """ medium_map = self.medium_map - medium_shapes = self._filter_plot_structures(x=x, y=y, z=z) + medium_shapes = self._filter_structures_plane(self.structures, x=x, y=y, z=z) for (medium, shape) in medium_shapes: params_updater = StructMediumParams(medium=medium, medium_map=medium_map) kwargs_struct = params_updater.update_params(**kwargs) @@ -412,7 +612,7 @@ def plot_structures_eps( # pylint: disable=too-many-arguments,too-many-locals eps_list.append(self.medium.eps_model(freq).real) eps_max = max(eps_list) eps_min = min(eps_list) - medium_shapes = self._filter_plot_structures(x=x, y=y, z=z) + medium_shapes = self._filter_structures_plane(self.structures, x=x, y=y, z=z) for (medium, shape) in medium_shapes: eps = medium.eps_model(freq).real params_updater = StructEpsParams(eps=eps, eps_max=eps_max) @@ -675,8 +875,9 @@ def _set_plot_bounds(self, ax: Ax, x: float = None, y: float = None, z: float = ax.set_ylim(ymin - pml_thick_y[0], ymax + pml_thick_y[1]) return ax - def _filter_plot_structures( - self, x: float = None, y: float = None, z: float = None + @staticmethod + def _filter_structures_plane( + structures: List[Structure], x: float = None, y: float = None, z: float = None ) -> List[Tuple[Medium, Shapely]]: """Compute list of shapes to plot on plane specified by {x,y,z}. Overlaps are removed or merged depending on medium. @@ -697,40 +898,18 @@ def _filter_plot_structures( """ shapes = [] - for struct in self.structures: + for structure in structures: # dont bother with geometries that dont intersect plane - if not struct.geometry.intersects_plane(x=x, y=y, z=z): + if not structure.geometry.intersects_plane(x=x, y=y, z=z): continue # get list of Shapely shapes that intersect at the plane - shapes_plane = struct.geometry.intersections(x=x, y=y, z=z) + shapes_plane = structure.geometry.intersections(x=x, y=y, z=z) # Append each of them and their medium information to the list of shapes for shape in shapes_plane: - shapes.append((struct.medium, shape)) - - # returned a merged list of mediums and shapes. - return self._merge_shapes(shapes) - - @staticmethod - def _merge_shapes(shapes: List[Tuple[Medium, Shapely]]) -> List[Tuple[Medium, Shapely]]: - # pylint:disable=line-too-long - """Merge list of shapes and mediums on plae by intersection with same medium. - Edit background shapes to remove volume by intersection. - - Parameters - ---------- - shapes : List[Tuple[:class:`Medium` or :class:`PoleResidue` or :class:`Lorentz` or :class:`Sellmeier` or :class:`Debye`, shapely.geometry.base.BaseGeometry]] - Ordered list of shapes and their mediums on a plane. - - Returns - ------- - List[Tuple[:AbstractMedium`, shapely.geometry.base.BaseGeometry]] - Shapes and their mediums on a plane - after merging and removing intersections with background. - """ - # pylint:enable=line-too-long + shapes.append((structure.medium, shape)) background_shapes = [] for medium, shape in shapes: