2D Anisotropic Grating#

Example of diffraction grating with trapezoidal ridges made from an anisotropic material.

from collections import OrderedDict

import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
from mpl_toolkits.axes_grid1 import make_axes_locatable

import gyptis as gy

We will study this benchmark and compare with results given in [PopovGratingBook].

fig, ax = plt.subplots(3, 2, figsize=(3.5, 5.5))


lambda0 = 633
period = 600

width_bottom, width_top = 500, 300
height = 600
eps_sub = 2.25
eps_rod = np.array([[2.592, 0.251, 0], [0.251, 2.592, 0], [0, 0, 2.829]])

pmesh = 10

thicknesses = OrderedDict(
    {
        "pml_bottom": 1 * lambda0,
        "substrate": 2 * lambda0,
        "groove": height * 1.5,
        "superstrate": 2 * lambda0,
        "pml_top": 1 * lambda0,
    }
)

mesh_param = dict(
    {
        "pml_bottom": pmesh * eps_sub**0.5,
        "substrate": pmesh * eps_sub**0.5,
        "groove": pmesh,
        "rod": pmesh * np.max(eps_rod) ** 0.5,
        "superstrate": pmesh,
        "pml_top": pmesh,
    }
)


geom = gy.Layered(2, period, thicknesses)
groove = geom.layers["groove"]
substrate = geom.layers["substrate"]
y0 = geom.y_position["groove"]
P = [geom.add_point(-width_bottom / 2, y0, 0)]
P.append(geom.add_point(width_bottom / 2, y0, 0))
P.append(geom.add_point(width_top / 2, y0 + height, 0))
P.append(geom.add_point(-width_top / 2, y0 + height, 0))
L = [
    geom.add_line(P[0], P[1]),
    geom.add_line(P[1], P[2]),
    geom.add_line(P[2], P[3]),
    geom.add_line(P[3], P[0]),
]
cl = geom.add_curve_loop(L)
rod = geom.add_plane_surface(geom.dimtag(cl, 1)[0])
substrate, groove, rod = geom.cut([substrate, groove], rod)
geom.add_physical(rod, "rod")
geom.add_physical(groove, "groove")
geom.add_physical(substrate, "substrate")
mesh_size = {d: lambda0 / param for d, param in mesh_param.items()}
geom.set_mesh_size(mesh_size)

geom.build(
    interactive=False,
    generate_mesh=True,
    write_mesh=True,
    read_info=True,
)
all_domains = geom.subdomains["surfaces"]
domains = [k for k in all_domains.keys() if k not in ["pml_bottom", "pml_top"]]

epsilon = {d: 1 for d in domains}
mu = {d: 1 for d in domains}

epsilon["substrate"] = eps_sub
epsilon["rod"] = eps_rod


nper = 8

for jangle, angle in enumerate([0, -20, -40]):
    angle_degree = angle * np.pi / 180

    pw = gy.PlaneWave(lambda0, angle_degree, dim=2)
    grating_TM = gy.Grating(geom, epsilon, mu, source=pw, polarization="TM", degree=2)
    grating_TM.solve()
    effs_TM = grating_TM.diffraction_efficiencies(2, orders=True)

    print(f"angle = {angle}, TM polarization")
    print("--------------------------------")
    print("R: ", effs_TM["R"])
    print("T: ", effs_TM["T"])

    ylim = geom.y_position["substrate"], geom.y_position["pml_top"]
    d = grating_TM.period
    vmin_TM, vmax_TM = -1.5, 1.7
    plt.sca(ax[jangle][0])
    per_plots, cb = grating_TM.plot_field(nper=nper)
    cb.remove()
    scatt_lines, layers_lines = grating_TM.plot_geometry(nper=nper, c="k")
    [layers_lines[i].remove() for i in [0, 1, 3, 4]]
    plt.ylim(ylim)
    plt.xlim(-d / 2, nper * d - d / 2)
    plt.axis("off")

    # TE

    grating_TE = gy.Grating(geom, epsilon, mu, source=pw, polarization="TE", degree=2)

    grating_TE.solve()
    effs_TE = grating_TE.diffraction_efficiencies(2, orders=True)

    H = grating_TE.solution["total"]
    print(f"angle = {angle}, TE polarization")
    print("--------------------------------")
    print("R: ", effs_TE["R"])
    print("T: ", effs_TE["T"])

    vmin_TE, vmax_TE = -2.5, 2.5
    plt.sca(ax[jangle][1])
    per_plots, cb = grating_TE.plot_field(nper=nper)
    cb.remove()
    scatt_lines, layers_lines = grating_TE.plot_geometry(nper=nper, c="k")
    [layers_lines[i].remove() for i in [0, 1, 3, 4]]
    plt.ylim(ylim)
    plt.xlim(-d / 2, nper * d - d / 2)
    plt.axis("off")

    ax[jangle][0].set_title(rf"$\theta = {angle}\degree$")
    ax[jangle][1].set_title(rf"$\theta = {angle}\degree$")


divider = make_axes_locatable(ax[0, 0])
cax = divider.new_vertical(size="5%", pad=0.5)
fig.add_axes(cax)
mTM = plt.cm.ScalarMappable(cmap="RdBu")
mTM.set_clim(vmin_TM, vmax_TM)

cbarTM = fig.colorbar(mTM, cax=cax, orientation="horizontal")
cax.set_title(r"${\rm Re}\, E_z$ (TM)")

divider = make_axes_locatable(ax[0, 1])
cax = divider.new_vertical(size="5%", pad=0.5)

mTE = plt.cm.ScalarMappable(cmap="RdBu")
mTE.set_clim(vmin_TE, vmax_TE)
fig.add_axes(cax)
cbarTE = fig.colorbar(mTE, cax=cax, orientation="horizontal")
cax.set_title(r"${\rm Re}\, H_z$ (TE)")

plt.tight_layout()
plt.subplots_adjust(wspace=-0.1, hspace=-0.3)
$\theta = 0\degree$, $\theta = 0\degree$, $\theta = -20\degree$, $\theta = -20\degree$, $\theta = -40\degree$, $\theta = -40\degree$, ${\rm Re}\, E_z$ (TM), ${\rm Re}\, H_z$ (TE)

Out:

angle = 0, TM polarization
--------------------------------
R:  [-4.800821302812278e-28, 1.302873843768947e-19, 0.008474283999092601, 1.3025887520542717e-19, -1.6323261920830313e-26]
T:  [2.57725693279246e-10, 0.20306992534769838, 0.5844347346981851, 0.20305413403395936, 2.5852272070126456e-10]
angle = 0, TE polarization
--------------------------------
R:  [-2.0894485013181953e-28, 5.522646067501839e-20, 0.01469264687176904, 7.194151935609829e-20, -6.306004034057696e-27]
T:  [1.3505410538920432e-11, 0.3225213754874706, 0.537712888377713, 0.12483912903990863, 9.772441142445927e-13]
angle = -20, TM polarization
--------------------------------
R:  [-7.437393998587917e-26, 0.00440145672406441, 0.015649696328782905, -1.0341928550174941e-23, -3.23075871794277e-26]
T:  [9.603605805621784e-14, 0.3994569273527532, 0.5750482867152704, 0.004633836044539676, 4.424726797287866e-08]
angle = -20, TE polarization
--------------------------------
R:  [-1.0887104429393313e-26, 0.005361690563464985, 0.01116230432450693, -5.129008063814091e-25, -2.2905800617653914e-26]
T:  [3.1591499562047014e-14, 0.5384830427561081, 0.4444623156789897, 0.00036518485135289064, -1.5422389655858709e-07]
angle = -40, TM polarization
--------------------------------
R:  [-1.0927670940551712e-23, 0.0025001016942677455, 0.05840998905701427, -3.0392602306671384e-24, -1.2123012245169451e-25]
T:  [0.026215511119240517, 0.4202511534997685, 0.49205443305936497, 3.5672226090250254e-14, -1.749884029070645e-05]
angle = -40, TE polarization
--------------------------------
R:  [-1.0560682830849532e-24, 0.0050336191740690565, 0.007695031429749402, -1.4753287186688208e-24, -1.7603721674863266e-26]
T:  [0.011869432010323797, 0.43403306218555954, 0.5411566793327887, 1.717031778233771e-14, -0.00026604411722746667]

Total running time of the script: ( 0 minutes 22.711 seconds)

Estimated memory usage: 55 MB

Gallery generated by Sphinx-Gallery