Tidy3D mode solver#

Tidy3D includes an open-source FDFD mode solver for simulating waveguide modes.

Waveguides#

Guided electromagnetic modes are those with an effective index larger than that of the waveguide cladding.

In the example below, we have a silicon waveguide (n = 3.4) surrounded by a SiO₂ cladding (n = 1.44). For a waveguide with dimensions of 220 nm in height and 450 nm in width, the effective index is calculated to be 2.466.

import gdsfactory as gf
import gplugins.tidy3d as gt
import matplotlib.pyplot as plt
import numpy as np
from gdsfactory.generic_tech import get_generic_pdk

gf.config.rich_output()
PDK = get_generic_pdk()
PDK.activate()

nm = 1e-3
strip = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=0.5,
    core_thickness=0.22,
    slab_thickness=0.0,
    core_material="si",
    clad_material="sio2",
)
strip.plot_index()

<matplotlib.collections.QuadMesh object at 0x7fbc27747110>

../_images/3ecca40ccc2ffd4b6370b3027f95c673b58598c6b6989bafa65137ed6f392b70.png
strip.plot_grid()

../_images/45bed9b85622ab57ff22b7d6318e4d2deaeca9c8409026a1579a0ac6426f4b55.png
strip.plot_field(field_name="Ex", mode_index=0)  # TE
01:02:55 UTC WARNING: The group index was not computed. To calculate group      
             index, pass 'group_index_step = True' in the 'ModeSpec'.           

<matplotlib.collections.QuadMesh object at 0x7fbc255ed5d0>

../_images/ad46e3ea51c9bff85aa12a31209a652eaa7c3ddb447b5c7f20b70d6181b57630.png
strip.plot_field(field_name="Ex", mode_index=0, value="dB")  # TE

<matplotlib.collections.QuadMesh object at 0x7fbc251638d0>

../_images/f74137dac1f39503fdd176e90e0e312f0d4a4f6f5951cf2330ab6eacee32d747.png
strip.plot_field(field_name="Ey", mode_index=1)  # TM

<matplotlib.collections.QuadMesh object at 0x7fbc2519eb10>

../_images/7ab6dd8a992e5702cacec6ab822d825a42265ce5775734e0d126f579c855146d.png
strip.n_eff

array([2.51134734+4.42777628e-05j, 1.86463643+2.09422901e-04j])
rib = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=0.5,
    core_thickness=0.22,
    slab_thickness=0.15,
    core_material="si",
    clad_material="sio2",
)
rib.plot_index()
rib.n_eff
01:02:56 UTC WARNING: Mode field at frequency index 0, mode index 1 does not    
             decay at the plane boundaries.                                     

array([2.67427418+3.10179974e-05j, 2.50854926+4.47861069e-05j])

../_images/1a1f2057e6cfae6e9bc6ebb2338cba4d68ec920c24f56f90b2273ee612abbe6c.png
rib.plot_field(field_name="Ex", mode_index=0)  # TE

<matplotlib.collections.QuadMesh object at 0x7fbc24e56690>

../_images/4a2e16b79f36c665ca60d3e6aceaff88f4a69756962778586bb4a00e34fba10f.png
nitride = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=1.0,
    core_thickness=400 * nm,
    slab_thickness=0.0,
    core_material="sin",
    clad_material="sio2",
)
nitride.plot_index()
nitride.n_eff

array([1.64461788+8.05914279e-05j, 1.57796343+1.42713228e-04j])

../_images/cc244328d9c5c0fca517db939161ee0f9f58cb1c69551c9776c4022321fa4ebc.png
nitride.plot_field(field_name="Ex", mode_index=0)  # TE

<matplotlib.collections.QuadMesh object at 0x7fbc24c4de10>

../_images/c4b847dd3b9da375d4d0cf7232f0fa309719ea2d7171dad58b0c5fa510e2440c.png

Sweep width#

You can sweep the waveguide width and compute the modes.

By increasing the waveguide width, the waveguide supports many more TE and TM modes. Where TE modes have a dominant Ex field and TM modes have larger Ey fields.

Notice that waveguides wider than 0.450 um support more than one TE mode. Therefore the maximum width for single mode operation is 0.450 um.

strip = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=1.0,
    slab_thickness=0.0,
    core_material="si",
    clad_material="sio2",
    core_thickness=220 * nm,
    num_modes=4,
)
w = np.linspace(400 * nm, 1000 * nm, 7)
n_eff = gt.modes.sweep_n_eff(strip, core_width=w)
fraction_te = gt.modes.sweep_fraction_te(strip, core_width=w)

for i in range(4):
    plt.plot(w, n_eff.sel(mode_index=i).real, c="k")
    plt.scatter(
        w, n_eff.sel(mode_index=i).real, c=fraction_te.sel(mode_index=i), vmin=0, vmax=1
    )
plt.axhline(y=1.44, color="k", ls="--")
plt.colorbar().set_label("TE fraction")
plt.xlabel("Width of waveguide µm")
plt.ylabel("Effective refractive index")
plt.title("Effective index sweep")
01:02:59 UTC WARNING: Mode field at frequency index 0, mode index 3 does not    
             decay at the plane boundaries.                                     
01:03:01 UTC WARNING: Mode field at frequency index 0, mode index 3 does not    
             decay at the plane boundaries.                                     

Text(0.5, 1.0, 'Effective index sweep')

../_images/18ea2c314a595e4b474ad2cb2a6d972481c0330389c9aced0e229916204e9436.png

Exercises

  • What is the maximum width to support a single TE mode at 1310 nm?

  • For a Silicon Nitride (n=2) 400nm thick waveguide surrounded by SiO2 (n=1.44), what is the maximum width to support a single TE mode at 1550 nm?

  • For two 500x220nm Silicon waveguides surrounded by SiO2, what is the coupling length (100% coupling) for 200 nm gap?

Group index#

You can also compute the group index for a waveguide.

nm = 1e-3

strip = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=500 * nm,
    slab_thickness=0.0,
    core_material="si",
    clad_material="sio2",
    core_thickness=220 * nm,
    num_modes=4,
    group_index_step=10 * nm,
)
print(strip.n_group)
01:03:13 UTC WARNING: Mode field at frequency index 0, mode index 3 does not    
             decay at the plane boundaries.                                     
             WARNING: Mode field at frequency index 1, mode index 3 does not    
             decay at the plane boundaries.                                     
01:03:14 UTC WARNING: Mode field at frequency index 2, mode index 3 does not    
             decay at the plane boundaries.                                     
[4.17803969 4.08299706 2.71577378 1.50332985]

Bend modes#

You can compute bend modes specifying the bend radius.

strip_bend = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=500 * nm,
    core_thickness=220 * nm,
    slab_thickness=0.0,
    bend_radius=4,
    core_material="si",
    clad_material="sio2",
)
strip_bend.plot_field(field_name="Ex", mode_index=0)  # TE
01:03:15 UTC WARNING: Mode field at frequency index 0, mode index 1 does not    
             decay at the plane boundaries.                                     

<matplotlib.collections.QuadMesh object at 0x7fbc23a091d0>

../_images/e24925126a60943ed2673a5953f6f3781862e7c58ea467f16b934c351860df52.png

Bend loss#

You can also compute the losses coming from the mode mismatch from the bend into a straight waveguide. To compute the bend loss due to mode mismatch you can calculate the mode overlap of the straight mode and the bent mode. Because there are two mode mismatch interfaces the total loss due to mode mismatch will be squared (from bend to straight and from straight to bend).

from paper

radii = np.arange(4, 7)
bend = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=500 * nm,
    core_thickness=220 * nm,
    core_material="si",
    clad_material="sio2",
    num_modes=1,
    bend_radius=radii.min(),
)
mismatch = gt.modes.sweep_bend_mismatch(bend, radii)

plt.plot(radii, 10 * np.log10(mismatch))
plt.title("Strip waveguide bend")
plt.xlabel("Radius (μm)")
plt.ylabel("Mismatch (dB)")

Text(0, 0.5, 'Mismatch (dB)')

../_images/32bdf5e2f250f17fbdec4ff7b515fc244538c39894ccf99c131a1b37a0dfaedf.png
dB_cm = 2  # dB/cm
length = 2 * np.pi * radii * 1e-6
propagation_loss = dB_cm * length * 1e2
print(propagation_loss)

plt.title("Bend90 loss for TE polarization")
plt.plot(radii, -10 * np.log10(mismatch), ".", label="mode loss")
plt.plot(radii, propagation_loss, ".", label="propagation loss")
plt.xlabel("bend radius (um)")
plt.ylabel("Loss (dB)")
plt.legend()
[0.00502655 0.00628319 0.00753982]

<matplotlib.legend.Legend object at 0x7fbc238de810>

../_images/d7d7b333ac71543cac5ab0e41eebb19587ecfc67e6c504bf84cdfa3461ee62d0.png
rib = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=1000 * nm,
    core_thickness=220 * nm,
    slab_thickness=110 * nm,
    bend_radius=15,
    core_material="si",
    clad_material="sio2",
)
rib.plot_field(field_name="Ex", mode_index=0)  # TE
01:03:20 UTC WARNING: Mode field at frequency index 0, mode index 1 does not    
             decay at the plane boundaries.                                     

<matplotlib.collections.QuadMesh object at 0x7fbc237a5a10>

../_images/66b23f95d8aead630b293fc081e84bb73649943e2b20eb3f4e57cbec30f45ba9.png
nitride_bend = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=1000 * nm,
    core_thickness=400 * nm,
    slab_thickness=0.0,
    bend_radius=30,
    core_material="sin",
    clad_material="sio2",
)
nitride_bend.plot_field(field_name="Ex", mode_index=0, value="abs")  # TE

<matplotlib.collections.QuadMesh object at 0x7fbc24cffb90>

../_images/1ce94a598c2cd49b1aefa16149df3d479a83fbb0112629b0743ce8e4c5b0f924.png
radii = np.array([30, 35, 40])
bend = gt.modes.Waveguide(
    wavelength=1.55,
    core_width=1000 * nm,
    core_thickness=400 * nm,
    core_material="sin",
    clad_material="sio2",
    num_modes=1,
    bend_radius=radii.min(),
)
mismatch = gt.modes.sweep_bend_mismatch(bend, radii)
dB_cm = 2  # dB/cm
length = 2 * np.pi * radii * 1e-6
propagation_loss = dB_cm * length * 1e2
print(propagation_loss)

plt.title("Bend90 loss for TE polarization")
plt.plot(radii, -10 * np.log10(mismatch), ".", label="mode loss")
plt.plot(radii, propagation_loss, ".", label="propagation loss")
plt.xlabel("bend radius (um)")
plt.ylabel("Loss (dB)")
plt.legend()
[0.03769911 0.0439823  0.05026548]

<matplotlib.legend.Legend object at 0x7fbc24fee6d0>

../_images/1771c6818973f0882c5a8364623ddca3d6d87b8c7dd54658251c4d13a609cb8f.png

Exercises

  • For a 500nm wide 220nm thick Silicon waveguide surrounded by SiO2, what is the minimum bend radius to have less than 0.04dB loss for TE polarization at 1550nm?

  • For a 500nm wide 220nm thick Silicon waveguide surrounded by SiO2, what is the minimum bend radius to have 99% power transmission for TM polarization at 1550nm?

Waveguide coupler#

You can also compute the modes of a waveguide coupler.

       ore_width[0]  core_width[1]
        <------->     <------->
         _______       _______   _
        |       |     |       | |
        |       |     |       |
        |       |_____|       | | core_thickness
        |slab_thickness       |
        |_____________________| |_
                <----->
                  gap


c = gt.modes.WaveguideCoupler(
    wavelength=1.55,
    core_width=(500 * nm, 500 * nm),
    gap=200 * nm,
    core_thickness=220 * nm,
    slab_thickness=100 * nm,
    core_material="si",
    clad_material="sio2",
)
c.plot_index()

<matplotlib.collections.QuadMesh object at 0x7fbc25021cd0>

../_images/129849ee6ebd61ca6748542adb5fff0f2aa87b314458d73c03b79764a5efdcc3.png
c.plot_field(field_name="Ex", mode_index=0)  # TE

<matplotlib.collections.QuadMesh object at 0x7fbc23810b90>

../_images/f7e785736c6bfffbb906ed5d0062483405b4b1e274d91f787b80c6d17a66fa28.png
c.plot_field(field_name="Ex", mode_index=1)  # TE

<matplotlib.collections.QuadMesh object at 0x7fbc231eb310>

../_images/da8d983867fe1a56947dad723052e2932a12dc585fe66c338febc6b6499f5819.png
coupler = gt.modes.WaveguideCoupler(
    wavelength=1.55,
    core_width=(450 * nm, 450 * nm),
    core_thickness=220 * nm,
    core_material="si",
    clad_material="sio2",
    num_modes=4,
    gap=0.1,
)

print("\nCoupler:", coupler)
print("Effective indices:", coupler.n_eff)
print("Mode areas:", coupler.mode_area)
print("Coupling length:", coupler.coupling_length())

gaps = np.linspace(0.05, 0.15, 11)
lengths = gt.modes.sweep_coupling_length(coupler, gaps)

_, ax = plt.subplots(1, 1)
ax.plot(gaps, lengths)
ax.set(xlabel="Gap (μm)", ylabel="Coupling length (μm)")
ax.legend(["TE", "TM"])
ax.grid()
Coupler: WaveguideCoupler(wavelength=array(1.55), core_width=['0.45', '0.45'], core_thickness='0.22', core_material='si', clad_material='sio2', box_material=None, slab_thickness='0.0', clad_thickness=None, box_thickness=None, side_margin=None, sidewall_angle='0.0', sidewall_thickness='0.0', sidewall_k='0.0', surface_thickness='0.0', surface_k='0.0', bend_radius=None, num_modes='4', group_index_step='False', precision='double', grid_resolution='20', max_grid_scaling='1.2', cache_path='/github/home/.gdsfactory/modes', overwrite='False', gap='0.1')
Effective indices: [2.4637647 +6.57552457e-05j 2.39007229+5.06214923e-05j
 1.9225165 +1.99036730e-04j 1.71420814+2.37015946e-04j]
Mode areas: [0.31003254 0.33258301 0.57286555 0.59002858]
Coupling length: [10.5166863   3.72044606]

../_images/6e5b17741bd2211f6602df1739112d11b90a14153f1c2e3a1309cb5f14561c19.png