From 0953062e68255c30a30c82d9a51fd42d804ff8af Mon Sep 17 00:00:00 2001 From: "Steven G. Johnson" Date: Sat, 25 Feb 2023 10:29:51 -0500 Subject: [PATCH 1/4] disable default decimation for nonlinear materials --- src/dft.cpp | 2 +- src/fields.cpp | 7 +++++++ src/meep.hpp | 9 ++++++++- src/structure.cpp | 10 ++++++++++ 4 files changed, 26 insertions(+), 2 deletions(-) diff --git a/src/dft.cpp b/src/dft.cpp index d230a5265..810f9c5d8 100644 --- a/src/dft.cpp +++ b/src/dft.cpp @@ -204,7 +204,7 @@ dft_chunk *fields::add_dft(component c, const volume &where, const double *freq, double freq_max = 0; for (size_t i = 0; i < Nfreq; ++i) freq_max = std::max(freq_max, std::abs(freq[i])); - if ((freq_max > 0) && (src_freq_max > 0)) + if ((freq_max > 0) && (src_freq_max > 0) && !has_nonlinearities()) decimation_factor = std::max(1, int(std::floor(1 / (dt * (freq_max + src_freq_max))))); else decimation_factor = 1; diff --git a/src/fields.cpp b/src/fields.cpp index d003148cd..8e6e422fa 100644 --- a/src/fields.cpp +++ b/src/fields.cpp @@ -663,6 +663,13 @@ void fields_chunk::use_real_fields() { } } +bool fields::has_nonlinearities() const { + bool nonlinear = false; + for (int i = 0; i < num_chunks; i++) + if (chunks[i]->is_mine()) nonlinear = nonlinear || chunks[i]->s->has_nonlinearities(); + return or_to_all(nonlinear); +} + int fields::phase_in_material(const structure *snew, double time) { if (snew->num_chunks != num_chunks) meep::abort("Can only phase in similar sets of chunks: %d vs %d\n", snew->num_chunks, diff --git a/src/meep.hpp b/src/meep.hpp index 8dee6b1b2..047a1393a 100644 --- a/src/meep.hpp +++ b/src/meep.hpp @@ -105,7 +105,9 @@ class susceptibility { int get_id() const { return id; } bool operator==(const susceptibility &s) const { return id == s.id; }; - // Returns the 1st order nonlinear susceptibility (generic) + virtual bool has_nonlinearities() const { return next ? next->has_nonlinearities() : false; } + + // Returns the 1st order (linear) susceptibility (generic) virtual std::complex chi1(realnum freq, realnum sigma = 1); // update all of the internal polarization state given the W field @@ -345,6 +347,8 @@ class multilevel_susceptibility : public susceptibility { virtual susceptibility *clone() const { return new multilevel_susceptibility(*this); } virtual ~multilevel_susceptibility(); + virtual bool has_nonlinearities() const { return true; } + virtual void update_P(realnum *W[NUM_FIELD_COMPONENTS][2], realnum *W_prev[NUM_FIELD_COMPONENTS][2], realnum dt, const grid_volume &gv, void *P_internal_data) const; @@ -616,6 +620,8 @@ class structure_chunk { pml_profile_func pml_profile, void *pml_profile_data, double pml_profile_integral, double pml_profile_integral_u); + bool has_nonlinearities() const; + void add_susceptibility(material_function &sigma, field_type ft, const susceptibility &sus); void mix_with(const structure_chunk *, double); @@ -1753,6 +1759,7 @@ class fields { void reset(); void log(const char *prefix = ""); void change_m(double new_m); + bool has_nonlinearities() const; // time.cpp std::vector time_spent_on(time_sink sink); diff --git a/src/structure.cpp b/src/structure.cpp index a4d905566..6ef42cb2a 100644 --- a/src/structure.cpp +++ b/src/structure.cpp @@ -545,6 +545,16 @@ bool structure_chunk::has_chi1inv(component c, direction d) const { return is_mine() && chi1inv[c][d] && !trivial_chi1inv[c][d]; } +bool structure_chunk::has_nonlinearities() const { + bool nonlinear = false; + if (!is_mine()) return false; + FOR_COMPONENTS(c) { nonlinear = nonlinear || chi2[c] || chi3[c]; } + FOR_FIELD_TYPES(ft) { + if (chiP[ft]) nonlinear = nonlinear || chiP[ft]->has_nonlinearities(); + } + return nonlinear; +} + void structure::mix_with(const structure *oth, double f) { if (num_chunks != oth->num_chunks) meep::abort("You can't phase materials with different chunk topologies...\n"); From ae9503d58999850a4c8c5a1ef80766f3d384368d Mon Sep 17 00:00:00 2001 From: "Steven G. Johnson" Date: Sun, 26 Feb 2023 08:10:44 -0500 Subject: [PATCH 2/4] move has_nonlinearities check --- src/dft.cpp | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/dft.cpp b/src/dft.cpp index 810f9c5d8..287a4f2b6 100644 --- a/src/dft.cpp +++ b/src/dft.cpp @@ -192,7 +192,9 @@ dft_chunk *fields::add_dft(component c, const volume &where, const double *freq, data.c = c; data.vc = vc; - if (decimation_factor == 0) { + if (decimation_factor == 0 && has_nonlinearities()) + decimation_factor = 1; + else if (decimation_factor == 0) { double src_freq_max = 0; for (src_time *s = sources; s; s = s->next) { if (s->get_fwidth() == 0) @@ -204,7 +206,7 @@ dft_chunk *fields::add_dft(component c, const volume &where, const double *freq, double freq_max = 0; for (size_t i = 0; i < Nfreq; ++i) freq_max = std::max(freq_max, std::abs(freq[i])); - if ((freq_max > 0) && (src_freq_max > 0) && !has_nonlinearities()) + if ((freq_max > 0) && (src_freq_max > 0)) decimation_factor = std::max(1, int(std::floor(1 / (dt * (freq_max + src_freq_max))))); else decimation_factor = 1; From 73223af69374d84ee85344476775509134f07afe Mon Sep 17 00:00:00 2001 From: "Steven G. Johnson" Date: Sun, 26 Feb 2023 14:32:46 -0500 Subject: [PATCH 3/4] Revert "move has_nonlinearities check" This reverts commit ae9503d58999850a4c8c5a1ef80766f3d384368d. --- src/dft.cpp | 6 ++---- 1 file changed, 2 insertions(+), 4 deletions(-) diff --git a/src/dft.cpp b/src/dft.cpp index 287a4f2b6..810f9c5d8 100644 --- a/src/dft.cpp +++ b/src/dft.cpp @@ -192,9 +192,7 @@ dft_chunk *fields::add_dft(component c, const volume &where, const double *freq, data.c = c; data.vc = vc; - if (decimation_factor == 0 && has_nonlinearities()) - decimation_factor = 1; - else if (decimation_factor == 0) { + if (decimation_factor == 0) { double src_freq_max = 0; for (src_time *s = sources; s; s = s->next) { if (s->get_fwidth() == 0) @@ -206,7 +204,7 @@ dft_chunk *fields::add_dft(component c, const volume &where, const double *freq, double freq_max = 0; for (size_t i = 0; i < Nfreq; ++i) freq_max = std::max(freq_max, std::abs(freq[i])); - if ((freq_max > 0) && (src_freq_max > 0)) + if ((freq_max > 0) && (src_freq_max > 0) && !has_nonlinearities()) decimation_factor = std::max(1, int(std::floor(1 / (dt * (freq_max + src_freq_max))))); else decimation_factor = 1; From 50118dbc8e54ac022e9e17400f43cc56e5bac560 Mon Sep 17 00:00:00 2001 From: "Steven G. Johnson" Date: Sun, 26 Feb 2023 14:37:22 -0500 Subject: [PATCH 4/4] globally sync decimation_factor --- src/dft.cpp | 7 ++++++- src/fields.cpp | 4 ++-- src/meep.hpp | 2 +- 3 files changed, 9 insertions(+), 4 deletions(-) diff --git a/src/dft.cpp b/src/dft.cpp index 810f9c5d8..c8256b1e5 100644 --- a/src/dft.cpp +++ b/src/dft.cpp @@ -204,10 +204,15 @@ dft_chunk *fields::add_dft(component c, const volume &where, const double *freq, double freq_max = 0; for (size_t i = 0; i < Nfreq; ++i) freq_max = std::max(freq_max, std::abs(freq[i])); - if ((freq_max > 0) && (src_freq_max > 0) && !has_nonlinearities()) + if ((freq_max > 0) && (src_freq_max > 0) && !has_nonlinearities(false)) decimation_factor = std::max(1, int(std::floor(1 / (dt * (freq_max + src_freq_max))))); else decimation_factor = 1; + + // with add_srcdata sources, it's possible that not all + // sources are present on all chunks, leading us to over-estimate + // the allowed decimation_factor -- take minimimum to be sure: + decimation_factor = min_to_all(decimation_factor); } data.decimation_factor = decimation_factor; diff --git a/src/fields.cpp b/src/fields.cpp index 8e6e422fa..9a5ae2b96 100644 --- a/src/fields.cpp +++ b/src/fields.cpp @@ -663,11 +663,11 @@ void fields_chunk::use_real_fields() { } } -bool fields::has_nonlinearities() const { +bool fields::has_nonlinearities(bool parallel) const { bool nonlinear = false; for (int i = 0; i < num_chunks; i++) if (chunks[i]->is_mine()) nonlinear = nonlinear || chunks[i]->s->has_nonlinearities(); - return or_to_all(nonlinear); + return parallel ? or_to_all(nonlinear) : nonlinear; } int fields::phase_in_material(const structure *snew, double time) { diff --git a/src/meep.hpp b/src/meep.hpp index 047a1393a..a1f728d40 100644 --- a/src/meep.hpp +++ b/src/meep.hpp @@ -1759,7 +1759,7 @@ class fields { void reset(); void log(const char *prefix = ""); void change_m(double new_m); - bool has_nonlinearities() const; + bool has_nonlinearities(bool parallel = true) const; // time.cpp std::vector time_spent_on(time_sink sink);