Mode solver simulations are typically fast because they are intrinsically 2D problems. However, in many cases, multiple mode solver simulations are required to understand mode behaviors in a waveguide, which can still take a significant amount of time if run sequentially. Tidy3D supports running a batch of mode solver tasks in parallel, very similar to running a batch of FDTD simulations.
This notebook demonstrates how to perform parallel mode solving using the example of calculating mode dispersion as a function of the lithium niobate (LN) strip waveguide width.
For a more focused mode solver tutorial, please refer to here.
import numpy as np
import matplotlib.pyplot as plt
import tidy3d as td
import tidy3d.web as web
from tidy3d.plugins import waveguide
First, we define the basic parameters and materials to be used. Specifically, we are solving for a lithium niobate (LN) strip waveguide with a sidewall angle of 20 degrees and a thickness of 400 nm at a wavelength of 1550 nm. The waveguide sits on a silicon dioxide BOX, and the top cladding is air.
For convenience, we will directly use the available materials from the material library. When defining the LN material, we need to specify the optical axis direction since LN is anisotropic, i.e. td.material_library['LiNbO3']['Zelmon1997'](optical axis)
, where optical axis = 1, 2, or 3
corresponds to the $x$-, $y$-, and $z$-axis of the simulation domain. In this case, the LN is X-cut with the mode propagating in the Y direction ($x$-axis of the simulation coordinate), so we set the optical axis to 1
, meaning the extraordinary refractive index is aligned with the $y$-axis of the simulation coordinates.
lda0 = 1.55 # wavelength of interest
h_LN = 0.4 # lithium niobate waveguide thickness
sidewall_angle = 20 # sidewall angle in degree
LN = td.material_library["LiNbO3"]["Zelmon1997"](1)
SiO2 = td.material_library["SiO2"]["Palik_Lossless"]
To perform the sweep of the waveguide width, we need to generate a dictionary of mode solver tasks each with a different waveguide width configuration. To do so conveniently, we will use the waveguide plugin as it eliminates the need to manually define the waveguide and BOX structures as well as the simulation domain and mode solving plane.
w_LN_range = np.arange(0.5, 4.1, 0.1) # waveguide width range
n_mode = 10 # number of modes to solve for
mode_solvers = {}
for w_LN in w_LN_range:
# using the waveguide plugin to create a mode solver for each waveguide width configuration
strip_waveguide = waveguide.RectangularDielectric(
wavelength=lda0,
core_width=w_LN,
core_thickness=h_LN,
sidewall_angle=np.deg2rad(sidewall_angle),
core_medium=LN,
clad_medium=td.Medium(),
box_medium=SiO2,
mode_spec=td.ModeSpec(
num_modes=n_mode, target_neff=LN.nk_model(td.C_0 / lda0)[0], num_pml=(12, 12)
),
)
# add the mode solver to the dictionary
mode_solvers[f"width={w_LN:.2f}"] = strip_waveguide.mode_solver
Once the list of mode solvers is created, we simply use the run_async
method to run the entire batch of tasks in parallel. Alternatively, one can create a Batch object first and run it with the run() method. Batch provides a better management of the simulation files since one can save the batch information to file and load the batch at a later time.
# run mode solvers in parallel
batch_results = web.run_async(
simulations=mode_solvers, # dictionary of mode solvers
path_dir="data", # path to store the result files
)
# alternative way to run the mode solvers in parallel is shown below
# batch = web.Batch(simulations=mode_solvers)
# batch_results = batch.run(path_dir="data")
Output()
18:39:41 UTC Started working on Batch containing 36 tasks.
18:39:56 UTC Maximum FlexCredit cost: 0.144 for the whole batch.
Use 'Batch.real_cost()' to get the billed FlexCredit cost after the Batch has completed.
Output()
18:40:03 UTC Batch complete.
Output()
After the batch of tasks is finished, we are ready to extract the results. We can plot all the mode profiles to inspect each mode at each waveguide width. As a demonstration, we plot all 10 modes with the widest waveguide.
# plot all mode profiles of the widest waveguide
fig, ax = plt.subplots(5, 2, figsize=(7, 10), tight_layout=True)
for i in range(n_mode):
mode_solvers[f"width={w_LN_range[-1]:.2f}"].plot_field(
"E", "abs", mode_index=i, ax=ax[i % 5][i // 5]
)
ax[i % 5][i // 5].set_title(f"mode_index={i}")

From the mode profiles, we can distinguish fundamental modes and higher order modes. To be more quantitative, we can extract other mode information. For example, we can extract the effective indices and polarization fractions of the modes.
# extract the effective index and polarization fraction
n_eff = np.array([result.n_eff[0] for _, result in batch_results.items()]).T
te_fraction = np.array([result.pol_fraction.te.data[0] for _, result in batch_results.items()]).T
Once we extract the information, we can plot them to visualize the mode dispersion. Here we plot the effective indices first as curves and then as scatterers with color, where the color indicates the polarization. Red corresponds to quasi-TE modes and blue corresponds to quasi-TM modes. We can see at a few places mode hybridization occurs. For example, at around 1.5 μm waveguide width, the TM0 and TE1 modes hybridize. This hybridization is often utilized to make mode converters for polarization control, as demonstrated in the polarization splitter rotator example.
# plot each mode dispersion curve and polarization as colored scatterers
plt.figure(figsize=(10, 5))
for i in range(n_mode):
plt.plot(w_LN_range, n_eff[i], c="black")
plt.scatter(w_LN_range, n_eff[i], c=te_fraction[i], cmap="coolwarm", s=40, vmin=0, vmax=1)
plt.ylim(1.5, 1.85)
plt.xlabel("Waveguide width (μm)")
plt.ylabel("Effective index")
plt.colorbar()
plt.show()
