# Mach-Zehnder lattice filter¶

This is an example of a Mach-Zehnder lattice filter with N stages.

import demolib.all as pdk
import ipkiss3.all as i3

import numpy as np
import pylab as plt


The non-dispersive symmetric directional coupler model:

$L_0 = \dfrac{\lambda_c arcsin(k_0)}{\pi \Delta n_{eff}}$

$\phi_0 = exp(i\dfrac{2\pi}{\lambda} n_{avg}L_{non-coupling})$

$\phi_{coupling} = exp(i\dfrac{2\pi}{\lambda} n_{avg} L_{coupling})$

$\Delta\beta L = \dfrac{2\pi}{\lambda}\Delta n_{eff} (L_{coupling}+ L_0)$

• Crossing transmission:

$T_{cross} = \dfrac{1}{2}(1 - exp(i\Delta\beta L))\phi_{coupling}\phi_0$

• Straight transmission:

$T_{straight} = \dfrac{1}{2}(1 + exp(i\Delta\beta L))\phi_{coupling}\phi_0$

class DirCoupModel(i3.CompactModel):
"""
Non-dispersive symmetric directional coupler model
"""
parameters = [
'delta_n',                 # Difference in effective index between the even and odd mode
'n_avg',                   # Average effective index
'coupler_length',          # Length of the coupling section
'non_coupling_length',     # Length of the directional coupler that is not part of the coupling
'coupling_at_zero_length', # Coupling for a zero length coupling, L_0 = zero length
'center_wavelength'        # Reference wavelength for all the parameters above
]

terms = [
i3.OpticalTerm(name='in1'),
i3.OpticalTerm(name='in2'),
i3.OpticalTerm(name='out1'),
i3.OpticalTerm(name='out2'),
]

def calculate_smatrix(parameters, env, S):
delta_n = parameters.delta_n
L_0 = parameters.center_wavelength * np.arcsin(parameters.coupling_at_zero_length ) / (np.pi * delta_n)
phi_0 = np.exp(1j * 2 * np.pi * parameters.n_avg * parameters.non_coupling_length / env.wavelength)
phi_coupler = np.exp(1j * 2 * np.pi * parameters.n_avg * parameters.coupler_length / env.wavelength)

delta_beta_x_length = 2 * np.pi / env.wavelength * delta_n * (parameters.coupler_length + L_0)

# Crossing transmission
cross = 0.5 * (1.0 - np.exp(1j * delta_beta_x_length)) * phi_coupler * phi_0
S['in1', 'out2'] = S['out2', 'in1'] = cross
S['in2', 'out1'] = S['out1', 'in2'] = cross

# Straight transmission
straight = 0.5 * (1.0 + np.exp(1j * delta_beta_x_length)) * phi_coupler * phi_0
S['in1', 'out1'] = S['out1', 'in1'] = straight
S['in2', 'out2'] = S['out2', 'in2'] = straight


The Directional Coupler described by the previous model:

class DirectionalCoupler(i3.PCell):
trace_template = i3.WaveguideTemplateProperty(doc="Trace template used for directional coupler")

def _default_trace_template(self):
return pdk.SiWireWaveguideTemplate()

class Layout(i3.LayoutView):
spacing = i3.PositiveNumberProperty(default=0.2, doc="Spacing between the waveguides", locked=True)
coupler_length = i3.PositiveNumberProperty(default=5.0, doc="Length of the straight section of the directional coupler")

def _generate_instances(self, insts):
delta = self.trace_template.core_width + self.spacing
coupler_length = self.coupler_length

shape = i3.Shape([
(-coupler_length / 2.0 - bend_radius, 0.0),
(-coupler_length / 2.0, 0.0),
(coupler_length / 2.0, 0.0),
(coupler_length / 2.0 + bend_radius, 0.0),
])

wg = i3.RoundedWaveguide(name=self.name + "_wg",
trace_template=self.trace_template)
wg_layout = wg.Layout(shape=shape)

insts += i3.SRef(reference=wg_layout,
name="wav_top",
position=(0, delta / 2.0))

insts += i3.SRef(reference=wg_layout,
name="wav_bot",
transformation=i3.VMirror() + i3.Translation(translation=(0, -delta / 2.0)))
return insts

def _generate_ports(self, ports):
ports += i3.expose_ports(self.instances, {
'wav_top:in': 'in2',
'wav_top:out': 'out2',
'wav_bot:in': 'in1',
'wav_bot:out': 'out1',
})
return ports

class CircuitModel(i3.CircuitModelView):
delta_n_eff = i3.NumberProperty(default=0.04, doc="Difference between even and odd supermode of the dc")
center_wavelength = i3.PositiveNumberProperty(default=1.31, doc="Reference wavelength for all parameters")
coupling_at_zero_length = i3.NonNegativeNumberProperty(default=0.03, doc="Field coupling for a zero length coupling")

def _generate_model(self):
trace_length = self.layout_view.instances['wav_top'].reference.trace_length()
coupler_length = self.layout_view.coupler_length
non_coupling_length = (trace_length - coupler_length)

return DirCoupModel(delta_n=self.delta_n_eff,
n_avg=self.trace_template.n_eff,
coupler_length=coupler_length,
non_coupling_length=non_coupling_length,
coupling_at_zero_length=self.coupling_at_zero_length,
center_wavelength=self.center_wavelength
)

class Netlist(i3.NetlistFromLayout):
pass


The default DirectionalCoupler:

dc = DirectionalCoupler()
dc_layout = dc.Layout()
print("coupler length: {0:.3f}".format(dc_layout.coupler_length))
dc_layout.visualize(annotate=True)


Out:

coupler length: 5.000


The MZI lattice filter that uses an undetermined number (n_stages + 1) of directional couplers.

class MZILattice(i3.PCell):
"""Mach-Zehnder based lattice filter with parameter n_stages.
"""
directional_couplers = i3.ChildCellListProperty(doc="list of directional couplers", locked=True)
delay_waveguides = i3.ChildCellListProperty(doc="List of waveguides between the directional couplers")
n_stages = i3.PositiveIntProperty(doc="Number of stages in the lattice.")
trace_template = i3.TraceTemplateProperty(doc="The trace template")

def _default_trace_template(self):
return pdk.SiWireWaveguideTemplate()

def _default_directional_couplers(self):
return [
DirectionalCoupler(name=self.name + "dc_{}".format(cnt),
trace_template=self.trace_template)
for cnt in range(self.n_stages + 1)
]

def _default_delay_waveguides(self):
north_wgs = [
i3.RoundedWaveguide(name=self.name + "wg_north_{}".format(ndx),
trace_template=self.trace_template)
for ndx in range(self.n_stages)]

south_wgs = [
i3.RoundedWaveguide(name=self.name + "wg_south_{}".format(ndx),
trace_template=self.trace_template)
for ndx in range(self.n_stages)]

return north_wgs + south_wgs

class Layout(i3.LayoutView):
delay_lengths = i3.ListProperty(default=[100.], doc="List of delay lengths")
dc_lengths = i3.ListProperty(default=[5., 5.], doc="List of directional coupler lengths")

def _default_directional_couplers(self):
dc_models = []
for idx, dc in enumerate(self.cell.directional_couplers):
dc_layout = dc.get_default_view(self)
dc_layout.set(coupler_length=self.dc_lengths[idx])
dc_models.append(dc_layout)
return dc_models

def _default_delay_waveguides(self):
north_delay_lengths = [max(dl, 0) for dl in self.delay_lengths]
south_delay_lengths = [abs(min(dl, 0)) for dl in self.delay_lengths]  # negative delay_length -> delay on south side

north_shapes = [i3.Shape([(0.0, 0.0), (0.0, min_dl + dl / 2.0), (2 * self.bend_radius, min_dl + dl / 2.0), (2 * self.bend_radius, 0.0)])
for dl in north_delay_lengths
]
south_shapes = [i3.Shape([(0.0, 0.0), (0.0, -min_dl - dl / 2.0), (2 * self.bend_radius, -min_dl - dl / 2.0), (2 * self.bend_radius, 0.0)])
for dl in south_delay_lengths
]

wgs = []
for wg, sh in zip(self.cell.delay_waveguides[:self.n_stages], north_shapes):
wg_lo = wg.get_default_view(i3.LayoutView)
wgs.append(wg_lo)

for wg, sh in zip(self.cell.delay_waveguides[self.n_stages:], south_shapes):
wg_lo = wg.get_default_view(i3.LayoutView)
wgs.append(wg_lo)

return wgs

def _generate_instances(self, insts):
"""Here we join the directional couplers and waveguides with the given delay length and
place the 'in1' port of the first directional coupler on (0.0, 0.0)."""
insts0 = {
"dc_{}".format(cnt): dc
for cnt, dc in enumerate(self.directional_couplers)
}

insts0.update({
'wg_north_{}'.format(cnt): wg
for cnt, wg in enumerate(self.delay_waveguides[:self.n_stages])
})

insts0.update({
'wg_south_{}'.format(cnt): wg
for cnt, wg in enumerate(self.delay_waveguides[self.n_stages:])
})

# Determine which ports of the directional couplers and waveguides should be joined:
south_in_joins = [("dc_{}:out1".format(cnt), "wg_south_{}:in".format(cnt))
for cnt in range(self.n_stages)]
south_out_joins = [("wg_south_{}:out".format(cnt), "dc_{}:in1".format(cnt + 1))
for cnt in range(self.n_stages)]
north_in_joins = [("dc_{}:out2".format(cnt), "wg_north_{}:in".format(cnt))
for cnt in range(self.n_stages)]
north_out_joins = [("wg_north_{}:out".format(cnt), "dc_{}:in2".format(cnt + 1))
for cnt in range(self.n_stages)]

joins = south_in_joins + south_out_joins + north_in_joins + north_out_joins

insts += i3.place_insts(
insts=insts0,
specs=[
i3.Place("dc_0:in1", (0.0, 0.0)),
i3.Join(joins)
]
)

return insts

def _generate_ports(self, ports):
ports = i3.expose_ports(
self.instances,
{"dc_0:in1": "in1",
"dc_0:in2": "in2",
"dc_{}:out1".format(self.n_stages): "out1",
"dc_{}:out2".format(self.n_stages): "out2",
}
)
return ports

class Netlist(i3.NetlistFromLayout):
pass

class CircuitModel(i3.CircuitModelView):
delta_n_eff = i3.NumberProperty(default=0.04, doc="Difference between even and odd mode of the dc")
center_wavelength = i3.PositiveNumberProperty(default=1.31, doc="Reference wavelength for all parameters")
coupling_at_zero_length = i3.NonNegativeNumberProperty(default=0.03, doc="Field coupling for a zero length coupling")

def _default_directional_couplers(self):
dc_models = []
for dc in self.cell.directional_couplers:
dc_cm = dc.get_default_view(self)
dc_cm.set(delta_n_eff=self.delta_n_eff,
center_wavelength=self.center_wavelength,
coupling_at_zero_length=self.coupling_at_zero_length)
dc_models.append(dc_cm)
return dc_models

def _generate_model(self):
return i3.HierarchicalModel.from_netlistview(self.netlist_view)


The model parameters:

tt = pdk.SiWireWaveguideTemplate()
tt.Layout().cross_section().visualize()

tt_cm = tt.CircuitModel(n_g=4.18920116149,
center_wavelength=1.31,
n_eff=2.4579358757)

delta_n_eff = 0.04
coupling_at_zero_length = 0.03


The filter coefficients below lead to a flat-band response.

fsr = 0.04
power_couplings = [0.5, 0.13, 0.12, 0.5, 0.25]
n_eff = tt_cm.n_eff
n_g = tt_cm.n_g
center_wavelength = tt_cm.center_wavelength

L = center_wavelength ** 2 / n_g / fsr
L_pi = center_wavelength / (2 * n_eff)

print("n_eff: {0:.2f}".format(n_eff))
print("n_g:   {0:.2f}".format(n_g))
print("center wavelength: {}".format(center_wavelength))
print("L:    {0:.2f}".format(L))
print("L_pi: {0:.2f}".format(L_pi))


Out:

n_eff: 2.46
n_g:   4.19
center wavelength: 1.31
L:    10.24
L_pi: 0.27


Calculate all layout parameters:

delay_lengths = [L, 2 * L, -(2 * L + L_pi), -2 * L]

def calculate_length_from_coupling(coupling):
"""Calculate the directional coupler length, required to have a certain amplitude coupling."""
L = center_wavelength * np.arcsin(coupling) / (np.pi * delta_n_eff)
L0 = center_wavelength * np.arcsin(coupling_at_zero_length) / (np.pi * delta_n_eff)
return L - L0

dc_lengths = [calculate_length_from_coupling(np.sqrt(kappa)) for kappa in power_couplings]

print("directional coupler lengths = {}".format(dc_lengths))
print("delay lengths = {}".format(delay_lengths))


Out:

directional coupler lengths = [7.8747136068958428, 3.5324806626416718, 3.374845604522557, 7.8747136068958428, 5.1455469402291767]
delay lengths = [10.241212667080566, 20.48242533416113, -20.748909096604063, -20.48242533416113]


Now we can make an MZI lattice with 4 stages, which uses the previously determined trace template. We can visualize it, and write it to a GDSII file.

cell = MZILattice(name="MziLattice", trace_template=tt, n_stages=4)

lv = cell.Layout(delay_lengths=delay_lengths, dc_lengths=dc_lengths)
lv.visualize(annotate=True)
lv.write_gdsii("mzi_lattice.gds")


We can get the S matrix of the circuit model and plot the transmission.

cm = cell.CircuitModel(delta_n_eff=delta_n_eff, center_wavelength=center_wavelength)
wavelengths = np.linspace(1.28, 1.35, 1001)

S = cm.get_smatrix(wavelengths=wavelengths)

plt.figure()
plt.plot(wavelengths, 20 * np.log10(np.abs(S["in1","out1"])), label='Transmission out 1')
plt.plot(wavelengths, 20 * np.log10(np.abs(S["in1","out2"])), label='Transmission out 2')
plt.ylabel("Transmission (dB)")
plt.xlabel("Wavelength ($\mu m$)")
plt.legend()
plt.ylim([-40, 0])
plt.show()


We can now check how sensitive our simulation is with respect to $\Delta n_{eff}$ . To do this we simulate the same circuit for different values of $\Delta n_{eff}$ and look at the responses.

plt.figure()

for delta_n_eff in np.linspace(0.03, 0.05, 5):
cm = cell.CircuitModel(delta_n_eff=delta_n_eff, center_wavelength=center_wavelength)
wavelengths = np.linspace(1.28, 1.35, 1001)

S = cm.get_smatrix(wavelengths=wavelengths)

plt.plot(wavelengths, 20 * np.log10(np.abs(S["in1","out1"])), label='$\Delta n_{eff}$ = '+'{0:.3f}'.format(delta_n_eff)+', $T_{out1}$')
plt.plot(wavelengths, 20 * np.log10(np.abs(S["in1","out2"])), label='$\Delta n_{eff}$ = '+'{0:.3f}'.format(delta_n_eff)+', $T_{out2}$')
plt.ylabel("Transmission (dB)")
plt.xlabel("Wavelength ($\mu m$)")
plt.legend()
plt.ylim([-40, 0])

plt.show()