Skip to content

Commit

Permalink
add needs_boundary_fix to add_srcdata
Browse files Browse the repository at this point in the history
  • Loading branch information
stevengj authored and mawc2019 committed Feb 20, 2022
1 parent 75521d2 commit ca688c6
Show file tree
Hide file tree
Showing 4 changed files with 20 additions and 18 deletions.
22 changes: 11 additions & 11 deletions src/fix_boundary_sources.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -45,28 +45,28 @@ void fields::fix_boundary_sources() {
// they are actually supposed to be located, storing info in boundarysources.
for (int i = 0; i < num_chunks; i++) {
FOR_FIELD_TYPES(ft) {
for (src_vol *src = chunks[i]->sources[ft]; src; src = src->next)
if (src->needs_boundary_fix) {
for (size_t ipt = 0; ipt < src->npts; ++ipt) {
component c = src->c;
ivec here = chunks[i]->gv.iloc(c, src->index[ipt]);
if (!chunks[i]->gv.owns(here) && src->A[ipt] != 0.0) {
if (src->t->id == 0) abort("bug: fix_boundary_sources called for non-registered source");
for (src_vol &src : chunks[i]->sources[ft])
if (src.needs_boundary_fix) {
for (size_t ipt = 0; ipt < src.num_points(); ++ipt) {
component c = src.c;
ivec here = chunks[i]->gv.iloc(c, src.index_at(ipt));
if (!chunks[i]->gv.owns(here) && src.amplitude(ipt) != 0.0) {
if (src.t()->id == 0) abort("bug: fix_boundary_sources called for non-registered source");

// find the chunk that owns this point, similar to logic in boundaries.cpp
std::complex<double> thephase;
if (locate_component_point(&c, &here, &thephase) && !on_metal_boundary(here)) {
for (int j = 0; j < num_chunks; j++)
if (chunks[j]->gv.owns(here)) {
srcpt_info s = { src->A[ipt]*conj(thephase), chunks[j]->gv.index(c, here), src->t->id, chunks[j]->chunk_idx, c };
srcpt_info s = { src.amplitude(ipt)*conj(thephase), chunks[j]->gv.index(c, here), src.t()->id, chunks[j]->chunk_idx, c };
boundarysources.push_back(s);
break;
}
}
src->A[ipt] = 0.0; // will no longer be needed
src.set_amplitude(ipt, 0.0); // will no longer be needed
}
}
src->needs_boundary_fix = false;
src.needs_boundary_fix = false;
}
}
}
Expand Down Expand Up @@ -141,7 +141,7 @@ for (int psrc = 0; psrc < P; ++psrc)
sourcedata srcdata = { (component) c, idx_arr, chunk_idx, amp_arr };
src_time *srctime = lookup_src_time(src_time_id);
if (srctime == NULL) abort("bug: unknown src_time_id (missing registration?)");
add_srcdata(srcdata, srctime);
add_srcdata(srcdata, srctime, size_t(0), NULL, false);
if (idx < N) {
chunk_idx = srcpts[idx].chunk_idx;
src_time_id = srcpts[idx].src_time_id;
Expand Down
2 changes: 1 addition & 1 deletion src/meep.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -1897,7 +1897,7 @@ class fields {
void require_source_components();
void _require_component(component c, bool aniso2d);
void require_component(component c) { _require_component(c, is_aniso2d()); sync_chunk_connections(); }
void add_srcdata(struct sourcedata cur_data, src_time *src, size_t n, std::complex<double>* amp_arr);
void add_srcdata(struct sourcedata cur_data, src_time *src, size_t n, std::complex<double>* amp_arr, bool needs_boundary_fix=false);
void register_src_time(src_time *src);
src_time *lookup_src_time(size_t id);

Expand Down
5 changes: 4 additions & 1 deletion src/meep_internals.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -55,13 +55,16 @@ class src_vol {
// Constructs a new source volume. Takes ownership of `ind` and `amps`.
// Requirement: ind.size() == amps.size()
src_vol(component cc, src_time *st, std::vector<ptrdiff_t> &&ind,
std::vector<std::complex<double> > &&amps);
std::vector<std::complex<double> > &&amps, bool fix_boundaries=false);

// Checks whether `a` and `b` are combinable, i.e. have the same indices and point to the same
// `src_time` instance, but have potentially different amplitudes.
static bool combinable(const src_vol &a, const src_vol &b);

ptrdiff_t index_at(size_t pos) const { return index[pos]; }
std::complex<double> amplitude(size_t j) const { return amp[j]; };
void set_amplitude(size_t j, std::complex<double> a) { amp[j] = a; };
void set_amplitude(size_t j, double a) { amp[j] = a; };
const std::complex<double> &amplitude_at(size_t pos) const { return amp[pos]; }
size_t num_points() const { return index.size(); };
const src_time *t() const { return src_t; };
Expand Down
9 changes: 4 additions & 5 deletions src/sources.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -166,13 +166,13 @@ bool custom_src_time::is_equal(const src_time &t) const {
/*********************************************************************/

src_vol::src_vol(component cc, src_time *st, std::vector<ptrdiff_t> &&ind,
std::vector<std::complex<double> > &&amps)
std::vector<std::complex<double> > &&amps, bool fix_boundaries)
: c([](component c) -> component {
if (is_D(c)) c = direction_component(Ex, component_direction(c));
if (is_B(c)) c = direction_component(Hx, component_direction(c));
return c;
}(cc)),
src_t(st), index(std::move(ind)), amp(std::move(amps)), needs_boundary_fix(false) {
src_t(st), index(std::move(ind)), amp(std::move(amps)), needs_boundary_fix(fix_boundaries) {
assert(index.size() == amp.size());
}

Expand Down Expand Up @@ -289,7 +289,7 @@ static void src_vol_chunkloop(fields_chunk *fc, int ichunk, component c, ivec is
fc->add_source(ft, src_vol(c, data->src, std::move(index_array), std::move(amps_array)));
}

void fields::add_srcdata(struct sourcedata cur_data, src_time *src, size_t n, std::complex<double>* amp_arr){
void fields::add_srcdata(struct sourcedata cur_data, src_time *src, size_t n, std::complex<double>* amp_arr, bool needs_boundary_fix){
if (n == 0) {
n = cur_data.idx_arr.size();
assert(amp_arr == NULL);
Expand All @@ -300,13 +300,12 @@ void fields::add_srcdata(struct sourcedata cur_data, src_time *src, size_t n, st
std::vector<ptrdiff_t> index_arr(cur_data.idx_arr);
std::vector<std::complex<double>> amplitudes(amp_arr, amp_arr+n);
component c = cur_data.near_fd_comp;

field_type ft = is_magnetic(c) ? B_stuff : D_stuff;
if (0 > cur_data.fc_idx or cur_data.fc_idx >= num_chunks) meep::abort("fields chunk index out of range");
fields_chunk *fc = chunks[cur_data.fc_idx];
if (!fc->is_mine()) meep::abort("wrong fields chunk");

fc->add_source(ft, src_vol(c, src, std::move(index_arr), std::move(amplitudes)));
fc->add_source(ft, src_vol(c, src, std::move(index_arr), std::move(amplitudes), needs_boundary_fix));
// We can't do require_component(c) since that only works if all processes are adding
// srcdata for the same components in the same order, which may not be true.
// ... instead, the caller should call fields::require_source_components()
Expand Down

0 comments on commit ca688c6

Please sign in to comment.