-
Notifications
You must be signed in to change notification settings - Fork 21
/
Copy pathmetagrating_meep.py
141 lines (110 loc) · 4.96 KB
/
metagrating_meep.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
# Computes the transmittance of the m=+1 order
# for a 3D metagrating with the 2D design
# imported from the file device1.csv
import meep as mp
import math
import numpy as np
def metagrating(P_pol: bool):
resolution = 100 # pixels/μm
nSi = 3.45
Si = mp.Medium(index=nSi)
nSiO2 = 1.45
SiO2 = mp.Medium(index=nSiO2)
theta_d = math.radians(50.0) # deflection angle
wvl = 1.05 # wavelength
fcen = 1/wvl
px = wvl/math.sin(theta_d) # period in x
py = 0.5*wvl # period in y
dpml = 1.0 # PML thickness
gh = 0.325 # grating height
dsub = 5.0 # substrate thickness
dair = 5.0 # air padding
sz = dpml+dsub+gh+dair+dpml
cell_size = mp.Vector3(px,py,sz)
boundary_layers = [mp.PML(thickness=dpml,direction=mp.Z)]
# periodic boundary conditions
k_point = mp.Vector3()
# plane of incidence is XZ
src_cmpt = mp.Ex if P_pol else mp.Ey
src_pt = mp.Vector3(0,0,-0.5*sz+dpml+0.5*dsub)
sources = [mp.Source(src=mp.GaussianSource(fcen,fwidth=0.1*fcen),
size=mp.Vector3(px,py,0),
center=src_pt,
component=src_cmpt)]
sim = mp.Simulation(resolution=resolution,
cell_size=cell_size,
sources=sources,
default_material=SiO2,
boundary_layers=boundary_layers,
k_point=k_point)
flux = sim.add_mode_monitor(fcen,
0,
1,
mp.ModeRegion(center=mp.Vector3(0,0,0.5*sz-dpml),
size=mp.Vector3(px,py,0)))
stop_cond = mp.stop_when_fields_decayed(10,src_cmpt,src_pt,1e-7)
sim.run(until_after_sources=stop_cond)
input_flux = mp.get_fluxes(flux)
sim.reset_meep()
# image resolution is ~340 pixels/μm and thus
# Meep resolution should not be larger than ~half this value
weights = np.genfromtxt('device1.csv',delimiter=',')
Nx, Ny = weights.shape
geometry = [mp.Block(size=mp.Vector3(mp.inf,mp.inf,dpml+dsub),
center=mp.Vector3(0,0,-0.5*sz+0.5*(dpml+dsub)),
material=SiO2),
mp.Block(size=mp.Vector3(px,py,gh),
center=mp.Vector3(0,0,-0.5*sz+dpml+dsub+0.5*gh),
material=mp.MaterialGrid(grid_size=mp.Vector3(Nx,Ny,1),
medium1=mp.air,
medium2=Si,
weights=weights))]
sim = mp.Simulation(resolution=resolution,
cell_size=cell_size,
sources=sources,
geometry=geometry,
boundary_layers=boundary_layers,
k_point=k_point,
eps_averaging=False)
flux = sim.add_mode_monitor(fcen,
0,
1,
mp.ModeRegion(center=mp.Vector3(0,0,0.5*sz-dpml),
size=mp.Vector3(px,py,0)))
sim.run(until_after_sources=stop_cond)
res = sim.get_eigenmode_coefficients(flux,
mp.DiffractedPlanewave((1,0,0),
mp.Vector3(1,0,0),
0 if P_pol else 1,
1 if P_pol else 0))
coeffs = res.alpha
tran = abs(coeffs[0,0,0])**2 / input_flux[0]
print("tran:, {}, {:.6f}".format('P' if P_pol else 'S',tran))
# for debugging:
# visualize three orthogonal cross sections of the 3D cell
# to ensure that structure matches expected design
if 0:
import matplotlib
matplotlib.use('agg')
import matplotlib.pyplot as plt
output_plane = mp.Volume(center=mp.Vector3(0,0,0.5*sz-dpml-dair-0.5*gh),
size=mp.Vector3(px,py,0))
plt.figure()
sim.plot2D(output_plane=output_plane,
eps_parameters={'resolution':100})
plt.savefig('cell_xy.png',dpi=150,bbox_inches='tight')
output_plane = mp.Volume(center=mp.Vector3(0,0,0),
size=mp.Vector3(0,py,sz))
plt.figure()
sim.plot2D(output_plane=output_plane,
eps_parameters={'resolution':100})
plt.savefig('cell_yz.png',dpi=150,bbox_inches='tight')
output_plane = mp.Volume(center=mp.Vector3(0,0,0),
size=mp.Vector3(px,0,sz))
plt.figure()
sim.plot2D(output_plane=output_plane,
eps_parameters={'resolution':100})
plt.savefig('cell_xz.png',dpi=150,bbox_inches='tight')
if __name__ == '__main__':
metagrating(False)
metagrating(True)