Leaky mode

Contents

Leaky mode#

Reproducing one example of [1]

Hide code cell source
from collections import OrderedDict

import matplotlib.pyplot as plt
import numpy as np
from shapely import LineString, box
from skfem import Basis, ElementDG, ElementTriP1
from skfem.io import from_meshio

from femwell.mesh import mesh_from_OrderedDict
from femwell.mode_solver_2d_periodic import plot_periodic, solve_periodic
height = 5.76 / 2 + 5
a = 0.010
b = 0.78
c = 0.2
slab = 0.920 + 5
pml = 3

wavelength = 1
k0 = 2 * np.pi / wavelength

left = LineString([(0, y) for y in np.linspace(-height, height, 2)])
right = LineString([(a, y) for y in np.linspace(-height, height, 2)])
top = LineString([(x, height) for x in np.linspace(0, a, 2)])
bottom = LineString([(x, -height) for x in np.linspace(0, a, 2)])

background = box(0, -height, a, height)
structure = box(0, -b / 2, a, b / 2)
structure1 = box(0, height - slab, a, height)
structure2 = box(0, -height + slab, a, -height)

resolutions = {
    "structure": {"resolution": 0.1, "distance": 0.1},
    "hole": {"resolution": 0.1, "distance": 0.1},
}

mesh = from_meshio(
    mesh_from_OrderedDict(
        OrderedDict(
            left=left,
            right=right,
            top=top,
            bottom=bottom,
            structure=structure,
            structure1=structure1,
            structure2=structure2,
            background=background,
        ),
        resolutions=resolutions,
        filename="mesh.msh",
        default_resolution_max=0.05,
        periodic_lines=[("left", "right")],
    )
)
basis_vec = Basis(mesh, ElementTriP1() * ElementTriP1())
basis_epsilon_r = basis_vec.with_element(ElementDG(ElementTriP1()))

epsilon_r = basis_epsilon_r.zeros(dtype=np.complex64) + 1.45
epsilon_r[basis_epsilon_r.get_dofs(elements="background")] = 1.39
epsilon_r **= 2
basis_epsilon_r.plot(np.real(epsilon_r), aspect=a, colorbar=True).show()

epsilon_r += basis_epsilon_r.project(
    lambda x: (0.5j) * (np.clip(np.abs(x[1]) - height + pml, 0, np.inf) / pml) ** 2,
    dtype=np.complex64,
)
basis_epsilon_r.plot(np.imag(epsilon_r), aspect=a, colorbar=True).show()

ks, basis_phi, phis = solve_periodic(basis_epsilon_r, epsilon_r, k0)

for i, k in enumerate(ks):
    if 0 < np.imag(k) < 0.002:
        fig, ax = plt.subplots(1, 1, figsize=(10, 4))
        plt.title(f"{k}")
        plot_periodic(k, a, basis_phi, phis[..., i], 100, ax)
        plt.show()
../../_images/ff7b98871c03d6b8021e99e57d73365f38f675df27476cb0866b28bd7373722f.png ../../_images/a3c5beffb2e20ff894108a8ee4497aad197eac8acf268dfbd456b937522e1a1f.png
/home/runner/miniconda3/lib/python3.12/site-packages/scipy/sparse/linalg/_dsolve/linsolve.py:603: SparseEfficiencyWarning: splu converted its input to CSC format
  return splu(A).solve
/home/runner/miniconda3/lib/python3.12/site-packages/scipy/sparse/linalg/_matfuncs.py:76: SparseEfficiencyWarning: spsolve is more efficient when sparse b is in the CSC matrix format
  Ainv = spsolve(A, I)
../../_images/8ae6d290c786f7e6a40ac4a955dd409761837eee292c75c7b7fef6f8dfb9a7a9.png

Bibliography#

[1]

Jonathan Hu and Curtis R. Menyuk. Understanding leaky modes: slab waveguide revisited. Advances in Optics and Photonics, 1(1):58, January 2009. URL: https://doi.org/10.1364/aop.1.000058, doi:10.1364/aop.1.000058.