MEEP modelling: Square Cavity – Part 1

Assuming that you have already installed MEEP, attempted some of the examples and followed the tutorials, let's try a simple problem of finding the frequencies and fields of the eigenmodes of a two-dimensional square of air/vacuum surrounded by perfect electric conductor (PEC) boundaries.

MEEP can use the Scheme programming language (a dialect of Lisp) as an interface to the nuts & bolts of the code. You essentially code the problem you want to solve by writing a program.

So let's get started. Open a new file called square.ctl in your favourite text editor then copy & paste the following snippets of code.

; Eigenmodes of a 2D square cavity

Anything after a colon is a remark. Always handy to remind yourself what the code does:

The first thing we need to do is define the geometry-lattice which defines the extent and dimensionality of our computational space. Here, we define a two-dimensional square with sides of length 1 in meep units (more on those later) in the \(x\)-\(y\) plane and with zero thickness i.e. no-size in the \(z\)-direction:

(set! geometry-lattice (make lattice (size 1 1 no-size)))

The next step – and sometimes the trickiest – is to define the geometry, the position, shape, size and material of structures within the computational space. For this simple case we just define a zero-thickness square of air/vacuum with the same dimensions:

(set! geometry (list
        (make block (center 0 0 0) (size 1 1 no-size) (material air))
))

Set the resolution, the number of cells per unit distance (a conservative value to keep simulation time low):

(set-param! resolution 20)

Define some parameters for the centre frequency and bandwidth of the point-source gaussian pulse using define-param:

(define-param fcen 1) ; pulse center frequency                            
(define-param df 2)  ; pulse width (in frequency)

Define the source; in this case a gaussian pulse with centre frequency fcen = 1.5 and bandwidth fwidth = 3.0, positioned at \(x\)=0.25, \(y\)=0.25 and exciting the \(E_x\) component of the electric field:

(set! sources (list
               (make source
                 (src (make gaussian-src (frequency fcen) (fwidth df)))
                 (component Ex) (center 0.33 0.33) )));

Run the FDTD simulation for the duration of the source pulse (exercise: work it out), then continue for a sufficient length of time. The following statement run-sources does just that. At the beginning it saves a permittivity file (which isn't much help for us at the moment), then after the source pulse, the after-sources command runs harminv which analyses the \(E_x\) time-domain data collected at point \(x\)=0.25, \(y\)=0.25 using harmonic inversion to estimate the resonant frequencies (and Q-factors) of any resonances. The command then runs the FDTD simulation for an additional 200 meep time units, which happen to be 200 periods of the pulse centre frequency fcen (more on those later).

(run-sources+ 200
    (after-sources (harminv Ex (vector3 0.033 0.033) fcen df)))

Once you've copy & pasted all of these sections into your text editor, save it as square.ctl or download it from here. Then get MEEP to process it using the command-line, something like meep square.ctl | tee square.out. The pipe (tee) echoes the output into a file.

You should (all being well) see some output text, showing that meep has executed the instructions in square.ctl and finished without errors. Something like:

-----------
Initializing structure...
Working in 2D dimensions.
Computational cell is 1 x 1 x 0 with resolution 21
     block, center = (0,0,0)
          size (1,1,1e-20)
          axes (1,0,0), (0,1,0), (0,0,1)
          dielectric constant epsilon diagonal = (1,1,1)
time for set_epsilon = 0.00232601 s
-----------
Meep progress: 1308.7380952380952/2003.3333332538605 = 65.3% done in 4.0s, 2.1s to go
on time step 54967 (time=1308.74), 7.27747e-05 s/step
harminv0:, frequency, imag. freq., Q, |amp|, amplitude, error
harminv0:, 0.49964982320796947, -6.387829272631824e-7, 391095.1607212498, 0.013285597622112183, -0.007052720387953484-0.011259051438998675i, 2.3287614080376595e-7+0.0i
harminv0:, 0.7067779190836112, 1.691710912831851e-6, -208894.41385126952, 0.005994921731086088, -0.002373006416742774+0.0055052635820590674i, 5.339765273626355e-7+0.0i
harminv0:, 0.9971987526065739, 1.6536690637644267e-7, -3015109.7775770864, 0.07795287534323195, 0.03595210583102817-0.06916716605870166i, 4.795095158586878e-8+0.0i
harminv0:, 1.1157887584886936, -1.399235272127685e-7, 3987137.7627295624, 0.06361903731428901, -0.044370834218134786-0.04559178631709551i, 3.2413084735784084e-8+0.0i
harminv0:, 1.4115657480682742, 4.0255510363580226e-8, -17532577.966584906, 0.0940011704004871, 0.06407402923620237+0.06878036648709931i, 1.0281322317081819e-8+0.0i
harminv0:, 1.4905289327310962, 9.500093703078488e-8, -7844811.742478356, 0.06351546181488635, 0.06351074111746435-7.743721771272548e-4i, 1.520738285961082e-8+0.0i
harminv0:, 1.5727074536238252, 2.96402781900332e-9, -265299037.26623282, 0.07776848493442345, 0.04957534656074576-0.059918463451407035i, 1.3776718980993644e-8+0.0i
harminv0:, 1.7956591793170913, -3.5930574726457174e-8, 24987899.483762957, 0.14467020137786865, -0.14410062932432377+0.012824811696352946i, 5.64638078279913e-9+0.0i
harminv0:, 1.9774963650323654, 5.985926891666783e-10, -1651787935.9546697, 0.32586317906949197, 0.13754120415113452+0.2954136568168976i, 2.5490377382288595e-9+0.0i
run 0 finished at t = 2003.3333333333333 (84140 timesteps)

Elapsed run time = 12.8566 s

So, what actually happened here? MEEP did a few things: setting up arrays representing the materials and boundary conditions of the problem, running the problem for an additional time of 200 after the sources and then listing the resonant frequencies detected by a probe monitoring the \(x\)-component of the electric field at position \((0.05, 0.05)\).

We have some data to analyze, onto the next post ...