Inverse taper edge coupler#

Note: The cost of running the entire notebook is higher than 1 FlexCredit.

Fiber-to-chip couplers are essential components in photonic integrated circuits (PICs) that need to couple in/out light using optical fibers. However, the huge mismatch between the mode sizes of a single mode fiber (\(\approx 10\mu m\) in diameter) and a silicon wire waveguide (220 x 500 nm) can cause significant coupling loss if they are coupled directly. For this reason, the design of the fiber-to-chip coupling devices is critical for the whole system’s performance and needs to be carefully accomplished. There are two main categories of fiber-to-chip optical coupling devices: grating couplers, where the optical fiber is positioned perpendicularly to the chip plane; and edge couplers, which allow in-plane light coupling by the direct alignment between the fiber and the integrated waveguide facets.

In this notebook, we will show an example of using Tidy3D to evaluate the performance of edge couplers built using inverted taper mode transformers of linear, quadratic, and exponential profiles. We will also see how to set up a Gaussian beam to simulate the field launched by a lensed fiber and the use of Batch simulations to perform parameter sweeps.

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 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.

We start by importing the necessary modules:

[1]:
# standard python imports
import numpy as np
import matplotlib.pyplot as plt

# Tidy3D imports
import tidy3d as td
import tidy3d.web as web

Simulation Setup#

The inverted taper structure is made using an SOI platform consisting of a silicon (Si) substrate, a 3 \(\mu m\) thick silicon dioxide (SiO2) layer, and a 220 nm thick silicon layer covered by 3 \(\mu m\) thick silicon dioxide. As shown in the image below, the inverted taper is defined at the 220 nm thick Si layer, and its main geometric parameters are also defined in the figure. An optical fiber mode with a Gaussian profile (2.5 \(\mu m\) spot-size) is focused at the taper tip and coupled into the output waveguide.

Schematic of the inverted taper edge coupler

Definition of the simulation wavelength and materials.

[2]:
wl_c = 1.55  # Central wavelength.
wl_bw = 0.100  # Wavelength bandwidth.
wl_n = 101  # Number of wavelength points to compute the output data.

mat_si = td.Medium(permittivity=3.48**2)  # Taper and substrate material.
mat_sio2 = td.Medium(permittivity=1.44**2)  # Box and cladding material.
mat_air = td.Medium(permittivity=1.00)  # External medium material.

Definition of the taper size, lensed fiber parameters, and SOI structure.

[3]:
taper_l = 25  # Taper length.
taper_w_in = 0.19  # Taper tip width.
taper_w_out = 0.50  # Taper output width.
taper_t = 0.22  # Taper thickness.

# Spot size of the gaussian mode launched by the lensed fiber at the taper tip.
spot_size = 2.5

box_thick = 3.0  # Silicon dioxide box layer.
clad_thick = 3.0  # Silicon dioxide layer covering the taper.

Calculated parameters and constants used in the simulation set up.

[4]:
# Extra space around the taper at x,y directions.
pad_x = 1 * wl_c
pad_y = 1.5 * wl_c

# Wavelength and frequency values.
wl_range = np.linspace(wl_c - wl_bw / 2, wl_c + wl_bw / 2, wl_n)
freq_c = td.C_0 / wl_c
freq_range = td.C_0 / wl_range
freq_width = 0.5 * (np.max(freq_range) - np.min(freq_range))
run_time = 30 / freq_width

# Large number to be used in replacement of td.inf when necessary.
_inf = 1e3

Function to create the inverted taper geometry based on the work of G. Ren et al. (2013).

[5]:
# List of inverted tapers names.
tap_names = ["linear", "quadratic", "exponential"]


def get_taper(
    taper_shape="linear",
    init_coord=[pad_x, taper_w_in / 2],
    end_coord=[pad_x + taper_l, taper_w_out / 2],
):
    """Return a polygon composed of the taper vertices given the desired shape."""

    x0 = init_coord[0]
    x1 = end_coord[0]
    y0 = init_coord[1]
    y1 = end_coord[1]
    x = np.linspace(x0, x1, 41)

    tap_name = "linear_tap"
    if taper_shape == "quadratic":
        alpha = ((x - x0) / (x1 - x0)) ** 2
        tap_name = "quadratic_tap"
    elif taper_shape == "exponential":
        alpha = np.expm1((x - x0) / (x1 - x0)) / np.expm1(1)
        tap_name = "exponential_tap"
    elif taper_shape == "linear":
        alpha = (x - x0) / (x1 - x0)
    else:
        print(
            "taper_shape is neither 'linear', 'quadratic', or 'exponential'!\n"
            + "Linear taper shape returned!"
        )
        alpha = (x - x0) / (x1 - x0)

    y = y0 + alpha * (y1 - y0)
    upper_bound = [[xv, yv] for xv, yv in zip(x, y)] + [[_inf, y1]]
    lower_bound = [[_inf, -y1]] + [[xv, -yv] for xv, yv in zip(x[::-1], y[::-1])]
    taper_polygon = upper_bound + lower_bound

    # Inverted taper structure using a PolySlab.
    taper = td.Structure(
        geometry=td.PolySlab(
            vertices=taper_polygon, axis=2, slab_bounds=(-taper_t / 2, taper_t / 2)
        ),
        medium=mat_si,
        name=tap_name,
    )
    return taper

Function to create the simulation object. Since the box and cladding oxide layers are considerably thick, the silicon substrate and the upper air layer will be removed from the simulation, so being possible to take advantage of the symmetry along the z-axis. In addition, we will insert the Gaussian source directly into the taper tip to reduce the simulation cost further.

[6]:
def get_simulations(
    tap_length=taper_l,
    tap_w_in=taper_w_in,
    tap_w_out=taper_w_out,
    tap_t=taper_t,
    tap_names=tap_names,
):
    """Return a dict of simulation objects for the content of 'tap_names'."""

    size_x = tap_length + 2 * pad_x
    size_y = tap_w_out + 2 * pad_y
    size_z = box_thick + clad_thick + tap_t

    # Gaussian source focused at the taper tip.
    # The Gaussian source must be placed at a plane within a uniform medium.
    gauss_source = td.GaussianBeam(
        center=(0.99 * pad_x, 0, 0),
        size=(0, size_y, size_z),
        source_time=td.GaussianPulse(freq0=freq_c, fwidth=freq_width),
        pol_angle=0,
        direction="+",
        num_freqs=7,
        waist_radius=spot_size / 2,
    )

    # Field monitor to visualize the fields throughout the length of the taper.
    field_monitor_xy = td.FieldMonitor(
        center=(size_x / 2, 0, 0),
        size=(size_x, size_y, 0),
        freqs=[freq_c],
        name="field_xy",
    )

    # Field monitor to visualize the Gaussian input fields.
    field_input = td.FieldMonitor(
        center=(pad_x, 0, 0),
        size=(0, size_y, size_z),
        freqs=[freq_c],
        name="field_input",
    )

    # Field monitor to visualize the fields at the output waveguide.
    field_output = td.FieldMonitor(
        center=(size_x - pad_x / 2, 0, 0),
        size=(0, size_y, size_z),
        freqs=[freq_c],
        name="field_output",
    )

    # Mode monitor to get the power coupled into the fundamental TE mode.
    mode_spec = td.ModeSpec(num_modes=1)
    mode_monitor = td.ModeMonitor(
        center=(size_x - pad_x / 2, 0, 0),
        size=(0, size_y, size_z),
        freqs=freq_range,
        mode_spec=mode_spec,
        name="mode_0",
    )

    # Silicon dioxide box + cladding layers
    sio2_medium = td.Structure(
        geometry=td.Box.from_bounds(
            rmin=(pad_x, -_inf, -_inf), rmax=(_inf, _inf, _inf)
        ),
        medium=mat_sio2,
    )

    # Simulation definition
    sim_tap = []
    for t_name in tap_names:
        taper_poly = get_taper(
            t_name,
            init_coord=[pad_x, tap_w_in / 2],
            end_coord=[size_x - pad_x, tap_w_out / 2],
        )

        sim_tap.append(
            td.Simulation(
                center=(size_x / 2, 0, 0),
                size=(size_x, size_y, size_z),
                medium=mat_air,
                grid_spec=td.GridSpec.auto(min_steps_per_wvl=15, wavelength=wl_c),
                structures=[sio2_medium, taper_poly],
                sources=[gauss_source],
                normalize_index=0,
                monitors=[field_monitor_xy, field_input, field_output, mode_monitor],
                boundary_spec=td.BoundarySpec.all_sides(boundary=td.PML(num_layers=20)),
                symmetry=(0, -1, 1),
                run_time=run_time,
            )
        )
    sims = {sim_name: sim for sim_name, sim in zip(tap_names, sim_tap)}
    return sims

Comparison between inverted tapers with different shapes.#

We start our analysis by comparing the response of inverted tapers with 3 different geometries: linear, quadratic, and exponential. Then, before running the simulations, we will visualize the shapes of the inverted tapers and examine the simulation set up to see if everything is correct.

[7]:
# Get the list of simulation objects.
sim_tap = get_simulations(tap_length=taper_l)

13:33:51 PDT WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
[8]:
# Structure visualization.
fig, axs = plt.subplots(3, 1, tight_layout=True, figsize=(8, 6))
for sim_s, tap_n, ax in zip(sim_tap.values(), sim_tap.keys(), axs):
    sim_s.plot_structures(z=0, ax=ax)
    ax.set_xlim(0, 2 * pad_x + taper_l)
    ax.set_ylim(-1.5 * taper_w_out / 2, 1.5 * taper_w_out / 2)
    ax.set_aspect("auto")  # Used to better visualize the shapes.
    ax.set_title(tap_n)
plt.show()

../_images/notebooks_EdgeCoupler_17_0.png
[9]:
# Simulation set up visualization.
fig = plt.figure(tight_layout=True, figsize=(8, 6))
gs = fig.add_gridspec(2, 2)
ax1 = fig.add_subplot(gs[0, :])
ax2 = fig.add_subplot(gs[1, 0])
ax3 = fig.add_subplot(gs[1, 1])
sim_tap["linear"].plot(z=0, ax=ax1)
sim_tap["linear"].plot(x=pad_x, ax=ax2)
sim_tap["linear"].plot(x=pad_x + taper_l + pad_x / 2, ax=ax3)
plt.show()

../_images/notebooks_EdgeCoupler_18_0.png

Now, we will run the simulations using a batch object.

[10]:
# Initialize a simulation batch.
batch = web.Batch(simulations=sim_tap, verbose=True)

# Get the estimated simulation cost.
tot_cost = 0
for bat in batch.get_info().values():
    sim_name = bat.taskName
    cost = web.estimate_cost(bat.taskId, verbose=False)
    tot_cost += cost
    print(f"Maximum FlexCredit cost for {sim_name} = {cost:.2f}")
print(f"Maximum FlexCredit cost for batch = {tot_cost:.2f}")

13:33:52 PDT Created task 'linear' with task_id
             'fdve-4df5f2ac-f4db-4730-bf85-f72a3bbb3bb2v1' and task_type 'FDTD'.
             Created task 'quadratic' with task_id
             'fdve-bad852b2-87cd-42ef-aabe-03ef55d1926dv1' and task_type 'FDTD'.
13:33:53 PDT Created task 'exponential' with task_id
             'fdve-ecb0e9ab-17b8-498e-bed9-30fd185706f0v1' and task_type 'FDTD'.
Maximum FlexCredit cost for linear = 0.12
Maximum FlexCredit cost for quadratic = 0.12
Maximum FlexCredit cost for exponential = 0.12
Maximum FlexCredit cost for batch = 0.36

Run the batch and store all of the data in the data_taper dir.

[11]:
batch_taper = batch.run(path_dir="data/data_taper")

13:34:10 PDT Started working on Batch.
13:34:21 PDT Maximum FlexCredit cost: 0.356 for the whole batch.
             Use 'Batch.real_cost()' to get the billed FlexCredit cost after the
             Batch has completed.
13:34:58 PDT Batch complete.

Inverted taper results#

Now we will see the field intensity and the modal coupling efficiency for the 3 inverted tapers.

[12]:
sim_tap_results = {tl: sn for tl, sn in batch_taper.items()}

             loading simulation from
             data/data_taper/fdve-4df5f2ac-f4db-4730-bf85-f72a3bbb3bb2v1.hdf5
13:34:59 PDT loading simulation from
             data/data_taper/fdve-bad852b2-87cd-42ef-aabe-03ef55d1926dv1.hdf5
13:35:00 PDT loading simulation from
             data/data_taper/fdve-ecb0e9ab-17b8-498e-bed9-30fd185706f0v1.hdf5
[13]:
fig, axs = plt.subplots(3, 1, tight_layout=True, figsize=(10, 8))
for tap_n, sim_data, ax in zip(sim_tap_results.keys(), sim_tap_results.values(), axs):
    sim_data.plot_field("field_xy", "Ey", f=freq_c, val="abs", ax=ax)
    ax.set_title(tap_n)
    ax.set_aspect("auto")  # Used to better visualize the shapes.
plt.show()

../_images/notebooks_EdgeCoupler_25_0.png

Bellow, we have the \(|E_{y}|\) field recorded by the monitors positioned at the taper tip and taper output of the quadratic structure. We can clearly visualize the Gaussian profile corresponding to the fields launched by a lensed fiber at the taper tip. In contrast, the output field corresponds to the transverse electric polarization of the fundamental output waveguide mode.

[14]:
fig, (ax1, ax2) = plt.subplots(1, 2, tight_layout=True, figsize=(8, 4))
sim_tap_results["quadratic"].plot_field(
    "field_input", "Ey", f=freq_c, val="abs", ax=ax1
)
ax1.set_title("|Ey| at taper tip")
ax1.collections[-1].set_clim(0, 20)
ax1.set_xlim(-2, 2)
ax1.set_ylim(-2, 2)

sim_tap_results["quadratic"].plot_field(
    "field_output", "Ey", f=freq_c, val="abs", ax=ax2
)
ax2.set_title("|Ey| at taper output")
ax2.collections[-1].set_clim(0, 50)
ax2.set_xlim(-2, 2)
ax2.set_ylim(-2, 2)
plt.show()

../_images/notebooks_EdgeCoupler_27_0.png

The modal coupling efficiencies shown below reveal that the quadratic taper profile has the lowest insertion losses.

[15]:
power_tap_25 = []
for sim_data in sim_tap_results.values():
    mode_amps = sim_data["mode_0"]
    coeffs_f = mode_amps.amps.sel(direction="+")
    power_0 = np.abs(coeffs_f.sel(mode_index=0)) ** 2
    power_0_db = 10 * np.log10(power_0)
    power_tap_25.append(power_0_db)

plot_colors = ("black", "red", "blue")
fig, ax1 = plt.subplots(1, figsize=(8, 4))
for data, color, label in zip(power_tap_25, plot_colors, sim_tap_results.keys()):
    ax1.plot(
        wl_range,
        data,
        color=color,
        linestyle="solid",
        linewidth=1.0,
        label=label,
    )
ax1.set_xlim([wl_range[0], wl_range[-1]])
ax1.set_xlabel("Wavelength (um)")
ax1.set_ylabel("Power (dB)")
ax1.set_title("Coupling Efficiency")
plt.legend()
plt.show()

../_images/notebooks_EdgeCoupler_29_0.png

Inverted Taper Length Sweep#

Before, we simulated the inverted tapers with a short length of only 25 \(\mu m\) to improve simulation time. Nevertheless, the designs can be improved by making the inverted taper larger in length so that its adiabaticity is increased. In the following sections, we increase the taper length to 50 \(\mu m\) and 100 \(\mu m\) (other values can be added to the sweep at a higher simulation cost).

[16]:
# Taper lengths on the sweep.
taper_sweep = [50, 100]
sim_sweep = [get_simulations(tap_length=t_l) for t_l in taper_sweep]

# Make a dictionary of {task name : simulation}
sims = {
    f"sim_{tap_l}_{sim_name}": sim
    for sim_l, tap_l in zip(sim_sweep, taper_sweep)
    for sim, sim_name in zip(sim_l.values(), sim_l.keys())
}
# Initialize a batch and run them all
batch = web.Batch(simulations=sims, verbose=True)

# Get the estimated simulation cost.
tot_cost = 0
for bat in batch.get_info().values():
    sim_name = bat.taskName
    cost = web.estimate_cost(bat.taskId, verbose=False)
    tot_cost += cost
    print(f"Maximum FlexCredit cost for {sim_name} = {cost:.2f}")
print(f"Maximum FlexCredit cost for batch = {tot_cost:.2f}")

13:35:04 PDT WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
             WARNING: Default value for the field monitor '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.                                    
             Created task 'sim_50_linear' with task_id
             'fdve-93460a36-f53c-44bb-af00-62758a847217v1' and task_type 'FDTD'.
13:35:05 PDT Created task 'sim_50_quadratic' with task_id
             'fdve-2ad824f0-912a-4d65-8f7d-94cd382870c1v1' and task_type 'FDTD'.
             Created task 'sim_50_exponential' with task_id
             'fdve-33777f58-6d7c-4b17-8e1d-d8aaa87419f2v1' and task_type 'FDTD'.
13:35:06 PDT Created task 'sim_100_linear' with task_id
             'fdve-59ab670e-e863-4c16-aaaf-68c7e5ffd5d7v1' and task_type 'FDTD'.
             Created task 'sim_100_quadratic' with task_id
             'fdve-ed1f7ff1-c31f-4ca1-8049-1a5bb827fe2av1' and task_type 'FDTD'.
13:35:07 PDT Created task 'sim_100_exponential' with task_id
             'fdve-46c356b0-f5a1-4cfa-9eb0-f5f697b7a27dv1' and task_type 'FDTD'.
Maximum FlexCredit cost for sim_50_linear = 0.22
Maximum FlexCredit cost for sim_50_quadratic = 0.22
Maximum FlexCredit cost for sim_50_exponential = 0.22
Maximum FlexCredit cost for sim_100_linear = 0.42
Maximum FlexCredit cost for sim_100_quadratic = 0.42
Maximum FlexCredit cost for sim_100_exponential = 0.42
Maximum FlexCredit cost for batch = 1.94

Run the batch and store all of the data in the data_sweep dir.

[17]:
batch_sweep = batch.run(path_dir="data/data_sweep")

13:35:40 PDT Started working on Batch.
13:36:19 PDT Maximum FlexCredit cost: 1.937 for the whole batch.
             Use 'Batch.real_cost()' to get the billed FlexCredit cost after the
             Batch has completed.
13:36:54 PDT Batch complete.

Print the simulation log corresponding to the longest inverted taper to inspection.

[18]:
print(batch_sweep["sim_100_quadratic"].log)

13:36:57 PDT loading simulation from
             data/data_sweep/fdve-ed1f7ff1-c31f-4ca1-8049-1a5bb827fe2av1.hdf5
Simulation domain Nx, Ny, Nz: [3477, 126, 136]
Applied symmetries: (0, -1, 1)
Number of computational grid points: 1.5829e+07.
Using subpixel averaging: True
Number of time steps: 8.9104e+04
Automatic shutoff factor: 1.00e-05
Time step (s): 5.3908e-17


Compute source modes time (s):     0.4150
Compute monitor modes time (s):    8.6281
Rest of setup time (s):            9.8532

Running solver for 89104 time steps...
- Time step   2364 / time 1.27e-13s (  2 % done), field decay: 1.00e+00
- Time step   3564 / time 1.92e-13s (  4 % done), field decay: 1.00e+00
- Time step   7128 / time 3.84e-13s (  8 % done), field decay: 7.76e-01
- Time step  10692 / time 5.76e-13s ( 12 % done), field decay: 5.40e-01
- Time step  14256 / time 7.69e-13s ( 16 % done), field decay: 4.18e-01
- Time step  17820 / time 9.61e-13s ( 20 % done), field decay: 3.55e-01
- Time step  21384 / time 1.15e-12s ( 24 % done), field decay: 3.02e-01
- Time step  24949 / time 1.34e-12s ( 28 % done), field decay: 1.31e-05
- Time step  28513 / time 1.54e-12s ( 32 % done), field decay: 1.39e-07
Field decay smaller than shutoff factor, exiting solver.

Solver time (s):                   30.2616
Data write time (s):               0.0695

Parameter sweep results#

Now, we will get the batch results and plot the modal coupling efficiency at the central wavelength as a function of the inverted taper length.

[19]:
coup_eff_tl = np.zeros((len(tap_names), len(taper_sweep) + 1))
for i in range(len(tap_names)):
    coup_eff_tl[i, 0] = power_tap_25[i].sel(f=freq_c)

sims = {
    f"sim_{tap_l}_{sim_name}": sim
    for sim_l, tap_l in zip(sim_sweep, taper_sweep)
    for sim, sim_name in zip(sim_l.values(), sim_l.keys())
}

for j, tl in enumerate(taper_sweep):
    for i, sn in enumerate(tap_names):
        data_sim = batch_sweep[f"sim_{tl}_{sn}"]
        mode_amps = data_sim["mode_0"]
        coeffs_f = mode_amps.amps.sel(direction="+", f=freq_c)
        power_0 = np.abs(coeffs_f.sel(mode_index=0)) ** 2
        power_0_db = 10 * np.log10(power_0)
        coup_eff_tl[i, j + 1] = power_0_db

13:36:58 PDT loading simulation from
             data/data_sweep/fdve-93460a36-f53c-44bb-af00-62758a847217v1.hdf5
13:36:59 PDT loading simulation from
             data/data_sweep/fdve-2ad824f0-912a-4d65-8f7d-94cd382870c1v1.hdf5
13:37:01 PDT loading simulation from
             data/data_sweep/fdve-33777f58-6d7c-4b17-8e1d-d8aaa87419f2v1.hdf5
13:37:02 PDT loading simulation from
             data/data_sweep/fdve-59ab670e-e863-4c16-aaaf-68c7e5ffd5d7v1.hdf5
             loading simulation from
             data/data_sweep/fdve-ed1f7ff1-c31f-4ca1-8049-1a5bb827fe2av1.hdf5
13:37:03 PDT loading simulation from
             data/data_sweep/fdve-46c356b0-f5a1-4cfa-9eb0-f5f697b7a27dv1.hdf5

As shown below, the longer the taper, the higher the modal coupling efficiency, as is expected due to the increased adiabaticity. However, the quadratic shape profile seems more advantageous for this case as it achieved a steady state coupling loss of \(\lt\) 0.6 dB for a 50 \(\mu m\) taper length. In contrast, the exponential shape needs double this value, and an even longer taper length would be necessary for the linear profile to obtain similar coupling loss. Detailed discussions about the effects of taper profile over the coupling efficiency, bandwidth, footprint size, or taper/fiber misalignments can be found in the literature [1,2,3].

[20]:
fig, ax = plt.subplots(1, figsize=(6, 3))
lin_style = ["ok", "sr", "hb"]
taper_lengths = [taper_l] + taper_sweep[:]
for ce, ls, lab in zip(coup_eff_tl, lin_style, tap_names):
    ax.plot(taper_lengths, ce, ls, label=lab)
ax.set_xlim([taper_lengths[0] - 2, taper_lengths[-1] + 2])
ax.set_xlabel("Taper Length ($\mu m$)")
ax.set_ylabel("Coupling Efficiency (dB)")
plt.legend()
plt.show()

../_images/notebooks_EdgeCoupler_39_0.png

Next, we will visualize the field distribution and the wavelength-dependent coupling coefficients for these longer inverted taper structures.

[21]:
final_sim_data = [
    batch_sweep["sim_100_linear"],
    batch_sweep["sim_50_quadratic"],
    batch_sweep["sim_100_exponential"],
]

fig, axs = plt.subplots(3, 1, tight_layout=True, figsize=(10, 8))
for sim_data, tap_n, ax in zip(final_sim_data, tap_names, axs):
    sim_data.plot_field("field_xy", "Ey", f=freq_c, val="abs", ax=ax)
    ax.set_title(tap_n)
    ax.set_aspect("auto")  # Used to better visualize the shapes.
plt.show()

             loading simulation from
             data/data_sweep/fdve-59ab670e-e863-4c16-aaaf-68c7e5ffd5d7v1.hdf5
13:37:04 PDT loading simulation from
             data/data_sweep/fdve-2ad824f0-912a-4d65-8f7d-94cd382870c1v1.hdf5
             loading simulation from
             data/data_sweep/fdve-46c356b0-f5a1-4cfa-9eb0-f5f697b7a27dv1.hdf5
../_images/notebooks_EdgeCoupler_41_3.png
[22]:
plot_colors = ("black", "red", "blue")
fig, ax1 = plt.subplots(1, figsize=(8, 4))
for sim_data, color, label in zip(final_sim_data, plot_colors, tap_names):
    mode_amps = sim_data["mode_0"]
    coeffs_f = mode_amps.amps.sel(direction="+")
    power_0 = np.abs(coeffs_f.sel(mode_index=0)) ** 2
    power_0_db = 10 * np.log10(power_0)

    ax1.plot(
        wl_range,
        power_0_db,
        color=color,
        linestyle="solid",
        linewidth=1.0,
        label=label,
    )
ax1.set_xlim([wl_range[0], wl_range[-1]])
ax1.set_xlabel("Wavelength (um)")
ax1.set_ylabel("Power (dB)")
ax1.set_title("Coupling Efficiency")
plt.legend()
plt.show()

../_images/notebooks_EdgeCoupler_42_0.png