1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201
|
#include "meep-ctl.hpp"
using namespace meep;
using namespace std;
/**************************************************************************/
/* The following are hook functions called from main() when
starting the program and just before exiting. */
static initialize *meep_init = 0;
void ctl_start_hook(int *argc, char ***argv) {
meep_init = new initialize(*argc, *argv);
#ifdef HAVE_LIBCTL_QUIET
extern int libctl_quiet;
libctl_quiet = !am_master();
#endif
}
void ctl_stop_hook(void) { delete meep_init; }
extern "C" void SWIG_init();
void ctl_export_hook(void) { SWIG_init(); }
/**************************************************************************/
ctlio::cvector3_list do_harminv(ctlio::cnumber_list vals, double dt, double fmin, double fmax,
int maxbands, double spectral_density, double Q_thresh,
double rel_err_thresh, double err_thresh, double rel_amp_thresh,
double amp_thresh) {
complex<double> *amp = new complex<double>[maxbands];
double *freq_re = new double[maxbands];
double *freq_im = new double[maxbands];
double *freq_err = new double[maxbands];
maxbands = do_harminv(reinterpret_cast<complex<double> *>(vals.items), vals.num_items, dt, fmin,
fmax, maxbands, amp, freq_re, freq_im, freq_err, spectral_density, Q_thresh,
rel_err_thresh, err_thresh, rel_amp_thresh, amp_thresh);
ctlio::cvector3_list res;
res.num_items = maxbands;
res.items = new cvector3[maxbands];
for (int i = 0; i < maxbands; ++i) {
res.items[i].x.re = freq_re[i];
res.items[i].x.im = freq_im[i];
res.items[i].y.re = real(amp[i]);
res.items[i].y.im = imag(amp[i]);
res.items[i].z.re = freq_err[i];
res.items[i].z.im = 0;
}
delete[] freq_err;
delete[] freq_im;
delete[] freq_re;
delete[] amp;
return res;
}
kpoint_list do_get_eigenmode_coefficients(fields *f, dft_flux flux, const volume &eig_vol,
int *bands, int num_bands, int parity,
std::complex<double> *coeffs, double *vgrp,
double eig_resolution, double eigensolver_tol,
meep::kpoint_func user_kpoint_func,
void *user_kpoint_data, int dir) {
size_t num_kpoints = num_bands * flux.freq.size();
meep::vec *kpoints = new meep::vec[num_kpoints];
meep::vec *kdom = new meep::vec[num_kpoints];
double *cscale = 0; // Not needed until adjoint calculation is implemented in scheme
f->get_eigenmode_coefficients(flux, eig_vol, bands, num_bands, parity, eig_resolution,
eigensolver_tol, coeffs, vgrp, user_kpoint_func, user_kpoint_data,
kpoints, kdom, cscale,
dir < 0 ? flux.normal_direction : direction(dir));
kpoint_list res = {kpoints, num_kpoints, kdom, num_kpoints};
return res;
}
/**************************************************************************/
/* This is a wrapper function to fool SWIG...since our list constructor
takes ownership of the next pointer, we have to make sure that SWIG
does not garbage-collect volume_list objects. We do
this by wrapping a "helper" function around the constructor which
does not have the %newobject SWIG attribute. Note that we then
need to deallocate the list explicitly in Scheme. */
volume_list *make_volume_list(const volume &v, int c, complex<double> weight, volume_list *next) {
return new volume_list(v, c, weight, next);
}
/***************************************************************************/
ctlio::number_list dft_flux_flux(dft_flux *f) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->flux();
return res;
}
ctlio::number_list dft_energy_electric(dft_energy *f) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->electric();
return res;
}
ctlio::number_list dft_energy_magnetic(dft_energy *f) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->magnetic();
return res;
}
ctlio::number_list dft_energy_total(dft_energy *f) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->total();
return res;
}
ctlio::number_list dft_force_force(dft_force *f) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->force();
return res;
}
ctlio::number_list dft_ldos_ldos(dft_ldos *f) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->ldos();
return res;
}
ctlio::cnumber_list dft_ldos_F(dft_ldos *f) {
ctlio::cnumber_list res;
res.num_items = f->freq.size();
res.items = (cnumber *)f->F();
return res;
}
ctlio::cnumber_list dft_ldos_J(dft_ldos *f) {
ctlio::cnumber_list res;
res.num_items = f->freq.size();
res.items = (cnumber *)f->J();
return res;
}
ctlio::cnumber_list dft_near2far_farfield(dft_near2far *f, const vec &x) {
ctlio::cnumber_list res;
res.num_items = f->freq.size() * 6;
res.items = (cnumber *)f->farfield(x);
return res;
}
ctlio::number_list dft_near2far_flux(dft_near2far *f, direction df, const volume &where,
double resolution) {
ctlio::number_list res;
res.num_items = f->freq.size();
res.items = f->flux(df, where, resolution);
return res;
}
/***************************************************************************/
ctlio::cnumber_list make_casimir_g(double T, double dt, double sigma, meep::field_type ft,
complex<double> (*eps_func)(complex<double> omega),
double Tfft) {
ctlio::cnumber_list res;
res.num_items = int(ceil(T / dt));
res.items = new cnumber[res.num_items];
complex<double> *g = meep::make_casimir_gfunc(T, dt, sigma, ft, eps_func, Tfft);
for (int i = 0; i < res.num_items; ++i) {
res.items[i].re = real(g[i]);
res.items[i].im = imag(g[i]);
}
delete[] g;
return res;
}
ctlio::cnumber_list make_casimir_g_kz(double T, double dt, double sigma, meep::field_type ft) {
ctlio::cnumber_list res;
res.num_items = int(ceil(T / dt));
res.items = new cnumber[res.num_items];
complex<double> *g = meep::make_casimir_gfunc_kz(T, dt, sigma, ft);
for (int i = 0; i < res.num_items; ++i) {
res.items[i].re = real(g[i]);
res.items[i].im = imag(g[i]);
}
delete[] g;
return res;
}
ctlio::number_list std_vector_double_to_scm(std::vector<double> *v) {
ctlio::number_list res;
res.num_items = int(v->size());
res.items = new number[res.num_items];
for (int i = 0; i < res.num_items; ++i)
res.items[i] = (*v)[i];
return res;
}
|