Open stevengj opened 4 years ago
In the following example of a ring resonator for a mode with Q of ~106 (computed using harminv
via a narrowband source with an error of ~10-10), the mode frequency/Q pair computed by solve_eigfreq
(with a guessfreq
from harminv
) has significantly worse accuracy (factor of 104!) and requires a substantially longer runtime than harminv
(1483.8517
s vs. 182.8247
s using two cores). Increasing the solve_eigfreq
convergence parameters maxiters
and cwmaxiters
improves accuracy at the expense of longer runtime.
From these results, it seems that the eigensolver is impractical for these applications.
harminv result
harminv0:, frequency, imag. freq., Q, |amp|, amplitude, error
harminv0:, 0.25763966078694694, -5.938000849363818e-07, 216941.41456257098, 0.06745820629360579, -0.06655029897341563+0.011030290245486258i, 2.0652221851419098e-10+0.0i
eigensolver result
Eigensolver step 30: 80 CG iters, freq = 0.256138-0.000940994i (change = 0.00169493+0.00010582i).
-- CONVERGENCE FAILURE (1) in solve_cw!
eigfreq:, 0.25613766961518836, 136.0995780261287
import meep as mp
n = 3.4 # index of waveguide
w = 1 # width of waveguide
r = 1 # inner radius of ring
pad = 4 # padding between waveguide and edge of PML
dpml = 2 # thickness of PML
sxy = 2*(r+w+pad+dpml) # cell size
c1 = mp.Cylinder(radius=r+w, material=mp.Medium(index=n))
c2 = mp.Cylinder(radius=r)
fcen = 0.257 # source center frequency
src = mp.Source(mp.ContinuousSource(fcen),
component=mp.Ez,
center=mp.Vector3(r+0.1))
sim = mp.Simulation(cell_size=mp.Vector3(sxy, sxy),
geometry=[c1, c2],
sources=[src],
resolution=40,
force_complex_fields=True,
symmetries=[mp.Mirror(mp.Y)],
boundary_layers=[mp.PML(dpml)])
sim.init_sim()
eigfreq = sim.solve_eigfreq(tol=1e-7, maxiters=30, guessfreq=fcen, cwtol=1e-6, cwmaxiters=80, L=10)
if mp.am_master():
print("eigfreq:, {}, {}".format(eigfreq.real,eigfreq.real/(-2*eigfreq.imag)))
By comparison, for a mode with smaller Q of ~104 the eigensolver produces results with slightly better accuracy but still far from adequate.
harminv result
harminv0:, frequency, imag. freq., Q, |amp|, amplitude, error
harminv0:, 0.20345453628737467, -1.2503433555924088e-05, 8135.946633274128, 0.040612299864003616, -0.02673202846745665-0.030573805034029296i, 6.3668396774220965e-12+0.0i
eigensolver result
Eigensolver step 30: 80 CG iters, freq = 0.178157-8.45776e-05i (change = 0.0159243+2.18228e-05i).
-- CONVERGENCE FAILURE (1) in solve_cw!
eigfreq:, 0.1781570997298217, 1053.2162229720261
Now that #1158 is merged, it might be nice to use
solve_eigfreq
in a tutorial.Typically, to find resonances we do multiple simulations: first we use
harminv
with a broadband pulse, to identify the resonant frequencies. Then we do additional calculations with narrow-band pulses around each resonant frequency in order to see individual modes and also to get more accurateharminv
data.Instead, you can now replace the second step by calling
solve_eigfreq
with aguessfreq
given by the broadbandharminv
result. This might be faster than a very narrowband calculation, I'm not sure yet (will probably depend on the problem). But it can probably give a more accurate result, especially for the Q of the mode (since high-Q calculations correspond to finding the complex ω to many decimal places).(For example, suppose you have a resonator with a Q=10⁵, then if you want 3 digits of accuracy in Q you need to compute the complex ω to about 8 digits. This might require a very long simulation for Harminv, if it is possible at all.)