From 8af9c099740f1648595c7f3d33d4c12375cfa8b1 Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Thu, 24 Feb 2022 10:53:20 +0100 Subject: [PATCH 1/8] started with tests --- nidn/tests/fdtd_test.py | 6 +++ nidn/tests/init_fdtd_grid_test.py | 83 +++++++++++++++++++++++++++++++ 2 files changed, 89 insertions(+) create mode 100644 nidn/tests/fdtd_test.py create mode 100644 nidn/tests/init_fdtd_grid_test.py diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py new file mode 100644 index 0000000..161c68a --- /dev/null +++ b/nidn/tests/fdtd_test.py @@ -0,0 +1,6 @@ +def test_fdtd_simulation(): + pass + + +if __name__ == "__main__": + test_fdtd_simulation() diff --git a/nidn/tests/init_fdtd_grid_test.py b/nidn/tests/init_fdtd_grid_test.py new file mode 100644 index 0000000..248d516 --- /dev/null +++ b/nidn/tests/init_fdtd_grid_test.py @@ -0,0 +1,83 @@ +from dotmap import DotMap +from numpy import source + +from ..fdtd_integration.init_fdtd import init_fdtd +import torch +from ..materials.layer_builder import LayerBuilder +from ..trcwa.compute_target_frequencies import compute_target_frequencies +from ..utils.load_default_cfg import load_default_cfg + + +def test_single_uniform_layer(): + # Create grid with uniform layer + cfg = load_default_cfg() + cfg.N_layers = 1 + cfg.N_freq = 1 + eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) + cfg.Nx = 1 + cfg.Ny = 1 + cfg.physical_wavelength_range[0] = 1e-5 + cfg.physical_wavelength_range[1] = 1e-5 + cfg.target_frequencies = compute_target_frequencies( + cfg.physical_wavelength_range[0], + cfg.physical_wavelength_range[1], + cfg.N_freq, + "linear", + ) + layer_builder = LayerBuilder(cfg) + eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") + grid, transmission_detector, reflection_detetctor = init_fdtd( + cfg, include_object=True, wavelength=1e-5, permittivity=eps_grid + ) + # Check that it was made properly + assert len(grid.objects) == 1 + assert grid.objects[0].permittivity == eps_grid[0, 0, 0, 0].real + assert len(grid.detectors) == 2 + assert len(grid.sources) == 1 + assert len(grid.boundaries) >= 2 + + +def test_multiple_uniform_layers(): + # Create grid with multiple uniform layer + cfg = DotMap() + cfg.N_layers = 4 + cfg.N_freq = 1 + eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) + cfg.Nx = 1 + cfg.Ny = 1 + cfg.physical_wavelength_range[0] = 1e-5 + cfg.physical_wavelength_range[1] = 1e-5 + cfg.target_frequencies = compute_target_frequencies( + cfg.physical_wavelength_range[0], + cfg.physical_wavelength_range[1], + cfg.N_freq, + "linear", + ) + cfg.PER_LAYER_THICKNESS = [1.0, 1.0, 1.0, 1.0] + cfg.FDTD_grid = [5.0, 2.0, 1.0] + layer_builder = LayerBuilder(cfg) + eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") + eps_grid[:, :, 1, :] = layer_builder.build_uniform_layer("zirconium") + eps_grid[:, :, 2, :] = layer_builder.build_uniform_layer("gallium_arsenide") + eps_grid[:, :, 3, :] = layer_builder.build_uniform_layer("silicon_nitride") + grid, transmission_detector, reflection_detetctor = init_fdtd( + cfg, include_object=True, wavelength=1e-5, permittivity=eps_grid + ) + # Check that it was made properly + assert len(grid.objects) == 4 + for i in range(4): + assert grid.objects[i].permittivity == eps_grid[0, 0, i, 0].real + assert len(grid.detectors) == 2 + assert len(grid.sources) == 1 + # If periodic boundaries in both x and y, it is two, if pml in x and periodic in y there is 3 and 4 if pml in both directions (I think) + assert len(grid.boundaries) >= 2 + + +def test_single_patterned_layer(): + # TODO: Test patterned layer, must be implemented first + pass + + +if __name__ == "__main__": + test_single_uniform_layer() + test_multiple_uniform_layers() From 6aab3c0461757c85eee9077bf96b32307c36fe47 Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Thu, 24 Feb 2022 13:24:26 +0100 Subject: [PATCH 2/8] Added tests for FDTD,a s well as validate added FDTD variables in cfg --- nidn/tests/fdtd_test.py | 32 +++++++++- nidn/tests/init_fdtd_grid_test.py | 19 +++--- nidn/training/utils/validate_config.py | 74 ++++++++++++++++++++---- nidn/utils/resources/default_config.toml | 4 +- 4 files changed, 104 insertions(+), 25 deletions(-) diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py index 161c68a..6fc8436 100644 --- a/nidn/tests/fdtd_test.py +++ b/nidn/tests/fdtd_test.py @@ -1,5 +1,35 @@ +import torch + +from ..fdtd_integration.init_fdtd import init_fdtd +from ..materials.layer_builder import LayerBuilder +from ..trcwa.compute_target_frequencies import compute_target_frequencies +from ..utils.load_default_cfg import load_default_cfg +from ..utils.compute_spectrum import compute_spectrum + + def test_fdtd_simulation(): - pass + # Create grid with uniform layer + cfg = load_default_cfg() + cfg.N_layers = 1 + cfg.N_freq = 1 + eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) + cfg.Nx = 1 + cfg.Ny = 1 + # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling + cfg.physical_wavelength_range[0] = 1e-6 + cfg.physical_wavelength_range[1] = 2e-6 + cfg.target_frequencies = compute_target_frequencies( + cfg.physical_wavelength_range[0], + cfg.physical_wavelength_range[1], + cfg.N_freq, + "linear", + ) + cfg.solver = "FDTD" + layer_builder = LayerBuilder(cfg) + eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") + transmission_spectrum, reflection_spectrum = compute_spectrum(eps_grid, cfg) + assert sum(transmission_spectrum) > 0 + assert sum(reflection_spectrum) > 0 if __name__ == "__main__": diff --git a/nidn/tests/init_fdtd_grid_test.py b/nidn/tests/init_fdtd_grid_test.py index 248d516..1444430 100644 --- a/nidn/tests/init_fdtd_grid_test.py +++ b/nidn/tests/init_fdtd_grid_test.py @@ -1,8 +1,6 @@ -from dotmap import DotMap -from numpy import source +import torch from ..fdtd_integration.init_fdtd import init_fdtd -import torch from ..materials.layer_builder import LayerBuilder from ..trcwa.compute_target_frequencies import compute_target_frequencies from ..utils.load_default_cfg import load_default_cfg @@ -16,8 +14,9 @@ def test_single_uniform_layer(): eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) cfg.Nx = 1 cfg.Ny = 1 - cfg.physical_wavelength_range[0] = 1e-5 - cfg.physical_wavelength_range[1] = 1e-5 + # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling + cfg.physical_wavelength_range[0] = 1e-6 + cfg.physical_wavelength_range[1] = 2e-6 cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], cfg.physical_wavelength_range[1], @@ -39,22 +38,22 @@ def test_single_uniform_layer(): def test_multiple_uniform_layers(): # Create grid with multiple uniform layer - cfg = DotMap() + cfg = load_default_cfg() cfg.N_layers = 4 cfg.N_freq = 1 eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) cfg.Nx = 1 cfg.Ny = 1 - cfg.physical_wavelength_range[0] = 1e-5 - cfg.physical_wavelength_range[1] = 1e-5 + cfg.physical_wavelength_range[ + 0 + ] = 1e-6 # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling + cfg.physical_wavelength_range[1] = 2e-6 cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], cfg.physical_wavelength_range[1], cfg.N_freq, "linear", ) - cfg.PER_LAYER_THICKNESS = [1.0, 1.0, 1.0, 1.0] - cfg.FDTD_grid = [5.0, 2.0, 1.0] layer_builder = LayerBuilder(cfg) eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") eps_grid[:, :, 1, :] = layer_builder.build_uniform_layer("zirconium") diff --git a/nidn/training/utils/validate_config.py b/nidn/training/utils/validate_config.py index 7684ae3..2bda4ce 100644 --- a/nidn/training/utils/validate_config.py +++ b/nidn/training/utils/validate_config.py @@ -36,8 +36,17 @@ def _validate_config(cfg: DotMap): "reg_loss_weight", "add_noise", "noise_scale", + "solver", "TRCWA_L_grid", "TRCWA_NG", + "FDTD_grid", + "FDTD_use_pointsource", + "FDTD_use_pulsesource", + "FDTD_pml_thickness", + "FDTD_source", + "FDTD_free_space_distance", + "FDTD_reflection_detector_x", + "FDTD_niter", "target_reflectance_spectrum", "target_transmittance_spectrum", "freq_distribution", @@ -60,6 +69,7 @@ def _validate_config(cfg: DotMap): "eps_oversampling", "seed", "TRCWA_NG", + "FDTD_niter", ] float_keys = [ "L", @@ -70,12 +80,24 @@ def _validate_config(cfg: DotMap): "siren_omega", "noise_scale", "reg_loss_weight", + "FDTD_free_space_distance", + "FDTD_reflection_detector_x", + "FDTD_pml_thickness", ] - boolean_keys = ["use_regularization_loss", "add_noise", "use_gpu", "avoid_zero_eps"] - string_keys = ["model_type", "type", "name", "freq_distribution"] + boolean_keys = [ + "use_regularization_loss", + "add_noise", + "use_gpu", + "avoid_zero_eps", + "FDTD_use_pulsesource", + "FDTD_use_pointsource", + ] + string_keys = ["model_type", "type", "name", "freq_distribution", "solver"] list_keys = [ "PER_LAYER_THICKNESS", "TRCWA_L_grid", + "FDTD_grid", + "FDTD_source", "target_reflectance_spectrum", "target_transmittance_spectrum", "physical_wavelength_range", @@ -111,6 +133,14 @@ def _validate_config(cfg: DotMap): if not (cfg.physical_wavelength_range[0] < cfg.physical_wavelength_range[1]): raise ValueError(f"physical_wavelength_range must be ordered from low to high") + if ( + cfg.FDTD_pml_thickness + cfg.FDTD_free_space_distance + < cfg.FDTD_reflection_detector_x + ) or (cfg.FDTD_reflection_detector_x < cfg.FDTD_pml_thickness): + raise ValueError( + f"Reflection detector must be placed in the free space before an eventual object, and after the pml layer" + ) + positive_value_keys = [ "L", "n_neurons", @@ -124,6 +154,10 @@ def _validate_config(cfg: DotMap): "noise_scale", "TRCWA_NG", "reg_loss_weight", + "FDTD_niter", + "FDTD_free_space_distance", + "FDTD_reflection_detector_x", + "FDTD_pml_thickness", ] for key in positive_value_keys: if not (cfg[key] > 0): @@ -142,14 +176,25 @@ def _validate_config(cfg: DotMap): if not cfg.TRCWA_L_grid[0][0] > cfg.TRCWA_L_grid[0][1]: raise ValueError(f"TRCWA_L_grid dim0 must be ordered from high to low") - if not all(cfg.TRCWA_L_grid) >= 0: - raise ValueError(f"TRCWA_L_grid must be positive") - - if not all(cfg.target_transmittance_spectrum) >= 0: - raise ValueError(f"target_transmittance_spectrum must be positive") + all_positive_list_keys = [ + "TRCWA_L_grid", + "PER_LAYER_THICKNESS", + "FDTD_grid", + "FDTD_source", + ] + all_positive_or_zero_list_keys = [ + "target_transmittance_spectrum", + "target_reflectance_spectrum", + ] - if not all(cfg.target_reflectance_spectrum) >= 0: - raise ValueError(f"target_reflectance_spectrum must be positive") + for key in all_positive_list_keys: + if not (all(cfg[key]) > 0.0): + raise ValueError(f"All elements in {key} must be a positive integer") + for key in all_positive_or_zero_list_keys: + if not (all(cfg[key]) >= 0.0): + raise ValueError( + f"All elements in {key} must be a positive integer or zero" + ) if not len(cfg.target_transmittance_spectrum) == cfg.N_freq: raise ValueError(f"target_transmittance_spectrum must have length N_freq") @@ -157,9 +202,6 @@ def _validate_config(cfg: DotMap): if not len(cfg.target_reflectance_spectrum) == cfg.N_freq: raise ValueError(f"target_reflectance_spectrum must have length N_freq") - if not (all(cfg.TRCWA_PER_LAYER_THICKNESS) > 0.0): - raise ValueError(f"thickness must a positive number.") - if not ( len(cfg.PER_LAYER_THICKNESS) == cfg.N_layers or len(cfg.PER_LAYER_THICKNESS) == 1 @@ -168,3 +210,11 @@ def _validate_config(cfg: DotMap): if not (cfg.freq_distribution == "linear" or cfg.freq_distribution == "log"): raise ValueError(f"freq_distribution must be either 'linear' or 'log'") + + if not len(cfg.FDTD_grid) == 3: + raise ValueError(f"FDTD_grid must me 3-dimentional") + + if not (len(cfg.FDTD_source) == 2 or len(cfg.FDTD_source) == 3): + raise ValueError( + f"The FDTD source needs either 2- or 3-dimensional coordinaets" + ) diff --git a/nidn/utils/resources/default_config.toml b/nidn/utils/resources/default_config.toml index e38e0d0..2fb35c7 100644 --- a/nidn/utils/resources/default_config.toml +++ b/nidn/utils/resources/default_config.toml @@ -29,7 +29,7 @@ imag_min_eps = 0.0 imag_max_eps = 3.0 # Simulation type -solver = "FDTD" # Options: FDTD, TRCWA +solver = "TRCWA" # Options: FDTD, TRCWA # Grid dimensions Nx = 1 @@ -48,7 +48,7 @@ FDTD_use_pointsource = false # If the source should be a point source in the FDT FDTD_use_pulsesource = false # If the source should be a single pulse in the FDTD simulations, otherwise set as continuous wave FDTD_pml_thickness = 1.5 # Thickness of PML layerin FDTD simulations, set in FDTD unit magnitudes. Perfectly Matched Layer are boundaries used in the x-direction, design to absrb all radiation FDTD_source = [1.5,1.0] # Coordinates of the source used in FDTD simulations, in FDTD unit magnitudes -FDTD_free_space_distance = 1 # The thickness of the free space layer before and after the material layers, in FDTD unit magnitudes +FDTD_free_space_distance = 1.0 # The thickness of the free space layer before and after the material layers, in FDTD unit magnitudes FDTD_reflection_detector_x = 2.4 # X-coordinates of the reflection detector for the FDTD simulation, in FDTD unit magnitudes FDTD_niter = 200 # Number of iterations / timesteps to run the FDTD for From 4f86c7eacde25ae2765b53c1aca89e8bfa8d398b Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Thu, 24 Feb 2022 15:11:54 +0100 Subject: [PATCH 3/8] Added plot_spectra updates --- nidn/plots/plot_spectra.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/nidn/plots/plot_spectra.py b/nidn/plots/plot_spectra.py index 74b13c8..8f72398 100644 --- a/nidn/plots/plot_spectra.py +++ b/nidn/plots/plot_spectra.py @@ -3,7 +3,7 @@ from matplotlib import pyplot as plt from matplotlib.ticker import FormatStrFormatter from ..utils.convert_units import freq_to_wl, wl_to_phys_wl -from ..trcwa.compute_spectrum import compute_spectrum +from ..utils.compute_spectrum import compute_spectrum from ..training.model.model_to_eps_grid import model_to_eps_grid From 1fa05cef58ab8bfd576bfda0cfd558391e3f43d1 Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Mon, 28 Feb 2022 16:36:10 +0100 Subject: [PATCH 4/8] Changed structure of tests and added more thourough testing, along with updating validate config and added some more logging. Also fixed unnoticed errors in transmission detector placement --- .../fdtd_integration/compute_spectrum_fdtd.py | 9 +- nidn/fdtd_integration/init_fdtd.py | 9 +- nidn/tests/calculate_coefficient_test.py | 25 +++ nidn/tests/fdtd_test.py | 195 +++++++++++++++++- nidn/tests/init_fdtd_grid_test.py | 82 -------- nidn/training/utils/validate_config.py | 10 +- nidn/utils/resources/default_config.toml | 4 +- 7 files changed, 226 insertions(+), 108 deletions(-) create mode 100644 nidn/tests/calculate_coefficient_test.py delete mode 100644 nidn/tests/init_fdtd_grid_test.py diff --git a/nidn/fdtd_integration/compute_spectrum_fdtd.py b/nidn/fdtd_integration/compute_spectrum_fdtd.py index 91feec9..96a27f2 100644 --- a/nidn/fdtd_integration/compute_spectrum_fdtd.py +++ b/nidn/fdtd_integration/compute_spectrum_fdtd.py @@ -23,9 +23,10 @@ def compute_spectrum_fdtd(permittivity, cfg: DotMap): transmission_spectrum = [] reflection_spectrum = [] physical_wavelengths, norm_freq = get_frequency_points(cfg) - logger.debug("Wavelenghts in spectrum") + logger.debug("Wavelenghts in spectrum : ") logger.debug(physical_wavelengths) - + logger.debug("Number of layers: ") + logger.debug(len(permittivity[0,0,:,0])) # For each wavelength, calculate transmission and reflection coefficents disable_progress_bar = logger._core.min_level >= 20 for i in tqdm(range(len(physical_wavelengths)), disable=disable_progress_bar): @@ -59,8 +60,6 @@ def compute_spectrum_fdtd(permittivity, cfg: DotMap): ) transmission_signal.append(transmission_material) reflection_signal.append(reflection_material) - time = [i for i in range(len(transmission_signal[0]))] - # Calculate transmission and reflection coefficients, # by using the signals from the free space simulation and the material simulation ( @@ -120,7 +119,7 @@ def _get_abs_value_from_3D_signal(signal): abs_value = [] for i in range(len(signal)): abs_value.append( - sqrt(tensor(signal[i][0] ** 2 + signal[i][1] ** 2 + signal[i][2] ** 2)) + sqrt(signal[i][0] ** 2 + signal[i][1] ** 2 + signal[i][2] ** 2) ) return abs_value diff --git a/nidn/fdtd_integration/init_fdtd.py b/nidn/fdtd_integration/init_fdtd.py index 506c8d9..55686e8 100644 --- a/nidn/fdtd_integration/init_fdtd.py +++ b/nidn/fdtd_integration/init_fdtd.py @@ -38,13 +38,14 @@ def init_fdtd(cfg: DotMap, include_object, wavelength, permittivity): cfg.FDTD_reflection_detector_x + cfg.N_layers * cfg.PER_LAYER_THICKNESS[0] ) + + 2 ), int(cfg.FDTD_reflection_detector_x * scaling), ) grid = _add_source( grid, - int(cfg.FDTD_source[0] * scaling), - int(cfg.FDTD_source[1] * scaling), + int(cfg.FDTD_source_position[0] * scaling), + int(cfg.FDTD_source_position[1] * scaling), wavelength / SPEED_OF_LIGHT, cfg.FDTD_use_pulsesource, cfg.FDTD_use_pointsource, @@ -116,7 +117,7 @@ def _add_source(grid, source_x, source_y, period, use_pulse_source, use_point_so if use_point_source: grid[source_x, source_y, 0] = fdtd.PointSource( period=period, - name="linesource", + name="pointsource", pulse=use_pulse_source, cycle=1, hanning_dt=1e-15, @@ -158,7 +159,7 @@ def _add_object(grid, object_start_x, object_end_x, permittivity, frequency): Returns: fdtd.Grid: The grid with the added object """ - # Not sure whether the conductivity should be relative or absolute, i.e. if it should be multiplied with EPS_0. + # Not sure whether the conductivity should be relative or absolute, i.e. if it should be multiplied with EPS_0. Multiplied with 2pi to get w(angular frequency)? # Since the permittivity is set to 1 for the free space grid, I'll leave it at an relative value for now. Also, permittivity for object is relative. grid[object_start_x:object_end_x, :, :] = fdtd.AbsorbingObject( permittivity=permittivity.real, diff --git a/nidn/tests/calculate_coefficient_test.py b/nidn/tests/calculate_coefficient_test.py new file mode 100644 index 0000000..14d5c78 --- /dev/null +++ b/nidn/tests/calculate_coefficient_test.py @@ -0,0 +1,25 @@ +from numpy import sin, pi +from ..fdtd_integration.calculate_transmission_reflection_coefficients import ( + calculate_transmission_reflection_coefficients, +) +from ..utils.load_default_cfg import load_default_cfg + + +def test_calculate_coefficient(): + cfg = load_default_cfg() + time_points = [0.002 * pi * i for i in range(1000)] + signal_a = 2 * sin(time_points) + signal_b = sin(time_points) + signal_array = [signal_a, signal_b] + ( + transmission_coefficient_ms, + reflection_coefficient_ms, + ) = calculate_transmission_reflection_coefficients( + signal_array, signal_array, "mean square", cfg + ) + # TODO: Add test for fdtd, when the method is complete + assert transmission_coefficient_ms - 0.25 == 0 + + +if __name__ == "__main__": + test_calculate_coefficient() diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py index 6fc8436..b4fac85 100644 --- a/nidn/tests/fdtd_test.py +++ b/nidn/tests/fdtd_test.py @@ -1,4 +1,7 @@ -import torch +from torch import zeros, tensor, cfloat +from numpy import subtract + +from nidn.utils.global_constants import SPEED_OF_LIGHT from ..fdtd_integration.init_fdtd import init_fdtd from ..materials.layer_builder import LayerBuilder @@ -7,17 +10,131 @@ from ..utils.compute_spectrum import compute_spectrum -def test_fdtd_simulation(): - # Create grid with uniform layer +def test_fdtd_grid_creation(): + """Test that the simulation is created in the right way, whith the correcto objects and correct realtive placement of them""" + # Create grid with multiple uniform layer cfg = load_default_cfg() - cfg.N_layers = 1 + cfg.N_layers = 4 cfg.N_freq = 1 - eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) + eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) cfg.Nx = 1 cfg.Ny = 1 # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling - cfg.physical_wavelength_range[0] = 1e-6 - cfg.physical_wavelength_range[1] = 2e-6 + cfg.target_frequencies = compute_target_frequencies( + cfg.physical_wavelength_range[0], + cfg.physical_wavelength_range[1], + cfg.N_freq, + "linear", + ) + layer_builder = LayerBuilder(cfg) + eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") + eps_grid[:, :, 1, :] = layer_builder.build_uniform_layer("zirconium") + eps_grid[:, :, 2, :] = layer_builder.build_uniform_layer("gallium_arsenide") + eps_grid[:, :, 3, :] = layer_builder.build_uniform_layer("silicon_nitride") + grid, transmission_detector, reflection_detetctor = init_fdtd( + cfg, + include_object=True, + wavelength=cfg.physical_wavelength_range[0], + permittivity=eps_grid, + ) + # Check that it was made properly + assert len(grid.objects) == 4 + for i in range(len(grid.objects)): + assert grid.objects[i].permittivity == eps_grid[0, 0, i, 0].real + assert ( + grid.objects[i].conductivity[0, 0, 0, 0] + - ( + eps_grid[0, 0, i, 0].imag + * SPEED_OF_LIGHT + / cfg.physical_wavelength_range[0] + ) + < 1e-8 + ) + + assert len(grid.detectors) == 2 + # Check that the reflection detector is placed before the first layer, and the transmission detector is placed after the last layer + assert transmission_detector.x[0] >= grid.objects[-1].x.stop + assert reflection_detetctor.x[0] <= grid.objects[0].x.start + + assert len(grid.sources) == 1 + # If periodic boundaries in both x and y, it is two, if pml in x and periodic in y there is 3 and 4 if pml in both directions (I think) + assert len(grid.boundaries) >= 2 + + +def test_fdtd_simulation_single_layer(): + """Test that checks that the calculate_spectrum function returns the correct spectrum for a single layer""" + # Create grid with uniform layer + cfg = load_default_cfg() + cfg.N_layers = 1 + eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) + # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling + cfg.target_frequencies = compute_target_frequencies( + cfg.physical_wavelength_range[0], + cfg.physical_wavelength_range[1], + cfg.N_freq, + "linear", + ) + cfg.solver = "FDTD" + layer_builder = LayerBuilder(cfg) + eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") + transmission_spectrum, reflection_spectrum = compute_spectrum(eps_grid, cfg) + validated_transmission_spectrum = [ + tensor(0.0), + tensor(0.0), + tensor(2.0448e-308), + tensor(2.5918e-277), + tensor(3.8102e-132), + tensor(0.6251), + tensor(0.5362), + tensor(0.6362), + tensor(0.7259), + tensor(0.6041), + tensor(0.5613), + tensor(0.6453), + tensor(0.6915), + tensor(0.5796), + tensor(0.4658), + tensor(0.5885), + tensor(0.5805), + tensor(0.4043), + tensor(0.4957), + tensor(0.4211), + ] + validated_reflection_spectrum = [ + tensor(0.9515), + tensor(0.8697), + tensor(0.8118), + tensor(0.7570), + tensor(0.6983), + tensor(0.2333), + tensor(0.3495), + tensor(0.2044), + tensor(0.0673), + tensor(0.2649), + tensor(0.3565), + tensor(0.2003), + tensor(0.0904), + tensor(0.2549), + tensor(0.3180), + tensor(0.1525), + tensor(0.1539), + tensor(0.3921), + tensor(0.2173), + tensor(0.3437), + ] + assert all( + abs(subtract(transmission_spectrum, validated_transmission_spectrum)) < 1e-4 + ) + assert all(abs(subtract(reflection_spectrum, validated_reflection_spectrum)) < 1e-4) + + +def test_fdtd_simulation_four_layers(): + """Test that checks that the calculate_spectrum function returns the correct spectrum for a simulation with four layers""" + # Create grid with four uniform layer + cfg = load_default_cfg() + cfg.N_layers = 4 + cfg.FDTD_niter = 400 + eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], cfg.physical_wavelength_range[1], @@ -27,10 +144,68 @@ def test_fdtd_simulation(): cfg.solver = "FDTD" layer_builder = LayerBuilder(cfg) eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") + eps_grid[:, :, 1, :] = layer_builder.build_uniform_layer("zirconium") + eps_grid[:, :, 2, :] = layer_builder.build_uniform_layer("gallium_arsenide") + eps_grid[:, :, 3, :] = layer_builder.build_uniform_layer("silicon_nitride") transmission_spectrum, reflection_spectrum = compute_spectrum(eps_grid, cfg) - assert sum(transmission_spectrum) > 0 - assert sum(reflection_spectrum) > 0 + validated_transmission_spectrum = [ + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + tensor(0.0), + ] + validated_reflection_spectrum = [ + tensor(0.9515), + tensor(0.8781), + tensor(0.8077), + tensor(0.7559), + tensor(0.7042), + tensor(0.5871), + tensor(0.5489), + tensor(0.6336), + tensor(0.4941), + tensor(0.6071), + tensor(0.4161), + tensor(0.6353), + tensor(0.5313), + tensor(0.6544), + tensor(0.6774), + tensor(0.4997), + tensor(0.6510), + tensor(0.4196), + tensor(0.5047), + tensor(0.2888), + ] + assert all( + abs(subtract(transmission_spectrum, validated_transmission_spectrum)) < 1e-4 + ) + assert all(abs(subtract(reflection_spectrum, validated_reflection_spectrum)) < 1e-4) + + +def test_single_patterned_layer(): + """Test that a pattern layer returns teh correct spectrum""" + # TODO: Test patterned layer, must be implemented first + pass if __name__ == "__main__": - test_fdtd_simulation() + test_fdtd_grid_creation() + test_fdtd_simulation_single_layer() + test_fdtd_simulation_four_layers() + test_single_patterned_layer() diff --git a/nidn/tests/init_fdtd_grid_test.py b/nidn/tests/init_fdtd_grid_test.py deleted file mode 100644 index 1444430..0000000 --- a/nidn/tests/init_fdtd_grid_test.py +++ /dev/null @@ -1,82 +0,0 @@ -import torch - -from ..fdtd_integration.init_fdtd import init_fdtd -from ..materials.layer_builder import LayerBuilder -from ..trcwa.compute_target_frequencies import compute_target_frequencies -from ..utils.load_default_cfg import load_default_cfg - - -def test_single_uniform_layer(): - # Create grid with uniform layer - cfg = load_default_cfg() - cfg.N_layers = 1 - cfg.N_freq = 1 - eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) - cfg.Nx = 1 - cfg.Ny = 1 - # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling - cfg.physical_wavelength_range[0] = 1e-6 - cfg.physical_wavelength_range[1] = 2e-6 - cfg.target_frequencies = compute_target_frequencies( - cfg.physical_wavelength_range[0], - cfg.physical_wavelength_range[1], - cfg.N_freq, - "linear", - ) - layer_builder = LayerBuilder(cfg) - eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") - grid, transmission_detector, reflection_detetctor = init_fdtd( - cfg, include_object=True, wavelength=1e-5, permittivity=eps_grid - ) - # Check that it was made properly - assert len(grid.objects) == 1 - assert grid.objects[0].permittivity == eps_grid[0, 0, 0, 0].real - assert len(grid.detectors) == 2 - assert len(grid.sources) == 1 - assert len(grid.boundaries) >= 2 - - -def test_multiple_uniform_layers(): - # Create grid with multiple uniform layer - cfg = load_default_cfg() - cfg.N_layers = 4 - cfg.N_freq = 1 - eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) - cfg.Nx = 1 - cfg.Ny = 1 - cfg.physical_wavelength_range[ - 0 - ] = 1e-6 # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling - cfg.physical_wavelength_range[1] = 2e-6 - cfg.target_frequencies = compute_target_frequencies( - cfg.physical_wavelength_range[0], - cfg.physical_wavelength_range[1], - cfg.N_freq, - "linear", - ) - layer_builder = LayerBuilder(cfg) - eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") - eps_grid[:, :, 1, :] = layer_builder.build_uniform_layer("zirconium") - eps_grid[:, :, 2, :] = layer_builder.build_uniform_layer("gallium_arsenide") - eps_grid[:, :, 3, :] = layer_builder.build_uniform_layer("silicon_nitride") - grid, transmission_detector, reflection_detetctor = init_fdtd( - cfg, include_object=True, wavelength=1e-5, permittivity=eps_grid - ) - # Check that it was made properly - assert len(grid.objects) == 4 - for i in range(4): - assert grid.objects[i].permittivity == eps_grid[0, 0, i, 0].real - assert len(grid.detectors) == 2 - assert len(grid.sources) == 1 - # If periodic boundaries in both x and y, it is two, if pml in x and periodic in y there is 3 and 4 if pml in both directions (I think) - assert len(grid.boundaries) >= 2 - - -def test_single_patterned_layer(): - # TODO: Test patterned layer, must be implemented first - pass - - -if __name__ == "__main__": - test_single_uniform_layer() - test_multiple_uniform_layers() diff --git a/nidn/training/utils/validate_config.py b/nidn/training/utils/validate_config.py index 2bda4ce..f6123db 100644 --- a/nidn/training/utils/validate_config.py +++ b/nidn/training/utils/validate_config.py @@ -43,7 +43,7 @@ def _validate_config(cfg: DotMap): "FDTD_use_pointsource", "FDTD_use_pulsesource", "FDTD_pml_thickness", - "FDTD_source", + "FDTD_source_position", "FDTD_free_space_distance", "FDTD_reflection_detector_x", "FDTD_niter", @@ -97,7 +97,7 @@ def _validate_config(cfg: DotMap): "PER_LAYER_THICKNESS", "TRCWA_L_grid", "FDTD_grid", - "FDTD_source", + "FDTD_source_position", "target_reflectance_spectrum", "target_transmittance_spectrum", "physical_wavelength_range", @@ -180,7 +180,7 @@ def _validate_config(cfg: DotMap): "TRCWA_L_grid", "PER_LAYER_THICKNESS", "FDTD_grid", - "FDTD_source", + "FDTD_source_position", ] all_positive_or_zero_list_keys = [ "target_transmittance_spectrum", @@ -214,7 +214,7 @@ def _validate_config(cfg: DotMap): if not len(cfg.FDTD_grid) == 3: raise ValueError(f"FDTD_grid must me 3-dimentional") - if not (len(cfg.FDTD_source) == 2 or len(cfg.FDTD_source) == 3): + if not (len(cfg.FDTD_source_position) == 2 or len(cfg.FDTD_source_position) == 3): raise ValueError( - f"The FDTD source needs either 2- or 3-dimensional coordinaets" + f"The FDTD source needs either 2- or 3-dimensional coordinates" ) diff --git a/nidn/utils/resources/default_config.toml b/nidn/utils/resources/default_config.toml index 2fb35c7..a7d3e9f 100644 --- a/nidn/utils/resources/default_config.toml +++ b/nidn/utils/resources/default_config.toml @@ -46,8 +46,8 @@ TRCWA_NG = 11 # Truncation order (actual number might be smaller) FDTD_grid = [5.0,2.0,1.0] # Dimensions of the grid used in FDTD simulations, in FDTD unit magnitudes FDTD_use_pointsource = false # If the source should be a point source in the FDTD simulaitons, otherwise set as linesource FDTD_use_pulsesource = false # If the source should be a single pulse in the FDTD simulations, otherwise set as continuous wave -FDTD_pml_thickness = 1.5 # Thickness of PML layerin FDTD simulations, set in FDTD unit magnitudes. Perfectly Matched Layer are boundaries used in the x-direction, design to absrb all radiation -FDTD_source = [1.5,1.0] # Coordinates of the source used in FDTD simulations, in FDTD unit magnitudes +FDTD_pml_thickness = 1.5 # Thickness of PML layerin FDTD simulations, set in FDTD unit magnitudes. Perfectly Matched Layer are boundaries used in the x-direction, design to absorb all radiation +FDTD_source_position = [1.5,1.0] # Coordinates of the source used in FDTD simulations, in FDTD unit magnitudes FDTD_free_space_distance = 1.0 # The thickness of the free space layer before and after the material layers, in FDTD unit magnitudes FDTD_reflection_detector_x = 2.4 # X-coordinates of the reflection detector for the FDTD simulation, in FDTD unit magnitudes FDTD_niter = 200 # Number of iterations / timesteps to run the FDTD for From 3b4cc39c225596b0794f429c4e38f5f886f061a1 Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Tue, 1 Mar 2022 00:13:33 +0100 Subject: [PATCH 5/8] Reduced number of frequencies in fdtd simulations used in the tests --- nidn/tests/fdtd_test.py | 80 ++++++----------------------------------- 1 file changed, 11 insertions(+), 69 deletions(-) diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py index b4fac85..5d9ebb1 100644 --- a/nidn/tests/fdtd_test.py +++ b/nidn/tests/fdtd_test.py @@ -65,6 +65,7 @@ def test_fdtd_simulation_single_layer(): """Test that checks that the calculate_spectrum function returns the correct spectrum for a single layer""" # Create grid with uniform layer cfg = load_default_cfg() + cfg.N_freq = 5 cfg.N_layers = 1 eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling @@ -80,46 +81,16 @@ def test_fdtd_simulation_single_layer(): transmission_spectrum, reflection_spectrum = compute_spectrum(eps_grid, cfg) validated_transmission_spectrum = [ tensor(0.0), - tensor(0.0), - tensor(2.0448e-308), - tensor(2.5918e-277), - tensor(3.8102e-132), - tensor(0.6251), - tensor(0.5362), - tensor(0.6362), - tensor(0.7259), - tensor(0.6041), - tensor(0.5613), - tensor(0.6453), - tensor(0.6915), - tensor(0.5796), - tensor(0.4658), - tensor(0.5885), - tensor(0.5805), - tensor(0.4043), - tensor(0.4957), + tensor(0.5564), + tensor(0.5902), + tensor(0.4664), tensor(0.4211), ] validated_reflection_spectrum = [ tensor(0.9515), - tensor(0.8697), - tensor(0.8118), - tensor(0.7570), - tensor(0.6983), - tensor(0.2333), - tensor(0.3495), - tensor(0.2044), - tensor(0.0673), - tensor(0.2649), - tensor(0.3565), - tensor(0.2003), - tensor(0.0904), - tensor(0.2549), - tensor(0.3180), - tensor(0.1525), - tensor(0.1539), - tensor(0.3921), - tensor(0.2173), + tensor(0.1605), + tensor(0.3508), + tensor(0.3171), tensor(0.3437), ] assert all( @@ -134,6 +105,7 @@ def test_fdtd_simulation_four_layers(): cfg = load_default_cfg() cfg.N_layers = 4 cfg.FDTD_niter = 400 + cfg.N_freq = 5 eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], @@ -154,42 +126,12 @@ def test_fdtd_simulation_four_layers(): tensor(0.0), tensor(0.0), tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), ] validated_reflection_spectrum = [ tensor(0.9515), - tensor(0.8781), - tensor(0.8077), - tensor(0.7559), - tensor(0.7042), - tensor(0.5871), - tensor(0.5489), - tensor(0.6336), - tensor(0.4941), - tensor(0.6071), - tensor(0.4161), - tensor(0.6353), - tensor(0.5313), - tensor(0.6544), - tensor(0.6774), - tensor(0.4997), - tensor(0.6510), - tensor(0.4196), - tensor(0.5047), + tensor(0.4005), + tensor(0.4919), + tensor(0.6812), tensor(0.2888), ] assert all( From 92875b6f557f60c39b62a37c70f91a7e17778f72 Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Wed, 2 Mar 2022 11:41:26 +0100 Subject: [PATCH 6/8] Removed external fdtd import, restructured customability of source, added more precision to tests and added validation of coefficents --- nidn/fdtd/detectors.py | 12 +-- ...te_transmission_reflection_coefficients.py | 18 ++++- .../fdtd_integration/compute_spectrum_fdtd.py | 10 +-- nidn/fdtd_integration/init_fdtd.py | 42 ++++++----- nidn/tests/calculate_coefficient_test.py | 11 +-- nidn/tests/fdtd_test.py | 73 +++++++++---------- nidn/training/utils/validate_config.py | 28 +++++-- nidn/utils/resources/default_config.toml | 6 +- 8 files changed, 118 insertions(+), 82 deletions(-) diff --git a/nidn/fdtd/detectors.py b/nidn/fdtd/detectors.py index 37304b3..f0f9998 100644 --- a/nidn/fdtd/detectors.py +++ b/nidn/fdtd/detectors.py @@ -14,7 +14,7 @@ # relative from .grid import Grid from .backend import backend as bd -from .constants import X, Y, Z +from .constants import FDTD_X, FDTD_Y, FDTD_Z ## Detector class LineDetector: @@ -441,20 +441,22 @@ def single_point_current(self, px, py, pz): # option for unit factor here? Units will get very complicated otherwise current_vector_1 = ( - (self.grid.H[px, py - 1, pz, X]) - (self.grid.H[px, py, pz, X]) + (self.grid.H[px, py - 1, pz, FDTD_X]) - (self.grid.H[px, py, pz, FDTD_X]) ) * self.grid.grid_spacing current_vector_2 = ( - (self.grid.H[px, py, pz, Y]) - (self.grid.H[px - 1, py, pz, Y]) + (self.grid.H[px, py, pz, FDTD_Y]) - (self.grid.H[px - 1, py, pz, FDTD_Y]) ) * self.grid.grid_spacing current_1 = current_vector_1 + current_vector_2 # current_1 = float(current.cpu()) current_vector_1 = ( - (self.grid.H[px, py - 1, pz - 1, X]) - (self.grid.H[px, py, pz - 1, X]) + (self.grid.H[px, py - 1, pz - 1, FDTD_X]) + - (self.grid.H[px, py, pz - 1, FDTD_X]) ) * self.grid.grid_spacing current_vector_2 += ( - (self.grid.H[px, py, pz - 1, Y]) - (self.grid.H[px - 1, py, pz - 1, Y]) + (self.grid.H[px, py, pz - 1, FDTD_Y]) + - (self.grid.H[px - 1, py, pz - 1, FDTD_Y]) ) * self.grid.grid_spacing # current_2 = float(current_2.cpu()) current_2 = current_vector_1 + current_vector_2 diff --git a/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py b/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py index f5205c6..b076f98 100644 --- a/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py +++ b/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py @@ -1,6 +1,8 @@ from dotmap import DotMap from torch.fft import rfft, rfftfreq -from torch import sqrt, tensor +import torch +from loguru import logger + from nidn.fdtd_integration.constants import FDTD_GRID_SCALE @@ -49,6 +51,20 @@ def calculate_transmission_reflection_coefficients( reflection_coefficient = _fft(true_reflection, cfg) / _fft( reflection_signals[0], cfg ) + + if not (transmission_coefficient >= 0 and transmission_coefficient <= 1): + raise ValueError( + f"The transmission coefficient is outsied of the possible range between 0 and 1" + ) + + if not (reflection_coefficient >= 0 and reflection_coefficient <= 1): + raise ValueError( + f"The reflection coefficient is outsied of the possible range between 0 and 1" + ) + if transmission_coefficient + reflection_coefficient > 1: + logger.warning( + f"The sum of the transmission and reflection coefficient is greater than 1, which is physically impossible" + ) return transmission_coefficient, reflection_coefficient diff --git a/nidn/fdtd_integration/compute_spectrum_fdtd.py b/nidn/fdtd_integration/compute_spectrum_fdtd.py index 96a27f2..618661a 100644 --- a/nidn/fdtd_integration/compute_spectrum_fdtd.py +++ b/nidn/fdtd_integration/compute_spectrum_fdtd.py @@ -1,6 +1,6 @@ from dotmap import DotMap from tqdm import tqdm -from torch import sqrt, tensor +import torch from loguru import logger from ..trcwa.get_frequency_points import get_frequency_points @@ -26,7 +26,7 @@ def compute_spectrum_fdtd(permittivity, cfg: DotMap): logger.debug("Wavelenghts in spectrum : ") logger.debug(physical_wavelengths) logger.debug("Number of layers: ") - logger.debug(len(permittivity[0,0,:,0])) + logger.debug(len(permittivity[0, 0, :, 0])) # For each wavelength, calculate transmission and reflection coefficents disable_progress_bar = logger._core.min_level >= 20 for i in tqdm(range(len(physical_wavelengths)), disable=disable_progress_bar): @@ -76,7 +76,7 @@ def compute_spectrum_fdtd(permittivity, cfg: DotMap): logger.debug("Reflection spectrum") logger.debug(reflection_spectrum) - return transmission_spectrum, reflection_spectrum + return torch.tensor(transmission_spectrum), torch.tensor(reflection_spectrum) def _get_detector_values(transmission_detector, reflection_detector): @@ -119,7 +119,7 @@ def _get_abs_value_from_3D_signal(signal): abs_value = [] for i in range(len(signal)): abs_value.append( - sqrt(signal[i][0] ** 2 + signal[i][1] ** 2 + signal[i][2] ** 2) + torch.sqrt(signal[i][0] ** 2 + signal[i][1] ** 2 + signal[i][2] ** 2) ) return abs_value @@ -135,7 +135,7 @@ def _average_along_detector(signal): """ avg = [] for e in signal: - s = [tensor(0.0), tensor(0.0), tensor(0.0)] + s = [torch.tensor(0.0), torch.tensor(0.0), torch.tensor(0.0)] for p in e: s[0] += p[0] / len(e) s[1] += p[1] / len(e) diff --git a/nidn/fdtd_integration/init_fdtd.py b/nidn/fdtd_integration/init_fdtd.py index 55686e8..1a2cdec 100644 --- a/nidn/fdtd_integration/init_fdtd.py +++ b/nidn/fdtd_integration/init_fdtd.py @@ -1,6 +1,6 @@ from dotmap import DotMap -import fdtd +from ..fdtd import * from ..utils.global_constants import EPS_0, SPEED_OF_LIGHT, UNIT_MAGNITUDE from .constants import FDTD_GRID_SCALE @@ -17,13 +17,13 @@ def init_fdtd(cfg: DotMap, include_object, wavelength, permittivity): Returns: fdtd:Grid: Grid with all the added object, ready to be run """ - fdtd.set_backend("torch") + set_backend("torch") scaling = UNIT_MAGNITUDE / (cfg.physical_wavelength_range[0] * FDTD_GRID_SCALE) x_grid_size = int( scaling * (cfg.FDTD_grid[0] + cfg.N_layers * cfg.PER_LAYER_THICKNESS[0]) ) y_grid_size = int(cfg.FDTD_grid[1] * scaling) - grid = fdtd.Grid( + grid = Grid( (x_grid_size, y_grid_size, 1), grid_spacing=cfg.physical_wavelength_range[0] * FDTD_GRID_SCALE, permittivity=1.0, @@ -38,17 +38,24 @@ def init_fdtd(cfg: DotMap, include_object, wavelength, permittivity): cfg.FDTD_reflection_detector_x + cfg.N_layers * cfg.PER_LAYER_THICKNESS[0] ) - + 2 + + cfg.FDTD_min_gridpoints_between_detectors ), int(cfg.FDTD_reflection_detector_x * scaling), ) + use_pulse = True + if cfg.FDTD_pulse_type == "pulse": + pass + elif cfg.FDTD_pulse_type == "continuous": + use_pulse = False + else: + raise ValueError(f' FDTD_pulse_type must either be "pulse" or "continuous"') grid = _add_source( grid, int(cfg.FDTD_source_position[0] * scaling), int(cfg.FDTD_source_position[1] * scaling), wavelength / SPEED_OF_LIGHT, - cfg.FDTD_use_pulsesource, - cfg.FDTD_use_pointsource, + use_pulse, + cfg.FDTD_source_type, ) if include_object: for i in range(cfg.N_layers): @@ -89,17 +96,17 @@ def _add_boundaries(grid, pml_thickness): fdtd.Grid: The grid with the added boundaries """ # Add PML boundary to the left side of the grid - grid[0:pml_thickness, :, :] = fdtd.PML(name="pml_xlow") + grid[0:pml_thickness, :, :] = PML(name="pml_xlow") # Add PML boundary to the right side of the grid - grid[-pml_thickness:, :, :] = fdtd.PML(name="pml_xhigh") + grid[-pml_thickness:, :, :] = PML(name="pml_xhigh") # Add periodic boundaries at both sides in y-direction - grid[:, 0, :] = fdtd.PeriodicBoundary(name="ybounds") + grid[:, 0, :] = PeriodicBoundary(name="ybounds") # Add periodic boundaries on both sides in z-direction. Only applicable for 3D grids # grid[:, :, 0] = fdtd.PeriodicBoundary(name="zbounds") return grid -def _add_source(grid, source_x, source_y, period, use_pulse_source, use_point_source): +def _add_source(grid, source_x, source_y, period, use_pulse_source, source_type): """Add a specified source to the fdtd grid Args: @@ -114,17 +121,18 @@ def _add_source(grid, source_x, source_y, period, use_pulse_source, use_point_so fdtd.Grid: The grid with the added source """ - if use_point_source: - grid[source_x, source_y, 0] = fdtd.PointSource( + if source_type == "point": + grid[source_x, source_y, 0] = PointSource( period=period, name="pointsource", pulse=use_pulse_source, cycle=1, hanning_dt=1e-15, ) + elif source_type == "line": + grid[source_x, :, 0] = LineSource(period=period, name="linesource") else: - grid[source_x, :, 0] = fdtd.LineSource(period=period, name="linesource") - + raise ValueError(f'FDTD_source_type must be either "line" or "point"') return grid @@ -140,8 +148,8 @@ def _add_detectors(grid, transmission_detector_x, reflection_detector_x): fdtd.Grid: The grid with the added detectors """ - transmission_detector = fdtd.LineDetector(name="t_detector") - reflection_detector = fdtd.LineDetector(name="r_detector") + transmission_detector = LineDetector(name="t_detector") + reflection_detector = LineDetector(name="r_detector") grid[transmission_detector_x, :, 0] = transmission_detector grid[reflection_detector_x, :, 0] = reflection_detector return grid, transmission_detector, reflection_detector @@ -161,7 +169,7 @@ def _add_object(grid, object_start_x, object_end_x, permittivity, frequency): """ # Not sure whether the conductivity should be relative or absolute, i.e. if it should be multiplied with EPS_0. Multiplied with 2pi to get w(angular frequency)? # Since the permittivity is set to 1 for the free space grid, I'll leave it at an relative value for now. Also, permittivity for object is relative. - grid[object_start_x:object_end_x, :, :] = fdtd.AbsorbingObject( + grid[object_start_x:object_end_x, :, :] = AbsorbingObject( permittivity=permittivity.real, conductivity=permittivity.imag * frequency, ) diff --git a/nidn/tests/calculate_coefficient_test.py b/nidn/tests/calculate_coefficient_test.py index 14d5c78..e60987b 100644 --- a/nidn/tests/calculate_coefficient_test.py +++ b/nidn/tests/calculate_coefficient_test.py @@ -1,4 +1,4 @@ -from numpy import sin, pi +import numpy as np from ..fdtd_integration.calculate_transmission_reflection_coefficients import ( calculate_transmission_reflection_coefficients, ) @@ -7,9 +7,9 @@ def test_calculate_coefficient(): cfg = load_default_cfg() - time_points = [0.002 * pi * i for i in range(1000)] - signal_a = 2 * sin(time_points) - signal_b = sin(time_points) + time_points = [0.002 * np.pi * i for i in range(1000)] + signal_a = 2 * np.sin(time_points) + signal_b = np.sin(time_points) signal_array = [signal_a, signal_b] ( transmission_coefficient_ms, @@ -17,8 +17,9 @@ def test_calculate_coefficient(): ) = calculate_transmission_reflection_coefficients( signal_array, signal_array, "mean square", cfg ) - # TODO: Add test for fdtd, when the method is complete + # TODO: Add test for fft, when the method is complete assert transmission_coefficient_ms - 0.25 == 0 + assert reflection_coefficient_ms - 0.25 == 0 if __name__ == "__main__": diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py index 5d9ebb1..587f202 100644 --- a/nidn/tests/fdtd_test.py +++ b/nidn/tests/fdtd_test.py @@ -1,5 +1,4 @@ -from torch import zeros, tensor, cfloat -from numpy import subtract +import torch from nidn.utils.global_constants import SPEED_OF_LIGHT @@ -11,15 +10,14 @@ def test_fdtd_grid_creation(): - """Test that the simulation is created in the right way, whith the correcto objects and correct realtive placement of them""" + """Test that the simulation is created in the right way, with the correct objects and correct relative placement of them""" # Create grid with multiple uniform layer cfg = load_default_cfg() cfg.N_layers = 4 cfg.N_freq = 1 - eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) + eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) cfg.Nx = 1 cfg.Ny = 1 - # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], cfg.physical_wavelength_range[1], @@ -67,8 +65,7 @@ def test_fdtd_simulation_single_layer(): cfg = load_default_cfg() cfg.N_freq = 5 cfg.N_layers = 1 - eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) - # Note: something went wrong when smalles wavelength was 1e-5, guess its the grid scaling + eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], cfg.physical_wavelength_range[1], @@ -79,24 +76,22 @@ def test_fdtd_simulation_single_layer(): layer_builder = LayerBuilder(cfg) eps_grid[:, :, 0, :] = layer_builder.build_uniform_layer("titanium_oxide") transmission_spectrum, reflection_spectrum = compute_spectrum(eps_grid, cfg) - validated_transmission_spectrum = [ - tensor(0.0), - tensor(0.5564), - tensor(0.5902), - tensor(0.4664), - tensor(0.4211), - ] - validated_reflection_spectrum = [ - tensor(0.9515), - tensor(0.1605), - tensor(0.3508), - tensor(0.3171), - tensor(0.3437), - ] + validated_transmission_spectrum = torch.tensor( + [ + 4.95440515e-176, + 6.98087684e-01, + 5.90159230e-01, + 4.66396898e-01, + 4.21121637e-01, + ] + ) + validated_reflection_spectrum = torch.tensor( + [0.95146948, 0.17358959, 0.35082144, 0.31708776, 0.34369091] + ) assert all( - abs(subtract(transmission_spectrum, validated_transmission_spectrum)) < 1e-4 + torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-4 ) - assert all(abs(subtract(reflection_spectrum, validated_reflection_spectrum)) < 1e-4) + assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-4) def test_fdtd_simulation_four_layers(): @@ -106,7 +101,7 @@ def test_fdtd_simulation_four_layers(): cfg.N_layers = 4 cfg.FDTD_niter = 400 cfg.N_freq = 5 - eps_grid = zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=cfloat) + eps_grid = torch.zeros(1, 1, cfg.N_layers, cfg.N_freq, dtype=torch.cfloat) cfg.target_frequencies = compute_target_frequencies( cfg.physical_wavelength_range[0], cfg.physical_wavelength_range[1], @@ -120,24 +115,22 @@ def test_fdtd_simulation_four_layers(): eps_grid[:, :, 2, :] = layer_builder.build_uniform_layer("gallium_arsenide") eps_grid[:, :, 3, :] = layer_builder.build_uniform_layer("silicon_nitride") transmission_spectrum, reflection_spectrum = compute_spectrum(eps_grid, cfg) - validated_transmission_spectrum = [ - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - tensor(0.0), - ] - validated_reflection_spectrum = [ - tensor(0.9515), - tensor(0.4005), - tensor(0.4919), - tensor(0.6812), - tensor(0.2888), - ] + validated_transmission_spectrum = torch.tensor( + [ + torch.tensor(0.0), + torch.tensor(0.0), + torch.tensor(0.0), + torch.tensor(0.0), + torch.tensor(0.0), + ] + ) + validated_reflection_spectrum = torch.tensor( + [0.95147296, 0.51595957, 0.49191941, 0.68122679, 0.28879608] + ) assert all( - abs(subtract(transmission_spectrum, validated_transmission_spectrum)) < 1e-4 + torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-4 ) - assert all(abs(subtract(reflection_spectrum, validated_reflection_spectrum)) < 1e-4) + assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-4) def test_single_patterned_layer(): diff --git a/nidn/training/utils/validate_config.py b/nidn/training/utils/validate_config.py index f6123db..732967d 100644 --- a/nidn/training/utils/validate_config.py +++ b/nidn/training/utils/validate_config.py @@ -1,4 +1,5 @@ from dotmap import DotMap +from ...utils.global_constants import UNIT_MAGNITUDE def _validate_config(cfg: DotMap): @@ -40,13 +41,14 @@ def _validate_config(cfg: DotMap): "TRCWA_L_grid", "TRCWA_NG", "FDTD_grid", - "FDTD_use_pointsource", - "FDTD_use_pulsesource", + "FDTD_source_type", + "FDTD_pulse_type", "FDTD_pml_thickness", "FDTD_source_position", "FDTD_free_space_distance", "FDTD_reflection_detector_x", "FDTD_niter", + "FDTD_min_gridpoints_between_detectors", "target_reflectance_spectrum", "target_transmittance_spectrum", "freq_distribution", @@ -70,6 +72,7 @@ def _validate_config(cfg: DotMap): "seed", "TRCWA_NG", "FDTD_niter", + "FDTD_min_gridpoints_between_detectors", ] float_keys = [ "L", @@ -89,10 +92,16 @@ def _validate_config(cfg: DotMap): "add_noise", "use_gpu", "avoid_zero_eps", - "FDTD_use_pulsesource", - "FDTD_use_pointsource", ] - string_keys = ["model_type", "type", "name", "freq_distribution", "solver"] + string_keys = [ + "model_type", + "type", + "name", + "freq_distribution", + "solver", + "FDTD_source_type", + "FDTD_pulse_type", + ] list_keys = [ "PER_LAYER_THICKNESS", "TRCWA_L_grid", @@ -158,6 +167,7 @@ def _validate_config(cfg: DotMap): "FDTD_free_space_distance", "FDTD_reflection_detector_x", "FDTD_pml_thickness", + "FDTD_min_gridpoints_between_detectors", ] for key in positive_value_keys: if not (cfg[key] > 0): @@ -212,9 +222,15 @@ def _validate_config(cfg: DotMap): raise ValueError(f"freq_distribution must be either 'linear' or 'log'") if not len(cfg.FDTD_grid) == 3: - raise ValueError(f"FDTD_grid must me 3-dimentional") + raise ValueError(f"FDTD_grid must be 3-dimentional") if not (len(cfg.FDTD_source_position) == 2 or len(cfg.FDTD_source_position) == 3): raise ValueError( f"The FDTD source needs either 2- or 3-dimensional coordinates" ) + if cfg.solver == "FDTD" and cfg.physical_wavelenght_range[0] >= ( + cfg.FDTD_grid[1] * UNIT_MAGNITUDE / 2 + ): + raise ValueError( + f" When using the FDTD solver, the smallest wavelength should be smaller than half the grid width in order for the scaling to work" + ) diff --git a/nidn/utils/resources/default_config.toml b/nidn/utils/resources/default_config.toml index a7d3e9f..407422f 100644 --- a/nidn/utils/resources/default_config.toml +++ b/nidn/utils/resources/default_config.toml @@ -44,14 +44,14 @@ TRCWA_NG = 11 # Truncation order (actual number might be smaller) # FDTD parameters FDTD_grid = [5.0,2.0,1.0] # Dimensions of the grid used in FDTD simulations, in FDTD unit magnitudes -FDTD_use_pointsource = false # If the source should be a point source in the FDTD simulaitons, otherwise set as linesource -FDTD_use_pulsesource = false # If the source should be a single pulse in the FDTD simulations, otherwise set as continuous wave +FDTD_source_type = "line" # Geometry of source, either "point" for pointsource or "line" for linesource +FDTD_pulse_type = "continuous" # If the source should be a single pulse in the FDTD simulations, or continuous wave. "pulse" or "continuous" accepted FDTD_pml_thickness = 1.5 # Thickness of PML layerin FDTD simulations, set in FDTD unit magnitudes. Perfectly Matched Layer are boundaries used in the x-direction, design to absorb all radiation FDTD_source_position = [1.5,1.0] # Coordinates of the source used in FDTD simulations, in FDTD unit magnitudes FDTD_free_space_distance = 1.0 # The thickness of the free space layer before and after the material layers, in FDTD unit magnitudes FDTD_reflection_detector_x = 2.4 # X-coordinates of the reflection detector for the FDTD simulation, in FDTD unit magnitudes FDTD_niter = 200 # Number of iterations / timesteps to run the FDTD for - +FDTD_min_gridpoints_between_detectors = 2 # Target spectra # To find out which frequencies below points correspond to use nidn.get_frequency_points(cfg) From c6b659a2cdf4e4a84a9dcaadaedbcc7031dd8521 Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Wed, 2 Mar 2022 11:44:18 +0100 Subject: [PATCH 7/8] Decreased error tolerance in test, according to the recently incresed precision --- nidn/tests/fdtd_test.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py index 587f202..7a6edf7 100644 --- a/nidn/tests/fdtd_test.py +++ b/nidn/tests/fdtd_test.py @@ -89,9 +89,9 @@ def test_fdtd_simulation_single_layer(): [0.95146948, 0.17358959, 0.35082144, 0.31708776, 0.34369091] ) assert all( - torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-4 + torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-8 ) - assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-4) + assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-8) def test_fdtd_simulation_four_layers(): @@ -128,9 +128,9 @@ def test_fdtd_simulation_four_layers(): [0.95147296, 0.51595957, 0.49191941, 0.68122679, 0.28879608] ) assert all( - torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-4 + torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-8 ) - assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-4) + assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-8) def test_single_patterned_layer(): From 53e851a449d4378c99a76dc912d91f137d036b3f Mon Sep 17 00:00:00 2001 From: torbjornstoro Date: Thu, 3 Mar 2022 09:59:58 +0100 Subject: [PATCH 8/8] Fixed final feedback --- ...te_transmission_reflection_coefficients.py | 14 +++++++------ .../fdtd_integration/compute_spectrum_fdtd.py | 2 +- nidn/fdtd_integration/init_fdtd.py | 21 ++++++++++++------- nidn/tests/fdtd_test.py | 8 +++++++ nidn/training/utils/validate_config.py | 5 +++++ nidn/utils/resources/default_config.toml | 2 +- 6 files changed, 36 insertions(+), 16 deletions(-) diff --git a/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py b/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py index b076f98..ed4fbde 100644 --- a/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py +++ b/nidn/fdtd_integration/calculate_transmission_reflection_coefficients.py @@ -52,14 +52,14 @@ def calculate_transmission_reflection_coefficients( reflection_signals[0], cfg ) - if not (transmission_coefficient >= 0 and transmission_coefficient <= 1): + if transmission_coefficient < 0 or transmission_coefficient > 1: raise ValueError( - f"The transmission coefficient is outsied of the possible range between 0 and 1" + f"The transmission coefficient is outside of the physical range between 0 and 1" ) - if not (reflection_coefficient >= 0 and reflection_coefficient <= 1): + if reflection_coefficient < 0 or reflection_coefficient > 1: raise ValueError( - f"The reflection coefficient is outsied of the possible range between 0 and 1" + f"The reflection coefficient is outside of the physical range between 0 and 1" ) if transmission_coefficient + reflection_coefficient > 1: logger.warning( @@ -90,9 +90,11 @@ def _fft(signal, cfg: DotMap): tuple[array,array]: fourier frequenices and their corresponding values """ sampling_frequencies = ( - cfg.physical_wavelength_range[0] * FDTD_GRID_SCALE / (sqrt(2) * SPEED_OF_LIGHT) + cfg.physical_wavelength_range[0] + * FDTD_GRID_SCALE + / (torch.sqrt(2) * SPEED_OF_LIGHT) ) - tensor_signal = tensor(signal) + tensor_signal = torch.tensor(signal) yf = rfft(tensor_signal) xf = rfftfreq(cfg.FDTD_niter, sampling_frequencies) diff --git a/nidn/fdtd_integration/compute_spectrum_fdtd.py b/nidn/fdtd_integration/compute_spectrum_fdtd.py index 618661a..56c6301 100644 --- a/nidn/fdtd_integration/compute_spectrum_fdtd.py +++ b/nidn/fdtd_integration/compute_spectrum_fdtd.py @@ -135,7 +135,7 @@ def _average_along_detector(signal): """ avg = [] for e in signal: - s = [torch.tensor(0.0), torch.tensor(0.0), torch.tensor(0.0)] + s = torch.zeros(3) for p in e: s[0] += p[0] / len(e) s[1] += p[1] / len(e) diff --git a/nidn/fdtd_integration/init_fdtd.py b/nidn/fdtd_integration/init_fdtd.py index 1a2cdec..b244e1a 100644 --- a/nidn/fdtd_integration/init_fdtd.py +++ b/nidn/fdtd_integration/init_fdtd.py @@ -1,6 +1,15 @@ from dotmap import DotMap -from ..fdtd import * +from ..fdtd import ( + AbsorbingObject, + set_backend, + Grid, + PML, + PeriodicBoundary, + LineSource, + LineDetector, + PointSource, +) from ..utils.global_constants import EPS_0, SPEED_OF_LIGHT, UNIT_MAGNITUDE from .constants import FDTD_GRID_SCALE @@ -42,13 +51,9 @@ def init_fdtd(cfg: DotMap, include_object, wavelength, permittivity): ), int(cfg.FDTD_reflection_detector_x * scaling), ) - use_pulse = True - if cfg.FDTD_pulse_type == "pulse": - pass - elif cfg.FDTD_pulse_type == "continuous": - use_pulse = False - else: - raise ValueError(f' FDTD_pulse_type must either be "pulse" or "continuous"') + assert cfg.FDTD_pulse_type in ["pulse", "continuous"] + use_pulse = cfg.FDTD_pulse_type == "pulse" + grid = _add_source( grid, int(cfg.FDTD_source_position[0] * scaling), diff --git a/nidn/tests/fdtd_test.py b/nidn/tests/fdtd_test.py index 7a6edf7..b78685c 100644 --- a/nidn/tests/fdtd_test.py +++ b/nidn/tests/fdtd_test.py @@ -92,6 +92,10 @@ def test_fdtd_simulation_single_layer(): torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-8 ) assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-8) + assert all(transmission_spectrum <= 1) + assert all(reflection_spectrum <= 1) + assert all(transmission_spectrum >= 0) + assert all(reflection_spectrum >= 0) def test_fdtd_simulation_four_layers(): @@ -131,6 +135,10 @@ def test_fdtd_simulation_four_layers(): torch.abs(transmission_spectrum - validated_transmission_spectrum) < 1e-8 ) assert all(torch.abs(reflection_spectrum - validated_reflection_spectrum) < 1e-8) + assert all(transmission_spectrum <= 1) + assert all(reflection_spectrum <= 1) + assert all(transmission_spectrum >= 0) + assert all(reflection_spectrum >= 0) def test_single_patterned_layer(): diff --git a/nidn/training/utils/validate_config.py b/nidn/training/utils/validate_config.py index 732967d..8b04faf 100644 --- a/nidn/training/utils/validate_config.py +++ b/nidn/training/utils/validate_config.py @@ -234,3 +234,8 @@ def _validate_config(cfg: DotMap): raise ValueError( f" When using the FDTD solver, the smallest wavelength should be smaller than half the grid width in order for the scaling to work" ) + if not cfg.FDTD_source_type in ["point", "line"]: + raise ValueError(f'The FDTD_source_type must either be "line" or "point"') + + if not cfg.FDTD_pulse_type in ["pulse", "continuous"]: + raise ValueError(f'The FDTD_pulse_type must either be "pulse" or "continuous"') diff --git a/nidn/utils/resources/default_config.toml b/nidn/utils/resources/default_config.toml index 407422f..1f86dfe 100644 --- a/nidn/utils/resources/default_config.toml +++ b/nidn/utils/resources/default_config.toml @@ -51,7 +51,7 @@ FDTD_source_position = [1.5,1.0] # Coordinates of the source used in FDTD simula FDTD_free_space_distance = 1.0 # The thickness of the free space layer before and after the material layers, in FDTD unit magnitudes FDTD_reflection_detector_x = 2.4 # X-coordinates of the reflection detector for the FDTD simulation, in FDTD unit magnitudes FDTD_niter = 200 # Number of iterations / timesteps to run the FDTD for -FDTD_min_gridpoints_between_detectors = 2 +FDTD_min_gridpoints_between_detectors = 2 # Minimal grid points between transmission detector and reflection detector in FDTD simulation. The true distance will be this value plus total layer thickness # Target spectra # To find out which frequencies below points correspond to use nidn.get_frequency_points(cfg)