View this project in Tidy3D Web App.

This is a simple example of using Tidy3D to perform a 2D simulation of a ring resonator side coupled to a dielectric waveguide.

diagram

With a center wavelength of 500 nm and 10 nm resolution, this is a challenging FDTD problem because of the large simulation size and the large number of time steps needed to capture the resonance of the ring. However, it can be easily handled with Tidy3D.

For more integrated photonic examples such as the 8-Channel mode and polarization de-multiplexer, the broadband bi-level taper polarization rotator-splitter, and the broadband directional coupler, please visit our examples page. If you are interested in ring resonator simulation in 3D, please check out the THz filter case study.

If you are new to the finite-difference time-domain (FDTD) method, we highly recommend going through our FDTD101 tutorials. FDTD simulations can diverge due to various reasons. If you run into any simulation divergence issues, please follow the steps outlined in our troubleshooting guide to resolve it.

[1]:

# standard python imports
import numpy as np
from numpy import random
import matplotlib.pyplot as plt

# tidy3D import
import tidy3d.web as web
import tidy3d as td

Initial setup¶

Our ring resonator will include a ring centered at (0,0) with a waveguide just above the ring spanning the x direction.

schematic

[2]:

# define geometry
wg_width = 0.25
couple_width = 0.13
ring_radius = 2.5
ring_wg_width = 0.25
wg_spacing = 0.5
buffer = 1.0

# compute quantities based on geometry parameters
x_span = 2 * wg_spacing + 2 * ring_radius + 2 * buffer
y_span = 2 * ring_radius + ring_wg_width + wg_width + couple_width + 2 * buffer
sim_center_y = (couple_width + wg_width) / 2
wg_insert_x = ring_radius + wg_spacing
wg_center_y = ring_radius + ring_wg_width / 2.0 + couple_width + wg_width / 2.0

Define frequency parameters for sources and monitors based on the wavelength range of interest. Since in this case study we are investigating the system's response over a wide range of frequencies we will use a broadband mode source that takes into account frequency-dependent character of waveguide modes. The Tidy3D's broadband feature is designed in such a way that it produces the most accurate results in the range (freq0 - 1.5 * fwidth, freq0 + 1.5 * fwidth), thus we define parameters freq0 and fwidth such that the frequency range of interest is exactly covered by (freq0 - 1.5 * fwidth, freq0 + 1.5 * fwidth).

[3]:

# wavelength range of interest
lambda_beg = 0.5
lambda_end = 0.6

# define pulse parameters
freq_beg = td.C_0 / lambda_end
freq_end = td.C_0 / lambda_beg
freq0 = (freq_beg + freq_end) / 2
fwidth = (freq_end - freq0) / 1.5

min_steps_per_wvl = 30
run_time = 15e-12

Define materials.

[4]:

n_bg = 1.0
n_solid = 1.5
background = td.Medium(permittivity=n_bg**2)
solid = td.Medium(permittivity=n_solid**2)

Define structures.

[5]:

# waveguide
waveguide = td.Structure(
    geometry=td.Box(
        center=[0, wg_center_y, 0],
        size=[td.inf, wg_width, td.inf],
    ),
    medium=solid,
    name="waveguide",
)

# outer ring boundary
outer = td.Cylinder(
    center=[0, 0, 0],
    axis=2,
    radius=ring_radius + ring_wg_width / 2.0,
    length=td.inf,
)

# inner ring boundary
inner = td.Cylinder(
    center=[0, 0, 0],
    axis=2,
    radius=ring_radius - ring_wg_width / 2.0,
    length=td.inf,
)

ring = td.Structure(
    geometry=outer - inner,
    medium=solid,
    name="outer_ring",
)

Compute and visualize the waveguide modes.

[6]:

from tidy3d.plugins.mode import ModeSolver
from tidy3d.plugins.mode.web import run as run_ms

mode_plane = td.Box(
    center=[-wg_insert_x, wg_center_y, 0],
    size=[0, 2, td.inf],
)

sim_modesolver = td.Simulation(
    center=[0, sim_center_y, 0],
    size=[x_span, y_span, 0],
    grid_spec=td.GridSpec.auto(
        min_steps_per_wvl=min_steps_per_wvl, wavelength=td.C_0 / freq0
    ),
    structures=[waveguide],
    run_time=1e-12,
    boundary_spec=td.BoundarySpec.all_sides(boundary=td.Periodic()),
    medium=background,
)

mode_spec = td.ModeSpec(num_modes=2)
mode_solver = ModeSolver(
    simulation=sim_modesolver, plane=mode_plane, mode_spec=mode_spec, freqs=[freq0]
)
mode_data = run_ms(mode_solver)

[15:53:03] Mode solver created with                                             
           task_id='fdve-6aeecfd9-acaa-4fb7-a055-db9cfecf5f17v1',               
           solver_id='mo-4c6d027a-be37-47c5-a370-2b2920daf9c2'.

Output()

Output()

[15:53:04] Mode solver status: queued

[15:53:06] Mode solver status: running

[15:53:17] Mode solver status: success

Output()

Tidy3D is warning us that the simulation does not contain a source. However, since this simulation is used to construct the mode solver and will not be run directly, we can ignore this warning.

[7]:

f, ((ax1, ax2, ax3), (ax4, ax5, ax6)) = plt.subplots(
    2, 3, tight_layout=True, figsize=(10, 6)
)
mode_data.Ex.sel(mode_index=0, f=freq0).abs.plot(ax=ax1)
mode_data.Ey.sel(mode_index=0, f=freq0).abs.plot(ax=ax2)
mode_data.Ez.sel(mode_index=0, f=freq0).abs.plot(ax=ax3)
mode_data.Ex.sel(mode_index=1, f=freq0).abs.plot(ax=ax4)
mode_data.Ey.sel(mode_index=1, f=freq0).abs.plot(ax=ax5)
mode_data.Ez.sel(mode_index=1, f=freq0).abs.plot(ax=ax6)
ax1.set_title("|Ex|: mode_index=0")
ax2.set_title("|Ey|: mode_index=0")
ax3.set_title("|Ez|: mode_index=0")
ax4.set_title("|Ex|: mode_index=1")
ax5.set_title("|Ey|: mode_index=1")
ax6.set_title("|Ez|: mode_index=1")

mode_data.to_dataframe()

[7]:

		wavelength	n eff	k eff	TE (Ey) fraction	wg TE fraction	wg TM fraction	mode area
f	mode_index
5.496195e+14	0	0.545455	1.348767	0.0	2.867437e-13	1.000000	0.882227	0.309920
5.496195e+14	1	0.545455	1.271590	0.0	1.000000e+00	0.789246	1.000000	0.454015

No description has been provided for this image

From the above plots, we see that

mode_index=0 corresponds to exciting 0-th order TM mode (E=Ez) and

mode_index=1 corresponds to exciting 0-th order TE mode (E=Ey).

We can therefore switch the mode index accordingly based on our polarization.

Let's select Ey and create the source for it. As mentioned before, we will use a broadband mode source that will approximate the chosen mode's frequency dependence in the range (freq0 - 1.5 * fwidth, freq0 + 1.5 * fwidth) using 11 Chebyshev polynomials.

[8]:

mode_source = td.ModeSource(
    size=mode_plane.size,
    center=mode_plane.center,
    source_time=td.GaussianPulse(freq0=freq0, fwidth=fwidth),
    mode_spec=td.ModeSpec(num_modes=2),
    mode_index=1,
    direction="+",
    num_freqs=11,  # using 11 (Chebyshev) points to approximate frequency dependence
)

In addition, let's monitor both the fields in plane as well as the output mode amplitudes into the fundamental TE mode.

[9]:

# monitor steady state fields at central frequency over whole domain
field_monitor = td.FieldMonitor(
    center=[0, 0, 0], size=[td.inf, td.inf, 0], freqs=[freq0], name="field"
)


# monitor the mode amps on the output waveguide
lambdas_measure = np.linspace(lambda_beg, lambda_end, 1001)
freqs_measure = td.C_0 / lambdas_measure[::-1]

mode_monitor = td.ModeMonitor(
    size=mode_plane.size,
    center=mode_plane.center,
    freqs=freqs_measure,
    mode_spec=td.ModeSpec(num_modes=2),
    name="mode",
)

# lets reset the center to the on the right hand side of the simulation though
mode_monitor = mode_monitor.copy(update=dict(center=[wg_insert_x, wg_center_y, 0]))

[16:16:36] WARNING: Default value for the field monitor           monitor.py:261
           'colocate' setting has changed to 'True' in Tidy3D                   
           2.4.0. All field components will be colocated to the                 
           grid boundaries. Set to 'False' to get the raw fields                
           on the Yee grid instead.

Define simulation. To make the simulation 2D, we can just set the simulation size in one of the dimensions to be 0. However, note that we still have to define a grid size in that direction. Additionally, we lower the shutoff factor to make sure the high-frequency response of this resonant system is accurately resolved.

[10]:

# create simulation
sim = td.Simulation(
    center=[0, sim_center_y, 0],
    size=[x_span, y_span, 0],
    grid_spec=td.GridSpec.auto(
        min_steps_per_wvl=min_steps_per_wvl, wavelength=td.C_0 / freq0
    ),
    structures=[waveguide, ring],
    sources=[mode_source],
    monitors=[field_monitor, mode_monitor],
    run_time=run_time,
    boundary_spec=td.BoundarySpec(
        x=td.Boundary.pml(), y=td.Boundary.pml(), z=td.Boundary.periodic()
    ),
    medium=background,
    shutoff=1e-9,
)

Visualize structure, source, and modes.

[11]:

# plot the two simulations
fig, ax = plt.subplots(1, 1, figsize=(6, 6))
sim.plot_eps(z=0.01, ax=ax)
plt.show()

Run Simulation¶

Run simulations on our server.

[12]:

# use function above to run simulation
sim_data = web.run(
    sim, task_name="ring_resonator", path="data/simulation_data.hdf5", verbose=True
)

           Created task 'ring_resonator' with task_id              webapi.py:188
           'fdve-8d63a91b-81e9-4fcb-bd8a-51e7fd091adev1'.

           View task using web UI at                               webapi.py:190
           'https://tidy3d.simulation.cloud/workbench?taskId=fdve-              
           8d63a91b-81e9-4fcb-bd8a-51e7fd091adev1'.

[16:16:37] status = queued                                         webapi.py:361

[16:16:46] status = preprocess                                     webapi.py:355

[16:16:51] Maximum FlexCredit cost: 0.098. Use                     webapi.py:341
           'web.real_cost(task_id)' to get the billed FlexCredit                
           cost after a simulation run.

           starting up solver                                      webapi.py:377

           running solver                                          webapi.py:386

           To cancel the simulation, use 'web.abort(task_id)' or   webapi.py:387
           'web.delete(task_id)' or abort/delete the task in the                
           web UI. Terminating the Python script will not stop the              
           job running on the cloud.

[16:20:08] early shutoff detected, exiting.                        webapi.py:404

           status = postprocess                                    webapi.py:419

[16:20:29] status = success                                        webapi.py:426

[16:20:34] loading SimulationData from data/simulation_data.hdf5   webapi.py:590

[13]:

# visualize normalization run
fig, (ax1, ax2) = plt.subplots(1, 2, tight_layout=True, figsize=(15, 6))

ax1 = sim_data.plot_field("field", "Hz", val="real", z=0, ax=ax1)
ax2 = sim_data.plot_field("field", "S", val="abs", z=0, ax=ax2)

plt.show()

Analyze Spectrum¶

Now let's analyze the mode amplitudes in the output waveguide.

First, let's grab the data to inspect it.

[14]:

sim_data["mode"].amps

As we see, the mode amplitude data is complex-valued with three 3 dimensions:

index into the mode order returned by solver (remember, we wanted mode_index=1 for fundamental TE).
direction of the propagation (for decomposition).
frequency.

Let's select into the first two dimensions to get mode amplitudes as a function of frequency.

[15]:

transmission_amps = sim_data["mode"].amps.sel(mode_index=1, direction="+")

Now let's plot the data.

[16]:

f, ax = plt.subplots(figsize=(10, 5))
transmission_amps.abs.plot.line(x="f", ax=ax, label="absolute value")
# transmission_amps.real.abs.plot.line(x='f', ax=ax, label='|real|')
ax.legend()
ax.set_title("transmission amplitudes into fundamental TE (forward)")
ax.set_ylim(0, None)
ax.set_xlim(freqs_measure[0], freqs_measure[-1])
plt.show()

Learning Center

2D ring resonator¶

Initial setup¶

Run Simulation¶

Analyze Spectrum¶

More Examples