diff --git a/examples/toroid25d.rs b/examples/toroid25d.rs index 84de66e..d5ca1d0 100644 --- a/examples/toroid25d.rs +++ b/examples/toroid25d.rs @@ -1,5 +1,6 @@ use coremem::{consts, Driver, Flt, mat, meas}; -use coremem::geom::{Coord, CylinderZ, Vec2}; +use coremem::geom::{Coord, CylinderZ, Vec2, Vec3}; +use coremem::stim::{Stimulus, Sinusoid}; fn main() { coremem::init_logging(); @@ -35,12 +36,12 @@ fn main() { m_to_um(ferro_inner_rad), m_to_um(ferro_outer_rad), )); + let conductor_region = CylinderZ::new( + Vec2::new(half_width, half_width), + conductor_outer_rad); // driver.add_term_renderer(); driver.add_measurement(meas::Label(format!("Conductivity: {}, Imax: {:.2e}", conductivity, peak_current))); - driver.add_measurement(meas::Current(CylinderZ::new( - Vec2::new(half_width, half_width), - conductor_outer_rad) - )); + driver.add_measurement(meas::Current(conductor_region.clone())); driver.add_measurement(meas::Magnetization( (half_width + ferro_inner_rad + 2.0*feat_size, half_width, half_depth).into() )); @@ -92,28 +93,11 @@ fn main() { let boundary = Coord::new(from_m(boundary_xy), from_m(boundary_xy), 20); driver.add_upml_boundary(boundary); + driver.add_stimulus(Stimulus::new( + conductor_region.clone(), + Sinusoid::new(Vec3::new(0.0, 0.0, peak_current * 1e-18), 1e9))); + loop { - // let drive_current = peak_current * match driver.state.step_no() { - // 0..=1000 => 1.0, - // 3000..=4000 => -1.0, - // _ => 0.0, - // }; - let drive_current = peak_current * 1e-18; - // J = \sigma*E = [Am^-2] - // I = \sigma*E*Area - // E = I / \sigma / Area - //let e = v/(2.0*feat_size); - let area = consts::PI*(conductor_outer_rad*conductor_outer_rad - conductor_inner_rad*conductor_inner_rad); - let e = drive_current/conductivity/area; - for y_px in from_m(half_width-conductor_outer_rad)..from_m(half_width+conductor_outer_rad) { - for x_px in from_m(half_width-conductor_outer_rad)..from_m(half_width+conductor_outer_rad) { - let loc = Coord::new(x_px, y_px, 0); - let d = Vec2::new(to_m(x_px), to_m(y_px)) - center; - if (conductor_inner_rad..conductor_outer_rad).contains(&d.mag()) { - driver.state.impulse_ez(loc, e); - } - } - } driver.step(); } } diff --git a/src/driver.rs b/src/driver.rs index 1eb0e7d..d3dcc4c 100644 --- a/src/driver.rs +++ b/src/driver.rs @@ -4,6 +4,7 @@ use crate::geom::{Coord, Vec3}; use crate::meas::{self, AbstractMeasurement}; use crate::render::{self, MultiRenderer, Renderer}; use crate::sim::{GenericSim as _, SimState}; +use crate::stim::AbstractStimulus; use log::{info, trace}; use std::path::PathBuf; @@ -14,25 +15,32 @@ pub struct Driver { renderer: MultiRenderer, steps_per_frame: u64, time_spent_stepping: Duration, + time_spent_on_stimuli: Duration, time_spent_rendering: Duration, measurements: Vec>, + stimuli: Vec>, start_time: SystemTime, } impl Driver { - // TODO: allow depth pub fn new(size: Coord, feature_size: Flt) -> Self { Driver { state: SimState::new(size, feature_size), renderer: Default::default(), steps_per_frame: 1, time_spent_stepping: Default::default(), + time_spent_on_stimuli: Default::default(), time_spent_rendering: Default::default(), measurements: vec![Box::new(meas::Time), Box::new(meas::Meta), Box::new(meas::Energy)], + stimuli: vec![], start_time: SystemTime::now(), } } + pub fn add_stimulus(&mut self, s: S) { + self.stimuli.push(Box::new(s)) + } + pub fn add_measurement(&mut self, m: M) { self.measurements.push(Box::new(m)); } @@ -126,6 +134,15 @@ impl Driver { self.time_spent_rendering += start_time.elapsed().unwrap(); trace!("render end"); } + { + trace!("stimuli begin"); + let start_time = SystemTime::now(); + for stim in &mut *self.stimuli { + stim.apply(&mut self.state); + } + self.time_spent_on_stimuli += start_time.elapsed().unwrap(); + } + trace!("step begin"); let start_time = SystemTime::now(); self.state.step(); @@ -133,19 +150,21 @@ impl Driver { trace!("step end"); let step = self.state.step_no(); if step % (10*self.steps_per_frame) == 0 { - let driver_time = self.time_spent_stepping.as_secs_f64() as Flt; - let render_time = self.time_spent_rendering.as_secs_f64() as Flt; - let overall_time = self.start_time.elapsed().unwrap().as_secs_f64() as Flt; - let fps = (self.state.step_no() as Flt) / overall_time; - let sim_time = self.state.time(); + let step_time = self.time_spent_stepping.as_secs_f64(); + let stim_time = self.time_spent_on_stimuli.as_secs_f64(); + let render_time = self.time_spent_rendering.as_secs_f64(); + let overall_time = self.start_time.elapsed().unwrap().as_secs_f64(); + let fps = (self.state.step_no() as f64) / overall_time; + let sim_time = self.state.time() as f64; info!( - "t={:.2e} frame {:06} fps: {:6.2} (sim: {:.1}s, render: {:.1}s, other: {:.1}s)", + "t={:.2e} frame {:06} fps: {:6.2} (sim: {:.1}s, stim: {:.1}s, render: {:.1}s, other: {:.1}s)", sim_time, step, fps, - driver_time, + step_time, + stim_time, render_time, - overall_time - driver_time - render_time + overall_time - step_time - stim_time - render_time ); } } diff --git a/src/geom/mod.rs b/src/geom/mod.rs index 0463eca..17d29ef 100644 --- a/src/geom/mod.rs +++ b/src/geom/mod.rs @@ -183,6 +183,7 @@ pub trait Region { fn contains(&self, p: Vec3) -> bool; } +#[derive(Copy, Clone)] pub struct CylinderZ { center: Vec2, radius: Real, diff --git a/src/lib.rs b/src/lib.rs index c45caf3..582cb17 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -13,6 +13,7 @@ pub mod mat; pub mod meas; pub mod render; pub mod sim; +pub mod stim; pub use driver::*; pub use mat::*; @@ -65,6 +66,7 @@ pub mod consts { // Vacuum Permittivity pub const EPS0: Flt = 8.854187812813e-12; // F⋅m−1 pub const PI: Flt = std::f64::consts::PI as Flt; + pub const TWO_PI: Flt = 2.0 * std::f64::consts::PI as Flt; pub mod real { use super::*; pub(crate) fn C() -> Real { diff --git a/src/stim.rs b/src/stim.rs new file mode 100644 index 0000000..5ac5ae7 --- /dev/null +++ b/src/stim.rs @@ -0,0 +1,56 @@ +use crate::consts; +use crate::flt::Flt; +use crate::geom::{Region, Vec3}; +use crate::sim::GenericSim; + +pub trait AbstractStimulus { + // /// Returns a Region over which to apply some net impulse (E). + // /// i.e. the caller will divide E by the volume of the Region and then + // /// apply that to each cell in the region. + // fn e_for(&mut self, sim: &dyn GenericSim) -> (Box, Vec3); + fn apply(&mut self, sim: &mut dyn GenericSim); +} + +pub trait TimeVarying { + /// Retrieve the E impulse to apply at the provided time (in seconds). + fn at(&mut self, t_sec: Flt) -> Vec3; +} + +pub struct Stimulus { + region: R, + time: T, +} + +impl Stimulus { + pub fn new(region: R, time: T) -> Self { + Self { + region, time + } + } +} + +impl AbstractStimulus for Stimulus { + fn apply(&mut self, sim: &mut dyn GenericSim) { + unimplemented!() + } +} + +pub struct Sinusoid { + amp: Vec3, + omega: Flt, +} + +impl Sinusoid { + pub fn new(amp: Vec3, freq: Flt) -> Self { + Self { + amp, + omega: freq * consts::TWO_PI, + } + } +} + +impl TimeVarying for Sinusoid { + fn at(&mut self, t_sec: Flt) -> Vec3 { + self.amp * (t_sec * self.omega).sin() + } +}