Trel725 / plasmon-meep

set of scripts for calculation of plasmon resonance/electric field enhancement on different structures
GNU General Public License v3.0
6 stars 4 forks source link

A specified frequency must be included in the frequency range #8

Closed yu-ting-py closed 2 years ago

yu-ting-py commented 3 years ago

Hi

I would like to see field enhancement at a specified frequency f=1/0.633

So I tune the frequency range:

cfreq = 1.5817
fwidth = 1.8

and corresponding frequencies

for i in range(nfreq):
 print(1 / flux_freqs[i])

wavlength [um] close to 1/cfreq:

0.6414010519621877
0.6377013525998186
0.6340440893231117
0.6304285361734508

after using the caculate2d_field_disk_cache.py file, the final freqs.h5 shows frequencies as above. How to get field enhancement exactly at frequency = 1/0.633?

Trel725 commented 3 years ago

Hi,

The easiest way is to interpolate the values. You could also try to increase the resolution in the simulation by increasing nfreq param: https://github.com/Trel725/plasmon-meep/blob/cfd371a15f05fe294829c32ad8d240e8ed142ba3/perform_fdtd.py#L57

although it will not be very different from the simple interpolation. If you want to measure the enhancement precisely at that frequency, you should rather use frequency domain solver or try to narrow the range of your source.

yu-ting-py commented 3 years ago

Hi

Thanks for your advice.

I also tried frequency-domain solver before. But it cannot converge somehow.

The example on MEEP tutorial by frequency-domain solver works on my computer. (https://nbviewer.jupyter.org/github/NanoComp/meep/blob/master/python/examples/solve-cw.ipynb)

I follow the steps from MEEP website and write a script of continuous-wave excitation by the frequency-domain solver.

For a quick try, I didn't compute the field with empty geometry.

Here is the code:

import meep as mp
import numpy as np
import matplotlib.pyplot as plt
from matplotlib import cm
from meep import materials

# radii
r = 50e-3

sizex = 240e-3
sizey = 240e-3

size_2d = (sizex, sizey)
# pixels / um
resolution = 800

pixel_size = 1/resolution
plane_th = 1*r
# original thickness = 0.2*sizex
pml_th = 0.040
fullx = sizex + 2*pml_th
fully = sizey + 2*pml_th
cell = mp.Vector3(fullx, fully, 0)

pml_layers = [mp.PML(thickness=pml_th, direction=mp.X), mp.PML(thickness=pml_th, direction=mp.Y)]

nonpml_vol = mp.Volume(mp.Vector3(), size=mp.Vector3(sizex, sizey))

wl = 0.633
cfreq = 1/wl

comp = mp.Ex

sources = [mp.Source(mp.ContinuousSource(frequency=cfreq, is_integrated=True),
                     size=mp.Vector3(fullx, 0, 0),
                     component=comp,
                     center=mp.Vector3(0, -0.5*sizey))]

geometry = [mp.Sphere(material=mp.Medium(epsilon=-11.685, D_conductivity=-2*np.pi*1.58*1.2672/11.685),
            center=mp.Vector3(0, 0),
            radius=r),
            ]
# Au permitivity at 0.633
#real part -11.685
#imaginary part 1.2672

sim = mp.Simulation(cell_size=cell,
                    boundary_layers=pml_layers,
                    geometry=geometry,
                    sources=sources,
                    resolution=resolution,
                    force_complex_fields=True,
                    split_chunks_evenly=False,)

f = plt.figure(dpi=150)
sim.plot2D(ax = f.gca())
plt.show()

sim.init_sim()
sim.solve_cw(1e-8, 10000, 10)

ex_norm =  np.real(sim.get_array(vol=nonpml_vol, component=mp.Ex))

ex_en = np.abs(ex_norm)**2

plt.figure()
plt.imshow(ex_en, cmap=cm.jet, interpolation='none')
plt.colorbar()

plt.show()

It appears "could not converge" at the end.

Could you please tell me the reason?

Trel725 commented 3 years ago

A quote from the MEEP docs:

The frequency-domain solver supports arbitrary geometries, PML, boundary conditions, symmetries, parallelism, conductors, and arbitrary nondispersive materials. Lorentz-Drude dispersive materials are not currently supported in the frequency-domain solver, but since you are solving at a known fixed frequency rather than timestepping, you should be able to pick conductivities etcetera in order to obtain any desired complex ε and μ at that frequency.

So, you should follow their suggestion and build nondispersive material. To be honest, I do not know how to do that.

Good luck!