Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

replace stability error with warning message in Lorentzian susceptibilities #150

Merged
merged 3 commits into from
Dec 22, 2017
Merged

replace stability error with warning message in Lorentzian susceptibilities #150

Changes from 1 commit
Commits
Show all changes
3 commits
Select commit Hold shift + click to select a range
058e515
fix bug due to rounding error in lorentzian_unstable
oskooi Dec 19, 2017
e8db852
replace the error with a warning message instead
oskooi Dec 20, 2017
59f4552
rm test
stevengj Dec 22, 2017
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions src/susceptibility.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -157,7 +157,7 @@ static bool lorentzian_unstable(double omega_0, double gamma, double dt) {
double w = 2*pi*omega_0, g = 2*pi*gamma;
double g2 = g*dt/2, w2 = (w*dt)*(w*dt);
double b = (1 - w2/2) / (1 + g2), c = (1 - g2) / (1 + g2);
return b*b > c && float(2*b*b - c + 2*fabs(b)*sqrt(b*b - c)) > 1;
return b*b > c && 2*b*b - c + 2*fabs(b)*sqrt(b*b - c) > 1;
}

#define SWAP(t,a,b) { t SWAP_temp = a; a = b; b = SWAP_temp; }
Expand All @@ -179,7 +179,7 @@ void lorentzian_susceptibility::update_P

if (!no_omega_0_denominator && gamma >= 0
&& lorentzian_unstable(omega_0, gamma, dt))
abort("Lorentzian pole at too high a frequency %g for stability with dt = %g: reduce the Courant factor, increase the resolution, or use a different dielectric model\n", omega_0, dt);
master_printf("WARNING: simulation is potentially unstable as the Lorentzian pole is at too high a frequency %g for stability with dt = %g. In the event the fields diverge: reduce the Courant factor, increase the resolution, or use a different dielectric model\n", omega_0, dt);
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I don't think we want to print out a warning on every timestep.

Copy link
Collaborator

@stevengj stevengj Dec 20, 2017 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe we should just omit the warning entirely. We don't warn about other potential instabilities, e.g. if you have a dielectric ε=0.1 somewhere.


FOR_COMPONENTS(c) DOCMP2 if (d->P[c][cmp]) {
const realnum *w = W[c][cmp], *s = sigma[c][component_direction(c)];
Expand Down