Calculation of effective mode area or nonlinear coefficent

[elizaleung830]

Does it support the calculation of effective mode area or nonlinear coefficent of waveguide? I find the equation effective area in the maxwell.md but i can't find it in the example.
Thanks in advance.

[lucasgrjn]

Hi, with the definition of the effective mode area $A_{\text{eff}}=\frac{\left( \int \left| \vec{\mathcal{E}} \right|^2 \mathrm{d}A \right)^2}{\int \left| \vec{\mathcal{E}} \right|^4 \mathrm{d}A}$, it should be possible to calculate it with the mode found and the integration of $\left|\vec{\mathcal{E}} \right|^2$ and $\left|\vec{\mathcal{E}} \right|^4$ over the surface!

Have you an example or a piece of code you are interested in? :)

[HelgeGehring]

Sounds good! Would be best to have some paper which reports a value which we can reproduce :)
@elizaleung830 do you know a paper which reports all the design values and a effective mode area?

[elizaleung830]

Hi I am trying to replicate the result of figure 1b and c from this paper

So far I have code the following scipt, which define the waveguide geometry and plot the n_eff vs width diagram. (I have tried to implement the method of calculate_effetive_mode myself but it doesn't work). It would be great if it can be implemented, Thanks! :)

import matplotlib.pyplot as plt
import shapely.affinity

from femwell.mesh import mesh_from_OrderedDict
from collections import OrderedDict

from skfem import Basis, ElementTriP0
from skfem.io.meshio import from_meshio

from femwell.maxwell.waveguide import compute_modes
wavelength = 1.55
capital_w = 1.4
width = 0.05 #sweep
capital_h = 0.3
h = 0.5
t = 0.1
n_silicon = 3.48
n_air = 1
n_nc = 1.6
n_silica = 1.45
w_list = [ x*1e-3 for x in range(50,700,10)]
neff_list = []
aeff_list = []
for width in w_list:
    air = shapely.geometry.box(-capital_w / 2, 0, +capital_w / 2, h + capital_h)
    nc = shapely.geometry.box(-width / 2, capital_h + (h - t) / 2, +width / 2, capital_h + (h - t) / 2 + t)
    silica = shapely.geometry.box(-capital_w / 2, 0, +capital_w / 2, capital_h)
    silicon = shapely.geometry.box(-width / 2, capital_h, +width / 2, capital_h + h)

    polygons = OrderedDict(
        core=nc,
        silicon=silicon,
        silica=silica,
        air=air
    )

    resolutions = dict(core={"resolution": 0.03, "distance": 0.5})

    mesh = from_meshio(mesh_from_OrderedDict(polygons, resolutions, default_resolution_max=10))

    basis0 = Basis(mesh, ElementTriP0())
    epsilon = basis0.zeros()
    for subdomain, n in {"core": n_nc, "silicon": n_silicon, "air": n_air, "silica": n_silica}.items():
        epsilon[basis0.get_dofs(elements=subdomain)] = n ** 2
    modes = compute_modes(basis0, epsilon, wavelength=wavelength, num_modes=1, order=2)

    for mode in modes:
        print(f"Effective refractive index: {mode.n_eff:.4f}")
        print(f"Effective mode area: {mode.calculate_effective_area():.4f}")
        neff_list.append(mode.n_eff)
        aeff_list.append(mode.calculate_effective_area())
plt.subplot(211)
plt.plot(w_list,aeff_list)
plt.ylabel('aeff')

plt.subplot(212)
plt.plot(w_list,neff_list)
plt.ylabel('neff')
plt.show()

[HelgeGehring]

Yay, looks great! :)

Could you make a pull-request that adds your example to the docs?
Like https://github.com/HelgeGehring/femwell/blob/main/docs/photonics/examples/waveguide_modes.py you'd need to add a header and some comments to split it into cells.
Then we can have a look on how to implement the calculate_effective_area().

Using https://automeris.io/WebPlotDigitizer/ we can also extract their data from the paper to have it in the docs too as a benchmark :)

All those definitions for the effective area area a bit confusing and https://opg.optica.org/ol/fulltext.cfm?uri=ol-39-15-4396&id=296216 gives even another formula for the effective area :/
Probably it's not a problem to add multiple ones and have the formulas in the docs? 🤔 as they all descrive slightly different effects? 🤔 @elizaleung830 what do you think?

[elizaleung830]

Yes sure! I will make the pull-request as soon as I add the documentation in the code.

For the effective areaa, it would be great to nhave mutiple formula for comparsion, but the main one that is used for my stimulation will probably be the formula mentioned by @lucasgrjn.

Hi, with the definition of the effective mode area Aeff=(∫|E→|2dA)2∫|E→|4dA, it should be possible to calculate it with the mode found and the integration of |E→|2 and |E→|4 over the surface!

Have you an example or a piece of code you are interested in? :)

[HelgeGehring]

👍 sounds good!

Hmm, that formula seems to me only adequate for fibers as it doesn't take into account what portion of the material is nonlinear 🤔 I'd be careful using that for nonlinear effects in waveguides

Whatever, having multiple formulas and comparing them to literature will probably be the best :)

[lucasgrjn]

The example is pretty good! We can also exploit the other data of the original paper?

I see the following that may add value to the doc:

• Add example for the eq.1 & 2 as pointed by @elizaleung830 in Implement effective_area + docs example #107
• Add data points for a waveguide width w < 250
• Reproduce Fig.2
• Reproduce Fig.3

Any thoughts @HelgeGehring @elizaleung830 ?

[HelgeGehring]

Sounds like a great plan!
For the widths below 250nm we probably need to figure out how they selected the modes, it seems not to be the TM-mode with the highest refractive index there :/
I'm also wondering how large they chose the surrounding air - for small, weakly guiding structures (-> super large effective areas) it has probably a large effect on the estimated parameters.

[elizaleung830]

Sounds great to me as well! I am not very familiar with implementing the equation, but i more than happy to help with reproducing figure 2 and 3 or data from other paper!

[lucasgrjn]

@elizaleung830 Perfect!

I will take a look in the equations while you investigate the figures then :)

[lucasgrjn]

Hey @elizaleung830
I implement the calculations for the Poynting vector which should be useful to calculate the equations!
We are close to have the article's equations implemented :)

[HelgeGehring]

@lucasgrjn thanks! I just merged the PR 🥳

@elizaleung830 I've also added the integrals needed for (1) and (3) - could you add the benchmarks to the paper again? 😃 We should now have everything to reproduce the whole paper, right? :)

[HelgeGehring]

@elizaleung830 Thanks, amazing PR!

I think what's left now is to understand better whats happening for widths < 250nm, my guess would be that two modes are crossing around there and we stick to the wrong one, that they have a different amount of air around 🤔 for the lowest confinement modes that would matter most, or that they just use different boundary conditions.

Could you have a look at that?
Easiest would probably be to reproduce the other mode images in Fig1, that should tell us directly if we have the right modes :)
(we could probably also use there a loop to plot all the modes)