step_h was using the wrong timestep; fixing it caused update issues based on the Courant number being 1.0, i.e. too large. Reducing it to \sqrt{3}/3 fixes stability issues. Unfortunately, a non-Unity Courant value prevents me from using the ideal Absorbing Boundary Condition I had just implemented, so that's been disabled.
52 lines
2.0 KiB
Rust
52 lines
2.0 KiB
Rust
use coremem::{GenericSim, Driver, mat};
|
|
use coremem::geom::Index;
|
|
|
|
fn main() {
|
|
coremem::init_logging();
|
|
let width = 401;
|
|
let height = 401;
|
|
let size = Index((width, height, 1).into());
|
|
let mut driver = Driver::new(size, 1e-6 /* feature size */);
|
|
driver.add_y4m_renderer("em_reflection.y4m");
|
|
// driver.add_term_renderer();
|
|
|
|
// driver.add_upml_boundary(Index((20, 40, 0).into()));
|
|
// for y in 75..100 {
|
|
// for x in 0..width {
|
|
// // from https://www.thoughtco.com/table-of-electrical-resistivity-conductivity-608499
|
|
// // NB: different sources give pretty different values for this
|
|
// // NB: Simulation misbehaves for values > 10... Proably this model isn't so great.
|
|
// // Maybe use \eps or \xi instead of conductivity.
|
|
// let loc = Index((x, y, 0).into());
|
|
// *driver.state.get_mut(loc).mat_mut() = mat::Static::conductor(2.17).into();
|
|
// }
|
|
// }
|
|
|
|
loop {
|
|
let imp = if driver.state.step_no() < 50 {
|
|
1e4 * ((driver.state.step_no() as f64)*0.04*std::f64::consts::PI).sin()
|
|
} else {
|
|
0.0
|
|
};
|
|
// let t = 4e12 * driver.state.time();
|
|
// let imp = 1e4 * (-t*t).exp();
|
|
// state.impulse_ex(50, 50, imp);
|
|
// state.impulse_ey(50, 50, imp);
|
|
// state.impulse_bz(20, 20, (imp / 3.0e8) as _);
|
|
// state.impulse_bz(80, 20, (imp / 3.0e8) as _);
|
|
for y in height/4..height*3/4 {
|
|
//for y in 0..height {
|
|
let loc = Index((200, y, 0).into());
|
|
driver.dyn_state().impulse_ez(loc, imp as _);
|
|
driver.dyn_state().impulse_hy(loc, (imp/376.730) as _);
|
|
}
|
|
// for x in 0..width {
|
|
// let loc = Index((x, 200, 0).into());
|
|
// driver.dyn_state().impulse_ez(loc, imp as _);
|
|
// driver.dyn_state().impulse_hx(loc, (imp/376.730) as _);
|
|
// }
|
|
driver.step();
|
|
//thread::sleep(time::Duration::from_millis(67));
|
|
}
|
|
}
|