基于Meep的伴随光子设计方法 的FDTD实现_python_代码_下载

import bayer_filter import math import meep as mp import interpolate as interp import time import numpy as np import sys # This should be pulled out to a different file that can be used for normalizations # So should the other type of normalization def normalize_by_frequency_gaussian_pulse(freq, fwidth, center_freq, t_offset): width = 1.0 / fwidth omega = 2 * math.pi * freq center_omega = 2 * math.pi * center_freq return np.abs(width) * (-omega / center_omega) * np.exp(-0.5 * width * width * (omega - center_omega)**2) * np.exp(np.complex(0, 1) * omega * t_offset) # return np.abs(width) * (-omega / center_omega) * np.exp(-0.5 * width * width * (omega - center_omega)**2) * np.exp(-np.complex(0, 1) * omega * t_offset) class BayerOptimization(): # Let's set up constants so that we can keep track of these quadrants and ranges in a more readable way RANGE_RED = 0 RANGE_GREEN = 1 RANGE_BLUE = 2 QUADRANT_RED = 2 QUADRANT_GREEN_XPOL = 1 QUADRANT_GREEN_YPOL = 3 QUADRANT_BLUE = 0 QUADRANT_TOP_RIGHT = 0 QUADRANT_TOP_LEFT = 1 QUADRANT_BOTTOM_LEFT = 2 QUADRANT_BOTTOM_RIGHT = 3 ################################### # # # # # # # P1 (G, X) # P0 (B) # # # # # # # ################################### # # # # # # # P2 (R) # P3 (G, Y) # # # # # # # ################################### def __init__(self, resolution): self.block_index = 1.53 self.block_permittivity = self.block_index**2 self.permittivity_bounds = [1.0, self.block_permittivity] self.length_scale_meters = 1e-6 self.length_scale_microns = self.length_scale_meters * 1e6 self.num_wavelengths_per_range = 1 self.num_ranges = 3 self.num_total_wavelengths = self.num_ranges * self.num_wavelengths_per_range self.num_total_frequencies = self.num_total_wavelengths self.blue_lower_lambda_um = .4 self.blue_upper_lambda_um = .5 self.green_lower_lambda_um = self.blue_upper_lambda_um self.green_upper_lambda_um = .6 self.red_lower_lambda_um = self.green_upper_lambda_um self.red_upper_lambda_um = .7 self.source_lower_lambda_um = self.blue_lower_lambda_um self.source_upper_lambda_um = self.red_upper_lambda_um self.source_upper_freq_unitless = self.length_scale_microns / self.source_lower_lambda_um; self.source_lower_freq_unitless = self.length_scale_microns / self.source_upper_lambda_um; self.source_freq_spread_unitless = self.source_upper_freq_unitless - self.source_lower_freq_unitless; self.source_freq_center_unitless = 0.5 * (self.source_lower_freq_unitless + self.source_upper_freq_unitless); self.dft_freq_min = self.source_lower_freq_unitless self.dft_freq_max = self.source_upper_freq_unitless self.dft_freq_spread = self.dft_freq_max - self.dft_freq_min self.dft_freq_center = self.dft_freq_min + 0.5 * self.dft_freq_spread self.frequencies_in_ranges = [[], [], []] self.max_intensity = [0] * self.num_total_frequencies for freq_idx in range(0, self.num_total_frequencies): frequency = ((float(freq_idx) / float(self.num_total_frequencies - 1)) * self.dft_freq_spread) + self.dft_freq_min wavelength_um = self.length_scale_microns / frequency self.max_intensity[freq_idx] = 1.0 / (wavelength_um * wavelength_um) self.frequencies_in_ranges[self.check_frequency_range(frequency)].append(frequency) self.max_frequencies_in_a_range = max(list(map(lambda y: len(y), self.frequencies_in_ranges))) # Maybe... set this by a real length scaled by the length scale? # self.spacing = 0.025 / self.length_scale_microns # self.resolution = int(1 / self.spacing) self.resolution = resolution self.spacing = 1 / resolution self.focal_length = 1.5 / self.length_scale_microns self.side_gap = 0.5 / self.length_scale_microns self.top_bottom_gap = .2 / self.length_scale_microns self.block_dim_x = 2 / self.length_scale_microns; self.block_dim_y = 2 / self.length_scale_microns; self.block_dim_z = 2 / self.length_scale_microns; self.block_size = mp.Vector3(self.block_dim_x, self.block_dim_y, self.block_dim_z) self.source_loc_z_from_pml = 0.0 / self.length_scale_microns; self.t_pml_x = 1.0 / self.length_scale_microns self.t_pml_y = 1.0 / self.length_scale_microns self.t_pml_z = 1.0 / self.length_scale_microns self.cell_dim_x = 2 * self.t_pml_x + 2 * self.side_gap + self.block_dim_x; self.cell_dim_y = 2 * self.t_pml_y + 2 * self.side_gap + self.block_dim_y; self.cell_dim_z = 2 * self.t_pml_z + 2 * self.top_bottom_gap + self.block_dim_z + self.focal_length; self.cell = mp.Vector3(self.cell_dim_x, self.cell_dim_y, self.cell_dim_z) self.block_center = mp.Vector3(0, 0, 0.5 * self.cell_dim_z - self.t_pml_z - self.top_bottom_gap - 0.5 * self.block_dim_z) self.top_source_center = mp.Vector3(0, 0, 0.5 * self.cell_dim_z - self.t_pml_z) self.device = bayer_filter.BayerFilter(self.block_center, self.resolution, self.block_size, self.permittivity_bounds, self.length_scale_meters); self.pml_layers = [mp.PML(thickness = self.t_pml_x, direction = mp.X), mp.PML(thickness = self.t_pml_y, direction = mp.Y), mp.PML(thickness = self.t_pml_z, direction = mp.Z)] # Run all the way up to and through the PML boundary layers self.top_bottom_source_size_x = self.cell_dim_x self.top_bottom_source_size_y = self.cell_dim_y # Surface currenets self.top_bottom_source_size_z = 0; self.top_bottom_source_size = mp.Vector3(self.top_bottom_source_size_x, self.top_bottom_source_size_y, self.top_bottom_source_size_z) # This is the default value for it, but let's be explicit self.gaussian_cutoff = 5.0 self.gaussian_width = 1.0 / self.source_freq_spread_unitless self.gaussian_peak_time = self.gaussian_cutoff * self.gaussian_width top_source_J = mp.Source(mp.GaussianSource(self.source_freq_center_unitless, fwidth=self.source_freq_spread_unitless, cutoff=self.gaussian_cutoff, start_time=0), component=mp.Ex, amplitude=1.0, center=self.top_source_center, size=self.top_bottom_source_size) self.forward_sources = [top_source_J] self.geometry = [mp.Block(size=self.block_size, center=self.block_center, epsilon_func=lambda loc: self.device.get_permittivity_value(loc))] self.z_measurement_point = self.block_center.z - 0.5 * self.block_dim_z - self.focal_length # Let's set up where we want to measure the figure of merit self.measure_blue_focal = mp.Vector3(self.block_dim_x / 4.0, self.block_dim_y / 4.0, self.z_measurement_point) self.measure_green_xpol_focal = mp.Vector3(-self.block_dim_x / 4.0, self.block_dim_y / 4.0, self.z_measurement_point) self.measure_red_focal = mp.Vector3(-self.block_dim_x / 4.0, -self.block_dim_y / 4.0, self.z_measurement_point) self.measure_green_ypol_focal = mp.Vector3(self.block_dim_x / 4.0, -self.block_dim_y / 4.0, self.z_measurement_point) self.dft_focal_plane_fields_center = mp.Vector3(0, 0, self.z_measurement_point) self.dft_focal_area_fields_size = mp.Vector3(self.block_dim_x, self.block_dim_y, 0); self.interpolate_focal_area_fields = interp.Interpolate(self.dft_focal_plane_fields_center, self.dft_focal_area_fields_size, self.resolution, 2) self.dft_block_fields_center = self.block_center self.dft_block_fields_size = self.block_size self.focal_plane_centers = [self.measure_blue_focal, self.measure_green_xpol_focal, self.measure_red_focal, self.measure_green_ypol_focal] self.xy_symmetry_transformation = [0] * 4 self.xy_symmetry_transformation[BayerOptimization.QUADRANT_TOP_RIGHT] = BayerOptimization.QUADRANT_TOP_RIGHT self.xy_symmetry_transformation[BayerOptimization.QUADRANT_TOP_LEFT] = BayerOptimization.QUADRANT_BOTTOM_RIGHT self.xy_symmetry_transformation[BayerOptimization.QUADRANT_BOTTOM_LEFT] = BayerOptimization.QUADRANT_BOTTOM_LEFT self.xy_symmetry_transformation[Baye