From e0af8d26fb520ce6a3525a681102171b7febe3e5 Mon Sep 17 00:00:00 2001 From: Weiliang Jin Date: Mon, 8 Apr 2024 16:46:30 -0700 Subject: [PATCH] Warning if a nonuniform custom medium is intersecting certain sources and monitors Add is_spatially_uniform property to materials Add is_uniform property to data array and unstructured dataset --- CHANGELOG.md | 1 + tests/test_components/test_custom.py | 108 ++++++++++++++++++++++++++- tests/test_data/test_data_arrays.py | 22 ++++++ tests/test_data/test_datasets.py | 31 +++++++- tidy3d/components/data/data_array.py | 6 ++ tidy3d/components/data/dataset.py | 5 ++ tidy3d/components/medium.py | 76 +++++++++++++++++++ tidy3d/components/simulation.py | 81 +++++++++++++------- 8 files changed, 300 insertions(+), 30 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index f9e07c633..6ea8d10c5 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added - Support for `.gz` files in `Simulation` version updater. +- Warning if a nonuniform custom medium is intersecting `PlaneWave`, `GaussianBeam`, `AstigmaticGaussianBeam`, `FieldProjectionCartesianMonitor`, `FieldProjectionAngleMonitor`, `FieldProjectionKSpaceMonitor`, and `DiffractionMonitor`. ### Changed diff --git a/tests/test_components/test_custom.py b/tests/test_components/test_custom.py index 9a5cbf54f..9f8dcdd99 100644 --- a/tests/test_components/test_custom.py +++ b/tests/test_components/test_custom.py @@ -64,12 +64,16 @@ def make_custom_current_source(): return td.CustomCurrentSource(size=SIZE, source_time=ST, current_dataset=current_dataset) -def make_spatial_data(value=0, dx=0, unstructured=False, seed=None): +def make_spatial_data(value=0, dx=0, unstructured=False, seed=None, uniform=False): """Makes a spatial data array.""" - data = np.random.random((Nx, Ny, Nz)) + value + if uniform: + data = value * np.ones((Nx, Ny, Nz)) + else: + data = np.random.random((Nx, Ny, Nz)) + value arr = td.SpatialDataArray(data, coords=dict(x=X + dx, y=Y, z=Z)) if unstructured: - return cartesian_to_unstructured(arr, seed=seed) + method = "direct" if uniform else "linear" + return cartesian_to_unstructured(arr, seed=seed, method=method) return arr @@ -634,12 +638,21 @@ def test_custom_isotropic_medium(unstructured): mat = CustomMedium(permittivity=permittivity, conductivity=sigmatmp) mat = CustomMedium(permittivity=permittivity, conductivity=sigmatmp, allow_gain=True) verify_custom_medium_methods(mat, ["permittivity", "conductivity"]) + assert not mat.is_spatially_uniform # inconsistent coords with pytest.raises(pydantic.ValidationError): sigmatmp = make_spatial_data(value=0, dx=1, unstructured=unstructured, seed=seed) mat = CustomMedium(permittivity=permittivity, conductivity=sigmatmp) + # uniform + permittivity = make_spatial_data(value=1, unstructured=unstructured, seed=seed, uniform=True) + mat = CustomMedium(permittivity=permittivity) + assert mat.is_spatially_uniform + + mat = CustomAnisotropicMedium(xx=mat, yy=mat, zz=mat) + assert mat.is_spatially_uniform + def verify_custom_dispersive_medium_methods(mat, reduced_fields=[]): """Verify that the methods in custom dispersive medium is producing expected results.""" @@ -709,6 +722,7 @@ def test_custom_pole_residue(unstructured): mat = CustomPoleResidue(eps_inf=eps_inf, poles=((a, c),)) verify_custom_dispersive_medium_methods(mat, ["eps_inf", "poles"]) assert mat.n_cfl > 1 + assert not mat.is_spatially_uniform # to custom non-dispersive medium # dispersive failure @@ -775,6 +789,7 @@ def test_custom_sellmeier(unstructured): mat = CustomSellmeier(coeffs=((b1, c1), (b2, c2))) verify_custom_dispersive_medium_methods(mat, ["coeffs"]) assert mat.n_cfl == 1 + assert not mat.is_spatially_uniform # from dispersion n = make_spatial_data(value=2, unstructured=unstructured, seed=seed) @@ -839,6 +854,7 @@ def test_custom_lorentz(unstructured): verify_custom_dispersive_medium_methods(mat, ["eps_inf", "coeffs"]) assert mat.n_cfl > 1 assert mat.pole_residue.subpixel + assert not mat.is_spatially_uniform @pytest.mark.parametrize("unstructured", [False, True]) @@ -876,6 +892,7 @@ def test_custom_drude(unstructured): mat = CustomDrude(eps_inf=eps_inf, coeffs=((f1, delta1), (f2, delta2))) verify_custom_dispersive_medium_methods(mat, ["eps_inf", "coeffs"]) assert mat.n_cfl > 1 + assert not mat.is_spatially_uniform @pytest.mark.parametrize("unstructured", [False, True]) @@ -926,6 +943,7 @@ def test_custom_debye(unstructured): mat = CustomDebye(eps_inf=eps_inf, coeffs=((eps1, tau1), (eps2, tau2))) verify_custom_dispersive_medium_methods(mat, ["eps_inf", "coeffs"]) assert mat.n_cfl > 1 + assert not mat.is_spatially_uniform @pytest.mark.parametrize("unstructured", [True]) @@ -957,6 +975,7 @@ def test_custom_anisotropic_medium(log_capture, unstructured): # anisotropic mat = CustomAnisotropicMedium(xx=mat_xx, yy=mat_yy, zz=mat_zz) verify_custom_medium_methods(mat) + assert not mat.is_spatially_uniform mat = CustomAnisotropicMedium(xx=mat_xx, yy=mat_yy, zz=mat_zz, subpixel=True) assert_log_level(log_capture, "WARNING") @@ -1075,3 +1094,86 @@ def test_io_dispersive(tmp_path, unstructured, z_custom): sim_load = td.Simulation.from_file(filename) assert sim_load == sim + + +def test_warn_planewave_intersection(log_capture): + """Warn that if a nonuniform custom medium is intersecting PlaneWave source.""" + src = td.PlaneWave( + source_time=td.GaussianPulse(freq0=3e14, fwidth=1e13), + center=(0, 0, 0), + size=(td.inf, td.inf, 0), + direction="+", + ) + + # uniform custom medium + permittivity = make_spatial_data(value=1, unstructured=False, seed=0, uniform=True) + mat = CustomMedium(permittivity=permittivity) + box = td.Structure( + geometry=td.Box(size=(td.inf, td.inf, 1)), + medium=mat, + ) + + sim = td.Simulation( + size=(1, 1, 2), + structures=[box], + grid_spec=td.GridSpec.auto(wavelength=1), + sources=[src], + run_time=1e-12, + boundary_spec=td.BoundarySpec.all_sides(boundary=td.Periodic()), + ) + assert_log_level(log_capture, None) + + # nonuniform custom medium + permittivity = make_spatial_data(value=1, unstructured=False, seed=0, uniform=False) + mat = CustomMedium(permittivity=permittivity) + box = td.Structure( + geometry=td.Box(size=(td.inf, td.inf, 1)), + medium=mat, + ) + sim.updated_copy(structures=[box]) + assert_log_level(log_capture, "WARNING") + + +def test_warn_diffraction_monitor_intersection(log_capture): + """Warn that if a nonuniform custom medium is intersecting Diffraction Monitor.""" + src = td.PointDipole( + source_time=td.GaussianPulse(freq0=2.5e14, fwidth=1e13), + center=(0, 0, 0.6), + polarization="Ex", + ) + monitor = td.DiffractionMonitor( + center=(0, 0, 0), + size=(td.inf, td.inf, 0), + freqs=[250e12], + name="monitor_diffraction", + normal_dir="+", + ) + + # uniform custom medium + permittivity = make_spatial_data(value=1, unstructured=False, seed=0, uniform=True) + mat = CustomMedium(permittivity=permittivity) + box = td.Structure( + geometry=td.Box(size=(td.inf, td.inf, 1)), + medium=mat, + ) + + sim = td.Simulation( + size=(1, 1, 2), + structures=[box], + grid_spec=td.GridSpec.auto(wavelength=1), + monitors=[monitor], + sources=[src], + run_time=1e-12, + boundary_spec=td.BoundarySpec.all_sides(boundary=td.Periodic()), + ) + assert_log_level(log_capture, None) + + # nonuniform custom medium + permittivity = make_spatial_data(value=1, unstructured=False, seed=0, uniform=False) + mat = CustomMedium(permittivity=permittivity) + box = td.Structure( + geometry=td.Box(size=(td.inf, td.inf, 1)), + medium=mat, + ) + sim.updated_copy(structures=[box]) + assert_log_level(log_capture, "WARNING") diff --git a/tests/test_data/test_data_arrays.py b/tests/test_data/test_data_arrays.py index 83c1378f0..8b0b3eb4c 100644 --- a/tests/test_data/test_data_arrays.py +++ b/tests/test_data/test_data_arrays.py @@ -381,3 +381,25 @@ def test_sel_inside(nx): with pytest.raises(DataError): _ = arr.does_cover([[0.1, 3, 2], [1, 2.5, 2]]) + + +def test_uniform_check(): + """check if each element in the array is of equal value.""" + arr = td.SpatialDataArray( + np.ones((2, 2, 2), dtype=np.complex128), + coords=dict(x=[0, 1], y=[1, 2], z=[2, 3]), + ) + assert arr.is_uniform + + # small variation is still considered as uniform + arr = td.SpatialDataArray( + np.ones((2, 2, 2)) + np.random.random((2, 2, 2)) * 1e-6, + coords=dict(x=[0, 1], y=[1, 2], z=[2, 3]), + ) + assert arr.is_uniform + + arr = td.SpatialDataArray( + np.ones((2, 2, 2)) + np.random.random((2, 2, 2)) * 1e-4, + coords=dict(x=[0, 1], y=[1, 2], z=[2, 3]), + ) + assert not arr.is_uniform diff --git a/tests/test_data/test_datasets.py b/tests/test_data/test_datasets.py index 9d160d5af..cfaa50aa7 100644 --- a/tests/test_data/test_datasets.py +++ b/tests/test_data/test_datasets.py @@ -39,7 +39,7 @@ def test_triangular_dataset(log_capture, tmp_path, ds_name, no_vtk=False): cells=tri_grid_cells, values=tri_grid_values, ) - + assert not tri_grid.is_uniform # test name redirect assert tri_grid.name == ds_name @@ -594,3 +594,32 @@ def test_cartesian_to_unstructured(nz, use_vtk, fill_value): assert sample_outside.values.item() == values[0, 0, 0] else: assert sample_outside.values.item() == fill_value + + +def test_triangular_dataset_uniform(): + import tidy3d as td + + # basic create + tri_grid_points = td.PointDataArray( + [[0.0, 0.0], [1.0, 0.0], [0.0, 1.0], [1.0, 1.0]], + dims=("index", "axis"), + ) + + tri_grid_cells = td.CellDataArray( + [[0, 1, 2], [1, 2, 3]], + dims=("cell_index", "vertex_index"), + ) + + tri_grid_values = td.IndexedDataArray( + [1.0, 1.0, 1.0, 1.0], + dims=("index"), + ) + + tri_grid = td.TriangularGridDataset( + normal_axis=1, + normal_pos=0, + points=tri_grid_points, + cells=tri_grid_cells, + values=tri_grid_values, + ) + assert tri_grid.is_uniform diff --git a/tidy3d/components/data/data_array.py b/tidy3d/components/data/data_array.py index f559c3eb8..86c304c05 100644 --- a/tidy3d/components/data/data_array.py +++ b/tidy3d/components/data/data_array.py @@ -144,6 +144,12 @@ def abs(self): """Absolute value of data array.""" return abs(self) + @property + def is_uniform(self): + """Whether each element is of equal value in the data array""" + raw_data = self.data.ravel() + return np.allclose(raw_data, raw_data[0]) + def to_hdf5(self, fname: Union[str, h5py.File], group_path: str) -> None: """Save an xr.DataArray to the hdf5 file or file handle with a given path to the group.""" diff --git a/tidy3d/components/data/dataset.py b/tidy3d/components/data/dataset.py index 94524dd7f..651cd1447 100644 --- a/tidy3d/components/data/dataset.py +++ b/tidy3d/components/data/dataset.py @@ -563,6 +563,11 @@ def points_right_dims(cls, val): ) return val + @property + def is_uniform(self): + """Whether each element is of equal value in ``values``.""" + return self.values.is_uniform + @pd.validator("points", always=True) def points_right_indexing(cls, val): """Check that points are indexed corrrectly.""" diff --git a/tidy3d/components/medium.py b/tidy3d/components/medium.py index 6194346e9..18a7e7167 100644 --- a/tidy3d/components/medium.py +++ b/tidy3d/components/medium.py @@ -680,6 +680,11 @@ def _validate_modulation_spec(cls, val, values): _name_validator = validate_name_str() + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + return True + @cached_property def is_time_modulated(self) -> bool: """Whether any component of the medium is time modulated.""" @@ -1484,6 +1489,13 @@ def _passivity_validation(cls, val, values): ) return val + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + if self.conductivity is None: + return self.permittivity.is_uniform + return self.permittivity.is_uniform and self.conductivity.is_uniform + @cached_property def n_cfl(self): """This property computes the index of refraction related to CFL condition, so that @@ -1841,6 +1853,11 @@ def is_isotropic(self) -> bool: return True return False + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + return self._medium.is_spatially_uniform + @cached_property def freqs(self) -> np.ndarray: """float array of frequencies. @@ -2920,6 +2937,18 @@ def _poles_correct_shape(cls, val, values): ) return val + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + if not self.eps_inf.is_uniform: + return False + + for coeffs in self.poles: + for coeff in coeffs: + if not coeff.is_uniform: + return False + return True + def eps_dataarray_freq( self, frequency: float ) -> Tuple[CustomSpatialDataType, CustomSpatialDataType, CustomSpatialDataType]: @@ -3273,6 +3302,15 @@ def _passivity_validation(cls, val, values): ) return val + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + for coeffs in self.coeffs: + for coeff in coeffs: + if not coeff.is_uniform: + return False + return True + def _pole_residue_dict(self) -> Dict: """Dict representation of Medium as a pole-residue model.""" poles_dict = Sellmeier._pole_residue_dict(self) @@ -3698,6 +3736,17 @@ def _passivity_validation(cls, val, values): ) return val + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + if not self.eps_inf.is_uniform: + return False + for coeffs in self.coeffs: + for coeff in coeffs: + if not coeff.is_uniform: + return False + return True + def eps_dataarray_freq( self, frequency: float ) -> Tuple[CustomSpatialDataType, CustomSpatialDataType, CustomSpatialDataType]: @@ -3946,6 +3995,17 @@ def _coeffs_correct_shape_and_sign(cls, val, values): raise SetupError("For stable medium, 'delta' must be positive.") return val + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + if not self.eps_inf.is_uniform: + return False + for coeffs in self.coeffs: + for coeff in coeffs: + if not coeff.is_uniform: + return False + return True + def eps_dataarray_freq( self, frequency: float ) -> Tuple[CustomSpatialDataType, CustomSpatialDataType, CustomSpatialDataType]: @@ -4207,6 +4267,17 @@ def _passivity_validation(cls, val, values): ) return val + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + if not self.eps_inf.is_uniform: + return False + for coeffs in self.coeffs: + for coeff in coeffs: + if not coeff.is_uniform: + return False + return True + def eps_dataarray_freq( self, frequency: float ) -> Tuple[CustomSpatialDataType, CustomSpatialDataType, CustomSpatialDataType]: @@ -4845,6 +4916,11 @@ def _ignored_fields(cls, values): ) return values + @cached_property + def is_spatially_uniform(self) -> bool: + """Whether the medium is spatially uniform.""" + return any(comp.is_spatially_uniform for comp in self.components.values()) + @cached_property def n_cfl(self): """This property computes the index of refraction related to CFL condition, so that diff --git a/tidy3d/components/simulation.py b/tidy3d/components/simulation.py index bc8cbf84b..888773308 100644 --- a/tidy3d/components/simulation.py +++ b/tidy3d/components/simulation.py @@ -2483,15 +2483,30 @@ def _projection_monitors_homogeneous(cls, val, values): structures = values.get("structures") or [] total_structures = [structure_bg] + list(structures) - for monitor in val: - if isinstance(monitor, (AbstractFieldProjectionMonitor, DiffractionMonitor)): - mediums = Scene.intersecting_media(monitor, total_structures) - # make sure there is no more than one medium in the returned list - if len(mediums) > 1: - raise SetupError( - f"{len(mediums)} different mediums detected on plane " - f"intersecting a {monitor.type}. Plane must be homogeneous." - ) + with log as consolidated_logger: + for monitor_ind, monitor in enumerate(val): + if isinstance(monitor, (AbstractFieldProjectionMonitor, DiffractionMonitor)): + mediums = Scene.intersecting_media(monitor, total_structures) + # make sure there is no more than one medium in the returned list + if len(mediums) > 1: + raise SetupError( + f"{len(mediums)} different mediums detected on plane " + f"intersecting a {monitor.type}. Plane must be homogeneous." + ) + # 0 medium, something is wrong + if len(mediums) < 1: + raise SetupError( + f"No medium detected on plane intersecting a {monitor.type}, " + "indicating an unexpected error. Please create a github issue so " + "that the problem can be investigated." + ) + # 1 medium, check if the medium is spatially uniform + if not list(mediums)[0].is_spatially_uniform: + consolidated_logger.warning( + f"Nonuniform custom medium detected on plane intersecting a {monitor.type}. " + "Plane must be homogeneous. Make sure custom medium is uniform on the plane.", + custom_loc=["monitors", monitor_ind], + ) return val @@ -2744,23 +2759,37 @@ def _source_homogeneous_isotropic(cls, val, values): total_structures = [structure_bg] + list(structures) # for each plane wave in the sources list - for source in val: - if isinstance(source, (PlaneWave, GaussianBeam, AstigmaticGaussianBeam)): - mediums = Scene.intersecting_media(source, total_structures) - # make sure there is no more than one medium in the returned list - if len(mediums) > 1: - raise SetupError( - f"{len(mediums)} different mediums detected on plane " - f"intersecting a {source.type} source. Plane must be homogeneous." - ) - if len(mediums) == 1 and isinstance( - list(mediums)[0], (AnisotropicMedium, FullyAnisotropicMedium) - ): - raise SetupError( - f"An anisotropic medium is detected on plane intersecting a {source.type} " - f"source. Injection of {source.type} into anisotropic media currently is " - "not supported." - ) + with log as consolidated_logger: + for source_id, source in enumerate(val): + if isinstance(source, (PlaneWave, GaussianBeam, AstigmaticGaussianBeam)): + mediums = Scene.intersecting_media(source, total_structures) + # make sure there is no more than one medium in the returned list + if len(mediums) > 1: + raise SetupError( + f"{len(mediums)} different mediums detected on plane " + f"intersecting a {source.type} source. Plane must be homogeneous." + ) + # 0 medium, something is wrong + if len(mediums) < 1: + raise SetupError( + f"No medium detected on plane intersecting a {source.type}, " + "indicating an unexpected error. Please create a github issue so " + "that the problem can be investigated." + ) + if isinstance(list(mediums)[0], (AnisotropicMedium, FullyAnisotropicMedium)): + raise SetupError( + f"An anisotropic medium is detected on plane intersecting a {source.type} " + f"source. Injection of {source.type} into anisotropic media currently is " + "not supported." + ) + + # check if the medium is spatially uniform + if not list(mediums)[0].is_spatially_uniform: + consolidated_logger.warning( + f"Nonuniform custom medium detected on plane intersecting a {source.type}. " + "Plane must be homogeneous. Make sure custom medium is uniform on the plane.", + custom_loc=["sources", source_id], + ) return val