Ybranch
Ybranch
Y branch
splitter and combiner
Source Code
1
2
3
4
5
6
7
# importing necessary packages
import gplugins.modes as gmode
import numpy as np
import matplotlib.pyplot as plt
import meep as mp
import gdsfactory as gf
import gplugins.gmeep as gm
1
2
Using MPI version 4.1, 1 processes
[32m2025-10-02 17:22:26.906[0m | [1mINFO [0m | [36mgplugins.gmeep[0m:[36m<module>[0m:[36m39[0m - [1mMeep '1.31.0' installed at ['/home/ramprakash/anaconda3/envs/si_photo/lib/python3.13/site-packages/meep'][0m
Y-branch from gdsfactory
Source Code
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
gf.clear_cache()
ybranch = gf.Component("ybranch")
inwg = gf.components.straight(length=5)
outwg = gf.components.straight(length=2)
splitter = gf.components.mmi1x2_with_sbend()
inwg_ref = ybranch.add_ref(inwg)
splitter_ref = ybranch.add_ref(splitter)
outwg_top_ref = ybranch.add_ref(outwg)
outwg_bot_ref = ybranch.add_ref(outwg)
# Connecting ports
inwg_ref.connect("o2", splitter_ref.ports["o1"])
outwg_top_ref.connect("o1", splitter_ref.ports["o2"])
outwg_bot_ref.connect("o1", splitter_ref.ports["o3"])
# adding port name
ybranch.add_port(name=f"o1", port=inwg_ref.ports["o1"])
ybranch.add_port(name=f"o2", port=outwg_top_ref.ports["o2"])
ybranch.add_port(name=f"o3", port=outwg_bot_ref.ports["o2"])
ybranch.auto_rename_ports()
# plotting
ybranch.draw_ports()
ybranch.plot()
Meep simulation
S-parameter calculation and steady-state fields
50-50 splitter
Input port 1, will not be exact 50-50 because of the losses at the split
Source Code
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
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
from gdsfactory.technology import LayerLevel, LayerStack
mp.verbosity(0)
# Set up frequency points for simulation
npoints = 100
lcen = 1.55
dlam = 0.02
wl = np.linspace(lcen - dlam / 2, lcen + dlam / 2, npoints)
fcen = 1 / lcen
fwidth = 3 * dlam / lcen**2
fpoints = 1 / wl # Center frequency
mode_parity = mp.EVEN_Y + mp.ODD_Z
dpml = 1
dpad = 1
resolution = 20
# Define materials
Si = mp.Medium(index=3.45)
SiO2 = mp.Medium(index=1.45)
cell_size = mp.Vector3(ybranch.xsize + 2 * dpml, ybranch.ysize + 2 * dpml + 2 * dpad, 0)
ybranch = gf.components.extend_ports(ybranch, port_names=["o1", "o2", "o3"], length=1)
ybranch = ybranch.copy()
ybranch.flatten()
ybranch.center = (0, 0)
geometry = gm.get_meep_geometry.get_meep_geometry_from_component(ybranch)
geometry = [
mp.Prism(geom.vertices, geom.height, geom.axis, geom.center, material=Si)
for geom in geometry
]
src = mp.GaussianSource(frequency=fcen, fwidth=fwidth)
source = [
mp.EigenModeSource(
src=src,
eig_band=1,
eig_parity=mode_parity,
size=mp.Vector3(0, 1),
center=mp.Vector3(ybranch.ports["o1"].x + dpml + 1, ybranch["o1"].y),
)
]
sim = mp.Simulation(
resolution=resolution,
cell_size=cell_size,
boundary_layers=[
mp.PML(dpml)
], # the boundary layers to absorb fields that leave the simulation
sources=source, # The sources
geometry=geometry, # The geometry
default_material=SiO2,
)
m1 = mp.Volume(
center=mp.Vector3(ybranch.ports["o1"].x + dpml + 1 + 0.5, ybranch["o1"].y),
size=mp.Vector3(0, 1),
)
m2 = mp.Volume(
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1 - 0.5, ybranch["o2"].y),
size=mp.Vector3(0, 1),
)
m3 = mp.Volume(
center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1 - 0.5, ybranch["o3"].y),
size=mp.Vector3(0, 1),
)
mode_monitor_1 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m1))
mode_monitor_2 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m2))
mode_monitor_3 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m3))
whole_dft = sim.add_dft_fields([mp.Ez], fcen, 0, 1, center=mp.Vector3(), size=cell_size)
sim.plot2D(labels=False)
sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
# Finds the S parameters
norm_mode_coeff = sim.get_eigenmode_coefficients(
mode_monitor_1, [1], eig_parity=mode_parity
).alpha[0, :, 0]
port1_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_1, [1], eig_parity=mode_parity).alpha[
0, :, 1
]
/ norm_mode_coeff
)
port2_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_2, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
port3_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_3, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
# Calculates the transmittance based off of the S parameters
port1_trans = abs(port1_coeff) ** 2
port2_trans = abs(port2_coeff) ** 2
port3_trans = abs(port3_coeff) ** 2
insertion_loss_dB = 10*np.log10(port2_trans)
fig = plt.figure(figsize=(15, 4))
ax_loss = fig.add_subplot(1,2,1)
ax_loss.set_title("Insertion Loss P2")
ax_loss.plot(wl,insertion_loss_dB)
ax_loss.set_xlabel(r"Wavelength (\mu m)")
ax_loss.set_ylabel("Loss (dB)")
# ax_trans1 = fig.add_subplot(1, 3, 1)
# ax_trans1.set_title("Transmission per Port")
# ax_trans1.plot(wl, port1_trans, label=r"s11^2")
# ax_trans1.plot(wl, port2_trans, label=r"s21^2")
# ax_trans1.plot(wl, port3_trans, label=r"s31^2")
# ax_trans1.set_xlabel(r"Wavelength (\mu m)")
# ax_trans1.set_ylabel(r"Transmission")
# ax_trans1.legend()
# ax_phase1 = fig.add_subplot(1, 3, 2)
# ax_phase1.set_title("Phase per Port")
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port1_coeff) * 180 / np.pi), label=r"phase of s11"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port2_coeff) * 180 / np.pi), label=r"phase of s21"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port3_coeff) * 180 / np.pi), label=r"phase of s31"
# )
# ax_phase1.set_xlabel(r"Wavelength (\mu m)")
# ax_phase1.set_ylabel("Angle (deg)")
# ax_phase1.legend()
# sim.reset_meep()
# sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
eps_data = sim.get_epsilon() #
ez_data = sim.get_dft_array(
whole_dft, mp.Ez, 0
) # Values for the component of the E-field in the z direction (in/out of screen)
# Creates the plot
ax_field = fig.add_subplot(1, 2, 2)
ax_field.set_title("Steady State Fields")
ax_field.imshow(np.transpose(eps_data), interpolation="spline36", cmap="binary")
ax_field.imshow(
np.flipud(np.transpose(np.real(ez_data))),
interpolation="spline36",
cmap="RdBu",
alpha=0.9,
)
ax_field.axis("off")
plt.show()
Input in single guide
Splitter
Source Code
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
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
from gdsfactory.technology import LayerLevel, LayerStack
mp.verbosity(0)
# Set up frequency points for simulation
npoints = 100
lcen = 1.55
dlam = 0.02
wl = np.linspace(lcen - dlam / 2, lcen + dlam / 2, npoints)
fcen = 1 / lcen
fwidth = 3 * dlam / lcen**2
fpoints = 1 / wl # Center frequency
mode_parity = mp.EVEN_Y + mp.ODD_Z
dpml = 1
dpad = 1
resolution = 20
# Define materials
Si = mp.Medium(index=3.45)
SiO2 = mp.Medium(index=1.45)
cell_size = mp.Vector3(ybranch.xsize + 2 * dpml, ybranch.ysize + 2 * dpml + 2 * dpad, 0)
ybranch = gf.components.extend_ports(ybranch, port_names=["o1", "o2", "o3"], length=1)
ybranch = ybranch.copy()
ybranch.flatten()
ybranch.center = (0, 0)
geometry = gm.get_meep_geometry.get_meep_geometry_from_component(ybranch)
geometry = [
mp.Prism(geom.vertices, geom.height, geom.axis, geom.center, material=Si)
for geom in geometry
]
src = mp.GaussianSource(frequency=fcen, fwidth=fwidth)
source = [
mp.EigenModeSource(
src=src,
eig_band=1,
eig_parity=mode_parity,
eig_kpoint = mp.Vector3(-1,0,0),
direction = mp.NO_DIRECTION,
size=mp.Vector3(0, 1),
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1, ybranch["o2"].y),
amplitude=1
),
# mp.EigenModeSource(
# src=src,
# eig_band=1,
# eig_parity=mode_parity,
# eig_kpoint = mp.Vector3(-1,0,0),
# direction = mp.NO_DIRECTION,
# size=mp.Vector3(0, 1),
# center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1, ybranch["o3"].y),
# amplitude=1
# ),
]
sim = mp.Simulation(
resolution=resolution,
cell_size=cell_size,
boundary_layers=[mp.PML(dpml)], #
sources=source, # The sources
geometry=geometry, # The geometry
default_material=SiO2,
)
m1 = mp.Volume(
center=mp.Vector3(ybranch.ports["o1"].x + dpml + 1 + 0.5, ybranch["o1"].y),
size=mp.Vector3(0, 1),
)
m2 = mp.Volume(
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1 - 0.5, ybranch["o2"].y),
size=mp.Vector3(0, 1),
)
m3 = mp.Volume(
center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1 - 0.5, ybranch["o3"].y),
size=mp.Vector3(0, 1),
)
mode_monitor_1 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m1))
mode_monitor_2 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m2))
mode_monitor_3 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m3))
whole_dft = sim.add_dft_fields([mp.Ez], fcen, 0, 1, center=mp.Vector3(), size=cell_size)
sim.plot2D(labels=False)
sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
# Finds the S parameters
norm_mode_coeff = sim.get_eigenmode_coefficients(
mode_monitor_2, [1], eig_parity=mode_parity
).alpha[0, :, 1] # Adding both the ports input values
port1_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_1, [1], eig_parity=mode_parity).alpha[
0, :, 1
]
/ norm_mode_coeff
)
port2_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_2, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
port3_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_3, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
# Calculates the transmittance based off of the S parameters
port1_trans = abs(port1_coeff) ** 2
port2_trans = abs(port2_coeff) ** 2
port3_trans = abs(port3_coeff) ** 2
insertion_loss_dB = 10*np.log10(port1_trans)
fig = plt.figure(figsize=(15, 4))
ax_loss = fig.add_subplot(1,2,1)
ax_loss.set_title("Insertion Loss P1")
ax_loss.plot(wl,insertion_loss_dB)
ax_loss.set_xlabel(r"Wavelength (\mu m)")
ax_loss.set_ylabel("Loss (dB)")
# ax_trans1 = fig.add_subplot(1, 3, 1)
# ax_trans1.set_title("Transmission per Port")
# ax_trans1.plot(wl, port1_trans, label=r"s1(23)^2") # Not the proper way to represent
# ax_trans1.plot(wl, port2_trans, label=r"s22^2")
# ax_trans1.plot(wl, port3_trans, label=r"s33^2")
# ax_trans1.set_xlabel(r"Wavelength (\mu m)")
# ax_trans1.set_ylabel(r"Transmission")
# ax_trans1.legend()
# ax_phase1 = fig.add_subplot(1, 3, 2)
# ax_phase1.set_title("Phase per Port")
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port1_coeff) * 180 / np.pi), label=r"phase of s1(23)"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port2_coeff) * 180 / np.pi), label=r"phase of s22"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port3_coeff) * 180 / np.pi), label=r"phase of s33"
# )
# ax_phase1.set_xlabel(r"Wavelength (\mu m)")
# ax_phase1.set_ylabel("Angle (deg)")
# ax_phase1.legend()
# sim.reset_meep()
# sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
eps_data = sim.get_epsilon() #
ez_data = sim.get_dft_array(
whole_dft, mp.Ez, 0
) # Values for the component of the E-field in the z direction (in/out of screen)
# Creates the plot
ax_field = fig.add_subplot(1, 2, 2)
ax_field.set_title("Steady State Fields")
ax_field.imshow(np.transpose(eps_data), interpolation="spline36", cmap="binary")
ax_field.imshow(
np.flipud(np.transpose(np.real(ez_data))),
interpolation="spline36",
cmap="RdBu",
alpha=0.9,
)
ax_field.axis("off")
plt.show()
Combiner
Out-phase
distructive interference
Source Code
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
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
from gdsfactory.technology import LayerLevel, LayerStack
mp.verbosity(0)
# Set up frequency points for simulation
npoints = 100
lcen = 1.55
dlam = 0.02
wl = np.linspace(lcen - dlam / 2, lcen + dlam / 2, npoints)
fcen = 1 / lcen
fwidth = 3 * dlam / lcen**2
fpoints = 1 / wl # Center frequency
mode_parity = mp.EVEN_Y + mp.ODD_Z
dpml = 1
dpad = 1
resolution = 20
# Define materials
Si = mp.Medium(index=3.45)
SiO2 = mp.Medium(index=1.45)
cell_size = mp.Vector3(ybranch.xsize + 2 * dpml, ybranch.ysize + 2 * dpml + 2 * dpad, 0)
ybranch = gf.components.extend_ports(ybranch, port_names=["o1", "o2", "o3"], length=1)
ybranch = ybranch.copy()
ybranch.flatten()
ybranch.center = (0, 0)
geometry = gm.get_meep_geometry.get_meep_geometry_from_component(ybranch)
geometry = [
mp.Prism(geom.vertices, geom.height, geom.axis, geom.center, material=Si)
for geom in geometry
]
src = mp.GaussianSource(frequency=fcen, fwidth=fwidth)
source = [
mp.EigenModeSource(
src=src,
eig_band=1,
eig_parity=mode_parity,
eig_kpoint = mp.Vector3(-1,0,0),
direction = mp.NO_DIRECTION,
size=mp.Vector3(0, 1),
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1, ybranch["o2"].y),
amplitude=1
),
mp.EigenModeSource(
src=src,
eig_band=1,
eig_parity=mode_parity,
eig_kpoint = mp.Vector3(-1,0,0),
direction = mp.NO_DIRECTION,
size=mp.Vector3(0, 1),
center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1, ybranch["o3"].y),
amplitude=-1
),
]
sim = mp.Simulation(
resolution=resolution,
cell_size=cell_size,
boundary_layers=[mp.PML(dpml)], #
sources=source, # The sources
geometry=geometry, # The geometry
default_material=SiO2,
)
m1 = mp.Volume(
center=mp.Vector3(ybranch.ports["o1"].x + dpml + 1 + 0.5, ybranch["o1"].y),
size=mp.Vector3(0, 1),
)
m2 = mp.Volume(
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1 - 0.5, ybranch["o2"].y),
size=mp.Vector3(0, 1),
)
m3 = mp.Volume(
center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1 - 0.5, ybranch["o3"].y),
size=mp.Vector3(0, 1),
)
mode_monitor_1 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m1))
mode_monitor_2 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m2))
mode_monitor_3 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m3))
whole_dft = sim.add_dft_fields([mp.Ez], fcen, 0, 1, center=mp.Vector3(), size=cell_size)
sim.plot2D(labels=False)
sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
# Finds the S parameters
norm_mode_coeff = sim.get_eigenmode_coefficients(
mode_monitor_2, [1], eig_parity=mode_parity
).alpha[0, :, 1] + sim.get_eigenmode_coefficients(
mode_monitor_3, [1], eig_parity=mode_parity
).alpha[0, :, 1] # Adding both the ports input values
port1_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_1, [1], eig_parity=mode_parity).alpha[
0, :, 1
]
/ norm_mode_coeff
)
port2_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_2, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
port3_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_3, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
# Calculates the transmittance based off of the S parameters
port1_trans = abs(port1_coeff) ** 2
port2_trans = abs(port2_coeff) ** 2
port3_trans = abs(port3_coeff) ** 2
insertion_loss_dB = 10*np.log10(port1_trans)
fig = plt.figure(figsize=(15, 4))
ax_loss = fig.add_subplot(1,2,1)
ax_loss.set_title("Insertion Loss P1")
ax_loss.plot(wl,insertion_loss_dB)
ax_loss.set_xlabel(r"Wavelength (\mu m)")
ax_loss.set_ylabel("Loss (dB)")
# ax_trans1 = fig.add_subplot(1, 3, 1)
# ax_trans1.set_title("Transmission per Port")
# ax_trans1.plot(wl, port1_trans, label=r"s1(23)^2") # Not the proper way to represent
# ax_trans1.plot(wl, port2_trans, label=r"s22^2")
# ax_trans1.plot(wl, port3_trans, label=r"s33^2")
# ax_trans1.set_xlabel(r"Wavelength (\mu m)")
# ax_trans1.set_ylabel(r"Transmission")
# ax_trans1.legend()
# ax_phase1 = fig.add_subplot(1, 3, 2)
# ax_phase1.set_title("Phase per Port")
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port1_coeff) * 180 / np.pi), label=r"phase of s1(23)"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port2_coeff) * 180 / np.pi), label=r"phase of s22"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port3_coeff) * 180 / np.pi), label=r"phase of s33"
# )
# ax_phase1.set_xlabel(r"Wavelength (\mu m)")
# ax_phase1.set_ylabel("Angle (deg)")
# ax_phase1.legend()
# sim.reset_meep()
# sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
eps_data = sim.get_epsilon() #
ez_data = sim.get_dft_array(
whole_dft, mp.Ez, 0
) # Values for the component of the E-field in the z direction (in/out of screen)
# Creates the plot
ax_field = fig.add_subplot(1, 2, 2)
ax_field.set_title("Steady State Fields")
ax_field.imshow(np.transpose(eps_data), interpolation="spline36", cmap="binary")
ax_field.imshow(
np.flipud(np.transpose(np.real(ez_data))),
interpolation="spline36",
cmap="RdBu",
alpha=0.9,
)
ax_field.axis("off")
plt.show()
In-phase
constructive interference
Source Code
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
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
from gdsfactory.technology import LayerLevel, LayerStack
mp.verbosity(0)
# Set up frequency points for simulation
npoints = 100
lcen = 1.55
dlam = 0.02
wl = np.linspace(lcen - dlam / 2, lcen + dlam / 2, npoints)
fcen = 1 / lcen
fwidth = 3 * dlam / lcen**2
fpoints = 1 / wl # Center frequency
mode_parity = mp.EVEN_Y + mp.ODD_Z
dpml = 1
dpad = 1
resolution = 20
# Define materials
Si = mp.Medium(index=3.45)
SiO2 = mp.Medium(index=1.45)
cell_size = mp.Vector3(ybranch.xsize + 2 * dpml, ybranch.ysize + 2 * dpml + 2 * dpad, 0)
ybranch = gf.components.extend_ports(ybranch, port_names=["o1", "o2", "o3"], length=1)
ybranch = ybranch.copy()
ybranch.flatten()
ybranch.center = (0, 0)
geometry = gm.get_meep_geometry.get_meep_geometry_from_component(ybranch)
geometry = [
mp.Prism(geom.vertices, geom.height, geom.axis, geom.center, material=Si)
for geom in geometry
]
src = mp.GaussianSource(frequency=fcen, fwidth=fwidth)
source = [
mp.EigenModeSource(
src=src,
eig_band=1,
eig_parity=mode_parity,
eig_kpoint = mp.Vector3(-1,0,0),
direction = mp.NO_DIRECTION,
size=mp.Vector3(0, 1),
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1, ybranch["o2"].y),
amplitude=1
),
mp.EigenModeSource(
src=src,
eig_band=1,
eig_parity=mode_parity,
eig_kpoint = mp.Vector3(-1,0,0),
direction = mp.NO_DIRECTION,
size=mp.Vector3(0, 1),
center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1, ybranch["o3"].y),
amplitude=1
),
]
sim = mp.Simulation(
resolution=resolution,
cell_size=cell_size,
boundary_layers=[mp.PML(dpml)], #
sources=source, # The sources
geometry=geometry, # The geometry
default_material=SiO2,
)
m1 = mp.Volume(
center=mp.Vector3(ybranch.ports["o1"].x + dpml + 1 + 0.5, ybranch["o1"].y),
size=mp.Vector3(0, 1),
)
m2 = mp.Volume(
center=mp.Vector3(ybranch.ports["o2"].x - dpml - 1 - 0.5, ybranch["o2"].y),
size=mp.Vector3(0, 1),
)
m3 = mp.Volume(
center=mp.Vector3(ybranch.ports["o3"].x - dpml - 1 - 0.5, ybranch["o3"].y),
size=mp.Vector3(0, 1),
)
mode_monitor_1 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m1))
mode_monitor_2 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m2))
mode_monitor_3 = sim.add_mode_monitor(fpoints, mp.ModeRegion(volume=m3))
whole_dft = sim.add_dft_fields([mp.Ez], fcen, 0, 1, center=mp.Vector3(), size=cell_size)
sim.plot2D(labels=False)
sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
# Finds the S parameters
norm_mode_coeff = sim.get_eigenmode_coefficients(
mode_monitor_2, [1], eig_parity=mode_parity
).alpha[0, :, 1] + sim.get_eigenmode_coefficients(
mode_monitor_3, [1], eig_parity=mode_parity
).alpha[0, :, 1] # Adding both the ports input values
port1_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_1, [1], eig_parity=mode_parity).alpha[
0, :, 1
]
/ norm_mode_coeff
)
port2_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_2, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
port3_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_3, [1], eig_parity=mode_parity).alpha[
0, :, 0
]
/ norm_mode_coeff
)
# Calculates the transmittance based off of the S parameters
port1_trans = abs(port1_coeff) ** 2
port2_trans = abs(port2_coeff) ** 2
port3_trans = abs(port3_coeff) ** 2
insertion_loss_dB = 10*np.log10(port1_trans)
fig = plt.figure(figsize=(15, 4))
ax_loss = fig.add_subplot(1,2,1)
ax_loss.set_title("Insertion Loss P1")
ax_loss.plot(wl,insertion_loss_dB)
ax_loss.set_xlabel(r"Wavelength (\mu m)")
ax_loss.set_ylabel("Loss (dB)")
# ax_trans1 = fig.add_subplot(1, 3, 1)
# ax_trans1.set_title("Transmission per Port")
# ax_trans1.plot(wl, port1_trans, label=r"s1(23)^2") # Not the proper way to represent
# ax_trans1.plot(wl, port2_trans, label=r"s22^2")
# ax_trans1.plot(wl, port3_trans, label=r"s33^2")
# ax_trans1.set_xlabel(r"Wavelength (\mu m)")
# ax_trans1.set_ylabel(r"Transmission")
# ax_trans1.legend()
# ax_phase1 = fig.add_subplot(1, 3, 2)
# ax_phase1.set_title("Phase per Port")
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port1_coeff) * 180 / np.pi), label=r"phase of s1(23)"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port2_coeff) * 180 / np.pi), label=r"phase of s22"
# )
# ax_phase1.plot(
# wl, np.unwrap(np.angle(port3_coeff) * 180 / np.pi), label=r"phase of s33"
# )
# ax_phase1.set_xlabel(r"Wavelength (\mu m)")
# ax_phase1.set_ylabel("Angle (deg)")
# ax_phase1.legend()
# sim.reset_meep()
# sim.run(until_after_sources=mp.stop_when_dft_decayed(tol=1e-4))
eps_data = sim.get_epsilon() #
ez_data = sim.get_dft_array(
whole_dft, mp.Ez, 0
) # Values for the component of the E-field in the z direction (in/out of screen)
# Creates the plot
ax_field = fig.add_subplot(1, 2, 2)
ax_field.set_title("Steady State Fields")
ax_field.imshow(np.transpose(eps_data), interpolation="spline36", cmap="binary")
ax_field.imshow(
np.flipud(np.transpose(np.real(ez_data))),
interpolation="spline36",
cmap="RdBu",
alpha=0.9,
)
ax_field.axis("off")
plt.show()
Source Code
1
2
3
4
5
6
7
8
9
10
11
12
norm_mode_coeff = sim.get_eigenmode_coefficients(
mode_monitor_3, [1], eig_parity=mode_parity
).alpha[0, :, 1]
port1_coeff = (
sim.get_eigenmode_coefficients(mode_monitor_1, [1], eig_parity=mode_parity).alpha[
0, :, 1
]
/ norm_mode_coeff
)
port1_trans = abs(port1_coeff) ** 2
insertion_loss_dB = 10*np.log10(port1_trans)
plt.plot(wl, insertion_loss_dB)
1
[<matplotlib.lines.Line2D at 0x7b118fc26fd0>]
Source Code
1
plt.plot(wl, port1_trans)
1
[<matplotlib.lines.Line2D at 0x7b118fca1d10>]
This post is licensed under CC BY 4.0 by the author.











