NanoComp / meep

free finite-difference time-domain (FDTD) software for electromagnetic simulations
GNU General Public License v2.0
1.16k stars 594 forks source link

Add support for `Absorber` in cylindrical coordinates #2824

Open oskooi opened 2 months ago

oskooi commented 2 months ago

The Absorber feature, which is a drop-in replacement for PML for cases involving lossy metals, is not currently supported in cylindrical coordinates. This would be a useful feature to add.

As a simple test case, a point-dipole emitter in vacuum in cylindrical coordinates should radiate the same power regardless of the boundary layers (PML or Absorber). (The radiated power can also be computed using the LDOS which is another way to validate the result obtained using the DFT flux monitors.) However, as shown below, this is not currently the case. For some reason, when the cell is terminated using an Absorber, the measured power is zero (!).

1. PML

dipole_power:, 5.942269

2. Absorber

dipole_power:, 0.000000

dipole_cyl

import matplotlib.pyplot as plt
import meep as mp

RESOLUTION_UM = 20

cell_size = mp.Vector3(5, 0, 5)

dipole_pos = mp.Vector3(1.5 / RESOLUTION_UM, 0, 0)
sources = [
    mp.Source(
        mp.GaussianSource(frequency=1.0, fwidth=0.1),
        component=mp.Er,
        center=dipole_pos,
    )
]

boundary_layers = [
    mp.Absorber(1.0, direction=mp.R),
    mp.Absorber(1.0, direction=mp.Z),
]

sim = mp.Simulation(
    resolution=RESOLUTION_UM,
    cell_size=cell_size,
    dimensions=mp.CYLINDRICAL,
    m=-1,
    sources=sources,
    boundary_layers=boundary_layers,
)

flux_mon = sim.add_flux(
    1.0,
    0,
    1,
    mp.FluxRegion(
        center=mp.Vector3(2.0, 0, 1.5),
        size=mp.Vector3(4.0, 0, 0),
    ),
    mp.FluxRegion(
        center=mp.Vector3(4.0, 0, 0),
        size=mp.Vector3(0, 0, 3.0),
    ),
    mp.FluxRegion(
        center=mp.Vector3(2.0, 0, -1.5),
        size=mp.Vector3(4.0, 0, 0),
    ),
)

fig, ax = plt.subplots()
sim.plot2D(ax=ax)
fig.savefig("dipole_cyl.png", dpi=150, bbox_inches="tight")

sim.run(
    until_after_sources=mp.stop_when_fields_decayed(
        25, mp.Er, dipole_pos, 1e-6
    )
)

dipole_power = mp.get_fluxes(flux_mon)[0]
print(f"dipole_power:, {dipole_power:.6f}")