Mach-Zehnder interferometer#

In this tutorial, we’ll define and simulate a simple circuit known as a Mach-Zender Interferometer (MZI). It assumes you’ve read the previous tutorial, Introduction to simphony.

Simphony uses SAX to define models and simulate circuits, which uses JAX as a computational engine. JAX can provide a nice speedup for larger circuits if their models are appropriately defined and you have a GPU. Otherwise, it will run perfectly fine on a CPU.

Note

We run the following command first to ensure that JAX, the library that allows Simphony and SAX to run calculations on GPU’s, uses double precision. Be aware that this setting must be set before JAX initializes, or the setting won’t take.

from jax import config
config.update("jax_enable_x64", True)

In this tutorial, we’ll use JAX as a drop in replacement for NumPy, matplotlib for visualizing our results, and SAX for constructing our circuit and running our simulations.

import jax.numpy as jnp
import matplotlib.pyplot as plt
import sax

The MZI from it’s components#

In an MZI, light entering the circuit is split and travels down two paths of differing lengths. When the light is recombined, it interferes, and the output is frequency-dependent.

Mach-Zehnder Interferometer (MZI) layout.

Fig. 1 Chip layout of a Mach-Zehnder Interferometer (MZI).#

The MZI we’ll create can be broken down into constituent parts. Simphony includes models for these building blocks below:

Mach-Zehnder Interferometer (MZI) block diagram

Fig. 2 Block diagram of a Mach-Zehnder Interferometer (MZI) with the components labeled.#

The grating couplers are the input and output for light in the circuit. The Y-branch can split and recombine light, and because the waveguides which carry light across the circuit are varying length, this produces interference when the light is recombined at the second Y-branch. We can now begin defining our circuit in Simphony using the components we have identified.

Models#

We’ll use models from the SiEPIC Ebeam PDK library (already included in simphony).

from simphony.libraries import siepic

As explained in the Introduction, models in SAX (and therefore in Simphony) are simply “callables” (functions) that return a dictionary of scattering parameters. The keys of that resulting dictionary are the port-to-port relationships. In Simphony, we follow the convention (output, input) for the keys of the dictionary, which is the same as the S-parameter matrix formulation (where \(S_{ij}\) is the scattering parameter representing the response at port \(j\) given a stimuli at port \(i\)).

Models in SAX must have default parameters in their function signatures; that is, no positional arguments are allowed.

siepic.grating_coupler?
Signature:
siepic.grating_coupler(
    wl: Union[float, jax.Array] = 1.55,
    pol: Literal['te', 'tm'] = 'te',
    thickness: float = 220.0,
    dwidth: float = 0,
) -> Dict[Tuple[str, str], jaxtyping.Complex[Array, '...']]
Docstring:
SiEPIC EBeam PDK grating coupler optimized for TE polarizations at
1550nm.

The grating coupler efficiently couples light from a fiber array positioned
above the chip into the circuit. For the TE mode, the angle is -25 degrees
[needs citation].

.. image:: /_static/images/ebeam_gc_te1550.png
    :alt: ebeam_bdc_te1550.png

Parameters
----------
wl : float or Array
    The wavelengths to evaluate at in microns.
pol : {"te", "tm"}
    Polarization of the input/output modes.
thickness : {210.0, 220.0, 230.0}
    Thickness of the grating coupler silicon in nm. Useful for simulating
    manufacturing variability.
dwidth : {-20.0, 0.0, 20.0}
    Change in width from nominal of the gratings. Representative of
    manufacturing variability. Must be one of -20, 0, or 20.

Raises
------
ValueError
    If `pol` is not 'te' or 'tm'.
ValueError
    If `thickness` is not one of 210, 220, or 230.
ValueError
    If `dwidth` is not one of -20, 0, or 20.

Notes
-----
See also the PDK documentation:
https://github.com/SiEPIC/SiEPIC_EBeam_PDK/blob/master/Documentation/SiEPIC_EBeam_PDK%20-%20Components%20with%20Models.docx
File:      ~/git/simphony/simphony/libraries/siepic/models.py
Type:      function

This means that they can always be called without arguments to inspect what their default return values are and to see what port names are provided.

siepic.grating_coupler()
{('o0', 'o0'): Array([-0.0307378-0.00345908j], dtype=complex128),
 ('o0', 'o1'): Array([0.75686856+0.02082852j], dtype=complex128),
 ('o1', 'o0'): Array([0.74360676+0.09760613j], dtype=complex128),
 ('o1', 'o1'): Array([0.0750638-0.02585451j], dtype=complex128)}

SAX also provides an explicit function to get the port names.

print(sax.get_ports(siepic.y_branch))
print(sax.get_ports(siepic.waveguide))
print(sax.get_ports(siepic.grating_coupler))
('port 1', 'port 2', 'port 3')
('o0', 'o1')
('o0', 'o1')

Note

Throughout all of Simphony and SAX, the "wl" argument is in microns, by convention. This is not strictly enforced, and various model libraries may follow different convention, so it’s always good to double check the documentation of the model you’re using (or, if you’re writing the model, create documentation for it)!

Writing the netlist#

A SAX circuit contains a netlist, which is a collection of component instances, their connections, and exposed ports; and a list of models, which makes it easy to swap out different models without rewriting your netlist to see how the circuit behavior changes.

The netlist is a dictionary with three fields:

  • "instances": A dictionary of instance names to human-readable component model strings (like “coupler” or “waveguide”). You will define the string-to-model mapping later.

  • "connections": A dictionary of ports to ports in the form "instance_name,port_name": "instance_name,port_name" (note there is no whitespace delimiting the instance from its port, just a comma).

  • "ports": Since a SAX circuit is basically a model itself, and could be used in other circuits, it has exposed ports. This field is a dictionary mapping port names of the composite object to the ports of its constituent instances.

We’ll create

mzi, info = sax.circuit(
    netlist={
        "instances": {
            "gc_in": "gc",
            "splitter": "ybranch",
            "long_wg": "waveguide",
            "short_wg": "waveguide",
            "combiner": "ybranch",
            "gc_out": "gc",
        },
        "connections": {
            "gc_in,o0": "splitter,port 1",
            "splitter,port 2": "long_wg,o0",
            "splitter,port 3": "short_wg,o0",
            "long_wg,o1": "combiner,port 2",
            "short_wg,o1": "combiner,port 3",
            "combiner,port 1": "gc_out,o0",
        },
        "ports": {
            "in": "gc_in,o1",
            "out": "gc_out,o1",
        },
    },
    models={
        "ybranch": siepic.y_branch,
        "waveguide": siepic.waveguide,
        "gc": siepic.grating_coupler,
    }
)

Simulation (using callables)#

sax.circuit() returns a tuple. The first element is another callable function. All parameters you call it with will be passed on to the models contained within the circuit, so long as they are named the same. In this way, a circuit itself can act as a model within another circuit.

Warning

It is important that your all models have common names for arguments. For example, all models that take a length parameter should all use the name length for that argument. Models that are wavelength dependent should all take the same wl keyword parameter (by convention). If you have your own model library, you can follow whatever convention you want, as long as it’s consistent.

The second element of the returned tuple is a information object that contains details about the circuit.

We can simulate the circuit by simply calling it with the appropriate arguments–in this case, the wavelengths we want to simulate at.

The circuit itself contains parameterized models. We can pass arguments targeting those models by passing a dictionary of keyword arguments and corresponding values when invoking the circuit function. The names you gave the instances in the netlist become the keyword arguments of the function. You can pass a dictionary which will be used to instantiate those components at simulation time.

If you’re unsure of the format of the settings dictionary you need to pass to the circuit, you can use SAX’s get_settings function.

sax.get_settings(mzi)
{'gc_in': {'wl': Array(1.55, dtype=float64),
  'pol': 'te',
  'thickness': Array(220., dtype=float64),
  'dwidth': Array(0., dtype=float64)},
 'splitter': {'wl': Array(1.55, dtype=float64),
  'pol': 'te',
  'thickness': Array(220., dtype=float64),
  'width': Array(500., dtype=float64)},
 'long_wg': {'wl': Array(1.55, dtype=float64),
  'pol': 'te',
  'length': Array(0., dtype=float64),
  'width': Array(500., dtype=float64),
  'height': Array(220., dtype=float64),
  'loss': Array(0., dtype=float64)},
 'short_wg': {'wl': Array(1.55, dtype=float64),
  'pol': 'te',
  'length': Array(0., dtype=float64),
  'width': Array(500., dtype=float64),
  'height': Array(220., dtype=float64),
  'loss': Array(0., dtype=float64)},
 'combiner': {'wl': Array(1.55, dtype=float64),
  'pol': 'te',
  'thickness': Array(220., dtype=float64),
  'width': Array(500., dtype=float64)},
 'gc_out': {'wl': Array(1.55, dtype=float64),
  'pol': 'te',
  'thickness': Array(220., dtype=float64),
  'dwidth': Array(0., dtype=float64)}}

This is slightly overkill in our case. We don’t want to update wl for every single object. We only want to change the length parameters for the long and short waveguides, and let everything else stay with the default or assume a global value.

We can set the "long_wg" and "short_wg" waveguides to different lengths alone. The waveguide model within the circuit takes a “length” argument, so that’s what we’ll use in our dictionary of parameters that we pass using keyword arguments. These keyword arguments correspond to the instance name in the netlist. By setting wl at the toplevel, it will trickle down to all models that require a wl parameter.

wl = jnp.linspace(1.5, 1.6, 1000)
S = mzi(wl=wl, long_wg={"length": 150.0}, short_wg={"length": 50.0})

S is now our evaluated s-parameter dictionary

We’re interested in the power transmitted from the input to the output, which is the magnitude squared of the s-parameter. We’ll square the magnitude of the "out,in" element of the resulting dictionary. Recall, too, that we renamed these “external” (unconnected) ports in the netlist when we created the circuit. It’s really easy to give the ports on your composite circuits meaningful names, and it makes your code much more readable.

Below, we plot both in linear and in log scale using matplotlib.

mag = jnp.abs(S["out", "in"])**2

fig, axs = plt.subplots(2, 1, sharex=True)
axs[0].plot(wl, mag)
axs[0].set_ylabel("Transmission")
axs[1].plot(wl, 10*jnp.log10(mag))
axs[1].set_ylabel("Transmission (dB)")
axs[1].set_xlabel("Wavelength (um)")
plt.suptitle("MZI Response")
plt.show()
../_images/mzi_31_0.png

Simulation (using convenience classes)#

Simphony also provides some handy classes for running simulations and doing further analysis on the results. We’ll use the ClassicalSim class to run a sweep simulation. This class takes a circuit and a list of wavelengths to simulate at. Additionally, we’ll add an ideal laser to the input and wavelength independent photodetectors to our output.

from simphony.classical import ClassicalSim

wl = jnp.linspace(1.5, 1.6, 1000)

sim = ClassicalSim(ckt=mzi, wl=wl, long_wg={"length": 150.0}, short_wg={"length": 50.0})
laser = sim.add_laser(ports=["in"], power=1.0)
detector = sim.add_detector(ports=["out"])

Now we simulate.

result = sim.run()

The dictionary of detectors in the simulation result have some handy functions for plotting the data.

result.detectors["out"].plot()
<Axes: xlabel='Wavelength (um)', ylabel='Power (mW)'>
../_images/mzi_37_1.png

Or, you can access the data directly from its .wl and .power attributes.

You can also access the underlying s-parameters directly. Note that once an input signal has been applied, the sdict is reduced to only contain the output ports as keys, instead of an input-output tuple for each relationship.

wl = result.wl
s = result.sdict["out"]

plt.plot(wl, jnp.abs(s)**2)
plt.title(r"MZI s-parameters $|S_{oi}|^2$")
plt.tight_layout()
plt.show()
../_images/mzi_39_0.png