Modes of a ring resonator

#include <meep.hpp>
using namespace meep;

const double n = 3.4;
const double w = 1.0;
const double r = 1.0;
const double pad = 4.0;
const double dpml = 2.0;
const double sxy = 2 * (r + w + pad + dpml);
const double origin = 0.5 * sxy;

const double fcen = 0.15;
const double df = 0.1;
const double freq_max = fcen + 0.5 * df;
const double freq_min = fcen - 0.5 * df;
const double resolution = 10.0;

const double harminv_time = 200.0;

inline double max(double a, double b) { return (a > b ? a : b); }

double eps (const vec &p)
{
  vec xx = p - vec(origin, origin);
  double pa = abs(xx);
  if (pa > r && pa <= r + w)
    return n * n;
  return 1.0;
}

int main(int argc, char **argv)
{
  initialize mpi(argc, argv);
  const volume v = vol2d(sxy, sxy, resolution);
  structure s(v, eps, pml(dpml), mirror(Y,v));
  s.set_epsilon(eps, true); // averaging
  fields f(&s);
  gaussian_src_time src(fcen, 1.0 / df, 0.0, 10.0 / df);
  f.add_point_source(Hx, src, vec(origin + r + 0.1, origin));
  f.output_hdf5(Dielectric, v.surroundings());
  h5file *ez_data = f.open_h5file("ez");
  int itot = f.last_source_time() / f.dt;
  int i = 0;
  while (f.time() < f.last_source_time()) {
    f.step();
    if (i%10 == 0) f.output_hdf5(Ez, v.surroundings(), ez_data, true);
    i++;
  }
  delete ez_data;
  int maxbands = 5;
  int ttot = int(harminv_time / f.dt)+1;
  complex<double> *mon_data = new complex<double>[ttot];
  complex<double> *amps = new complex<double>[maxbands];
  double *freq_re = new double[maxbands];
  double *freq_im = new double[maxbands];
  double *errors = new double[maxbands];
  i = 0;
  master_printf("ttot = %d\n", ttot);

  while (f.time() < f.last_source_time() + harminv_time) {
    f.step();
    mon_data[i++] = f.get_field(Ez, vec(origin + r + 0.1, origin));
  }

  int nfreq = do_harminv(mon_data, ttot, f.dt, freq_min, freq_max,
		 maxbands, amps, freq_re, freq_im, errors); 
  for (int i = 0; i < nfreq; ++i)
    master_printf("harminv0: %g, %g, %g, %g, %g+%gi, %g \n",
	  freq_re[i], freq_im[i], - 0.5 * freq_re[i] / freq_im[i],
		  abs(amps[i]), real(amps[i]), imag(amps[i]), errors[i]);
 
  delete[] errors;
  delete[] freq_im;
  delete[] freq_re;
  delete[] amps;
  delete[] mon_data;
  return 0;
}

アニメーションの作成

$ mkdir figs
$ cd figs
$ h5topng -t 0:199 -R -Zc dkbluered -a yarg -A ../eps-000000.00.h5 -d ez.r ../ez.h5
$ convert -loop 0 ez.t*.png ../rr.gif

戻る