\begin{comment}
/*
\end{comment}
\section{A simple 2D system.}
\begin{comment}
*/
\end{comment}

This example is intended to let you quickly get started using meep to run
a simple calculation.  As such, it will include within it the complete code
of the example itself.  Meep is a C++ library, so your control file is a
small C++ program.

At the beginning of your control file, you have to include the ``meep.h''
header and use the ``meep'' namespace...
\begin{verbatim}
#include <meep.h>
using namespace meep;
\end{verbatim}

Next we create a function to define epsilon.  This function accepts a
``vec'' argument, and returns a double, which is the value of epsilon.  For
this example, we use an index-guided waveguide with some air slits cut in
it.  You can choose whatever units you like in which to define your
structure.  In this case we choose the width of the waveguide as our unit,
which is also equal to 1 micron.

\begin{verbatim}
const double half_cavity_width = 0.5*0.68, air_slit_width = 0.38,
             grating_periodicity = 0.48,
             half_waveguide_width = 1.0,
             num_air_slits = 15.0,
             high_dielectric = 12.0, low_dielectric = 11.5;
const double pml_thickness = 1.0;
const double x_center = 7.7 + pml_thickness;
const double y_center = 10.5 + pml_thickness;
double eps(const vec &rr) {
  // Displacement from center of cavity is r:
  const vec r = rr - vec2d(x_center, y_center);
  // First the air slits:
  double dx = fabs(r.x()) - half_cavity_width;
  if (dx < num_air_slits*grating_periodicity && dx > 0.0) {
    while (dx > grating_periodicity) dx -= grating_periodicity;
    if (dx < air_slit_width) return 1.0;
  }
  // Now check if the y value is within the waveguide:
  if (fabs(r.y()) < half_waveguide_width) return high_dielectric;
  // Otherwise we must be in the surrounding low dielectric:
  return low_dielectric;
}
\end{verbatim}
\begin{figure}
\label{simple_figure}
\caption{$E_z$}
\begin{center}
\includegraphics[width=4.8cm,clip=true]{simple-out/ez-000200-00}
\end{center}
\end{figure}
The main function should always start by creating an initialize object.
This object is responsible for setting up MPI if we are running on multiple
processors, and cleaning up properly when it is deleted (which means we are
done).
\begin{verbatim}
int main(int argc, char *argv[]) {
  initialize mpi(argc, argv);
\end{verbatim}
The ``mat'' structure defines the contents of the unit cell.  It needs a
volume, which includes the size of the grid and the resolution, as well as
the epsilon function we defined earlier.  Here we also choose to use PML
absorbing boundary conditions in all directions, since we are interested in
the high Q mode in the cavity.
\begin{verbatim}
  const double amicron = 10; // a micron is this many grid points.
  const volume vol = voltwo(2*x_center, 2*y_center, amicron);
  const symmetry S = mirror(Y, vol) + rotate2(Y, vol);
  structure s(vol, eps, 0, S);
  s.use_pml_everywhere(pml_thickness);
\end{verbatim}
To avoid clutter, we'll create a directory to hold the output.  The
function \verb!make_output_directory! creates a directory based on the name
of the example program along with an extension.  It also backs up the C++
source file if it can find it.  If the directory already exists, then it
reuses it, unless the C++ control file has changed, in which case it
creates a new one.
\begin{verbatim}
  const char *dirname = make_output_directory(__FILE__);
  s.set_output_directory(dirname);
\end{verbatim}
The structure only holds the epsilon.  We will also need a ``fields''
object to hold our electric and magnetic fields.  We add a point source
oriented in the $E_z$ direction, located in the center of our cavity.
\begin{verbatim}
  fields f(&s);
  const double wavelength = 1.72;
  const double freq = 1.0/wavelength;
  f.add_point_source(Hy, freq, 10.0, 0.0, 5.0, vec2d(x_center,y_center));
\end{verbatim}
I'm not interested in seeing the source itself, so I'll keep time stepping
until the current time is greater than the last time at which the source is
running.
\begin{verbatim}
  while (f.time() < f.find_last_source()) f.step();
\end{verbatim}
Now we'll wait a bit (to let the low-Q modes die off) and then take a
snapshot of the fields in encapsulated postscript format.
\begin{verbatim}
  while (f.time() < 200.0) f.step();
  f.eps_slices();
\end{verbatim}
And now we're done, although you might wonder if we've done anything
worthwhile, since all we got out of this was a picture... All that is left
is (as a matter of principle) to delete the string containing the directory
name of our output directory.
\begin{verbatim}
  delete[] dirname;
}
\end{verbatim}