colin/fdtd-coremem
1
0
Fork 0
Code Issues Pull Requests Projects Releases Wiki Activity

Abstract coordinates behind a Coord type

Browse Source 
This will allow adding another dimension to the coordinate, more easily.
This commit is contained in:
Colin
2020-09-22 11:23:58 -07:00
parent b96f9c5dcc
commit 6d77c85179
7 changed files with 95 additions and 59 deletions
Download Patch File Download Diff File Expand all files Collapse all files

49
examples/toroid25d.rs
View File

@@ -1,4 +1,5 @@
use coremem::{consts, Driver, Flt, mat, meas};
use coremem::coord::Coord;
use coremem::geom::{Circle, Point};
fn main() {
@@ -33,26 +34,45 @@ fn main() {
driver.add_measurement(meas::Current(Circle::new(Point::new(
to_m(half_width), to_m(half_width)),
to_m(conductor_outer_rad))));
driver.add_measurement(meas::Magnetization(half_width + ferro_inner_rad + 2, half_width));
driver.add_measurement(meas::MagneticFlux(half_width + ferro_inner_rad + 2, half_width));
driver.add_measurement(meas::MagneticStrength(half_width + ferro_inner_rad + 2, half_width));
driver.add_measurement(meas::Magnetization(half_width + ferro_inner_rad + 1, half_width));
driver.add_measurement(meas::MagneticFlux(half_width + ferro_inner_rad + 1, half_width));
driver.add_measurement(meas::MagneticStrength(half_width + ferro_inner_rad + 1, half_width));
driver.add_measurement(meas::Magnetization(half_width + (ferro_inner_rad + ferro_outer_rad) / 2, half_width));
driver.add_measurement(meas::MagneticFlux(half_width + (ferro_inner_rad + ferro_outer_rad) / 2, half_width));
driver.add_measurement(meas::MagneticStrength(half_width + (ferro_inner_rad + ferro_outer_rad) / 2, half_width));
driver.add_measurement(meas::Magnetization(
(half_width + ferro_inner_rad + 2, half_width).into()
));
driver.add_measurement(meas::MagneticFlux(
(half_width + ferro_inner_rad + 2, half_width).into()
));
driver.add_measurement(meas::MagneticStrength(
(half_width + ferro_inner_rad + 2, half_width).into()
));
driver.add_measurement(meas::Magnetization(
(half_width + ferro_inner_rad + 1, half_width).into()
));
driver.add_measurement(meas::MagneticFlux(
(half_width + ferro_inner_rad + 1, half_width).into()
));
driver.add_measurement(meas::MagneticStrength(
(half_width + ferro_inner_rad + 1, half_width).into()
));
driver.add_measurement(meas::Magnetization(
(half_width + (ferro_inner_rad + ferro_outer_rad) / 2, half_width).into()
));
driver.add_measurement(meas::MagneticFlux(
(half_width + (ferro_inner_rad + ferro_outer_rad) / 2, half_width).into()
));
driver.add_measurement(meas::MagneticStrength(
(half_width + (ferro_inner_rad + ferro_outer_rad) / 2, half_width).into()
));
let center = Point::new(half_width as _, half_width as _);
for y in 0..width {
for x in 0..width {
let d = Point::new(x as _, y as _) - center;
let loc = Coord::new(x, y);
let d = Point::from(loc) - center;
let r = d.mag();
if (conductor_inner_rad as _..conductor_outer_rad as _).contains(&r) {
*driver.state.get_mut(x, y).mat_mut() = mat::Static::conductor(conductivity).into();
*driver.state.get_mut(loc).mat_mut() = mat::Static::conductor(conductivity).into();
} else if (ferro_inner_rad as _..ferro_outer_rad as _).contains(&r) {
*driver.state.get_mut(x, y).mat_mut() = mat::db::ferroxcube_3r1();
*driver.state.get_mut(loc).mat_mut() = mat::db::ferroxcube_3r1();
}
}
}
@@ -75,9 +95,10 @@ fn main() {
let e = drive_current/conductivity/area;
for y in half_width-conductor_outer_rad..half_width+conductor_outer_rad {
for x in half_width-conductor_outer_rad..half_width+conductor_outer_rad {
let d = Point::new(x as _, y as _) - center;
let loc = Coord::new(x, y);
let d = Point::from(loc) - center;
if (conductor_inner_rad as _..conductor_outer_rad as _).contains(&d.mag()) {
driver.state.impulse_ez(x, y, e);
driver.state.impulse_ez(loc, e);
}
}
}

13
src/driver.rs
View File

@@ -1,5 +1,6 @@
use crate::{flt::Flt, mat};
use crate::consts;
use crate::coord::Coord;
use crate::meas::{self, AbstractMeasurement};
use crate::render::{self, MultiRenderer, Renderer};
use crate::sim::SimState;
@@ -22,7 +23,7 @@ pub struct Driver {
impl Driver {
pub fn new(width: u32, height: u32, feature_size: Flt) -> Self {
Driver {
state: SimState::new(width, height, feature_size),
state: SimState::new(Coord::new(width, height), feature_size),
renderer: Default::default(),
steps_per_frame: 1,
time_spent_stepping: Default::default(),
@@ -61,15 +62,15 @@ impl Driver {
let conductivity = base_conductivity * (depth*depth) as Flt;
for x in inset..self.state.width() - inset {
// left
*self.state.get_mut(x, inset).mat_mut() = mat::Static::conductor(conductivity).into();
*self.state.get_mut((x, inset).into()).mat_mut() = mat::Static::conductor(conductivity).into();
// right
*self.state.get_mut(x, self.state.height() - 1 - inset).mat_mut() = mat::Static::conductor(conductivity).into();
*self.state.get_mut((x, self.state.height() - 1 - inset).into()).mat_mut() = mat::Static::conductor(conductivity).into();
}
for y in inset..self.state.height() - inset {
// top
*self.state.get_mut(inset, y).mat_mut() = mat::Static::conductor(conductivity).into();
*self.state.get_mut((inset, y).into()).mat_mut() = mat::Static::conductor(conductivity).into();
// bottom
*self.state.get_mut(self.state.width() - 1 - inset, y).mat_mut() = mat::Static::conductor(conductivity).into();
*self.state.get_mut((self.state.width() - 1 - inset, y).into()).mat_mut() = mat::Static::conductor(conductivity).into();
}
}
}
@@ -100,7 +101,7 @@ impl Driver {
let cond_x = scale * (depth_x as Flt/thickness as Flt).powf(3.0);
let cond_y = scale * (depth_y as Flt/thickness as Flt).powf(3.0);
let conductor = mat::Static::anisotropic_conductor(Vec3::new(cond_x, cond_y, 0.0));
*self.state.get_mut(x, y).mat_mut() = conductor.into();
*self.state.get_mut((x, y).into()).mat_mut() = conductor.into();
}
}
}

8
src/geom.rs
View File

@@ -1,4 +1,6 @@
use crate::coord::Coord;
use crate::flt::{Flt, Real};
use std::convert::From;
use std::fmt::{self, Display};
use std::ops::{Add, AddAssign, Mul, Neg, Sub};
@@ -16,6 +18,12 @@ pub struct Point {
pub y: Real,
}
impl From<Coord> for Point {
fn from(c: Coord) -> Self {
Self::new(c.x() as _, c.y() as _)
}
}
impl Add for Point {
type Output = Point;
fn add(self, other: Point) -> Point {

1
src/lib.rs
View File

@@ -7,6 +7,7 @@
use log::info;
pub mod coord;
pub mod driver;
pub mod geom;
pub mod mat;

33
src/meas.rs
View File

@@ -1,3 +1,4 @@
use crate::coord::Coord;
use crate::geom::{Point, Region};
use crate::mat::Material as _;
use crate::sim::{Cell, SimSnapshot};
@@ -38,9 +39,9 @@ impl AbstractMeasurement for Label {
}
}
fn sum_over_region<T: Default + Send + Sum<T>, R: Region + Sync, F: Fn(u32, u32, &Cell) -> T + Sync>(state: &SimSnapshot, r: &R, f: F) -> T {
state.map_sum_enumerated(|x, y, cell| if r.contains(Point::new(x as _, y as _)*state.feature_size()) {
f(x, y, cell)
fn sum_over_region<T: Default + Send + Sum<T>, R: Region + Sync, F: Fn(Coord, &Cell) -> T + Sync>(state: &SimSnapshot, r: &R, f: F) -> T {
state.map_sum_enumerated(|coord, cell| if r.contains(Point::from(coord)*state.feature_size()) {
f(coord, cell)
} else {
Default::default()
})
@@ -50,47 +51,47 @@ pub struct Current<R>(pub R);
impl<R: Region + Display + Sync> AbstractMeasurement for Current<R> {
fn eval(&self, state: &SimSnapshot) -> String {
let current = sum_over_region(state, &self.0, |x, y, _cell| state.current(x, y));
let current = sum_over_region(state, &self.0, |coord, _cell| state.current(coord));
format!("I({}): ({:.2e}, {:.2e}, {:.2e})", self.0, current.x(), current.y(), current.z())
}
}
/// M
pub struct Magnetization(pub u32, pub u32);
pub struct Magnetization(pub Coord);
impl AbstractMeasurement for Magnetization {
fn eval(&self, state: &SimSnapshot) -> String {
let m = state.get(self.0, self.1).mat().m();
format!("M({}, {}): ({:.2e}, {:.2e}, {:.2e})", self.0, self.1, m.x(), m.y(), m.z())
let m = state.get(self.0).mat().m();
format!("M{}: ({:.2e}, {:.2e}, {:.2e})", self.0, m.x(), m.y(), m.z())
}
}
/// B
pub struct MagneticFlux(pub u32, pub u32);
pub struct MagneticFlux(pub Coord);
impl AbstractMeasurement for MagneticFlux {
fn eval(&self, state: &SimSnapshot) -> String {
let b = state.get(self.0, self.1).b();
format!("B({}, {}): ({:.2e}, {:.2e}, {:.2e})", self.0, self.1, b.x(), b.y(), b.z())
let b = state.get(self.0).b();
format!("B{}: ({:.2e}, {:.2e}, {:.2e})", self.0, b.x(), b.y(), b.z())
}
}
/// H
pub struct MagneticStrength(pub u32, pub u32);
pub struct MagneticStrength(pub Coord);
impl AbstractMeasurement for MagneticStrength {
fn eval(&self, state: &SimSnapshot) -> String {
let h = state.get(self.0, self.1).h();
format!("H({}, {}): ({:.2e}, {:.2e}, {:.2e})", self.0, self.1, h.x(), h.y(), h.z())
let h = state.get(self.0).h();
format!("H{}: ({:.2e}, {:.2e}, {:.2e})", self.0, h.x(), h.y(), h.z())
}
}
pub struct ElectricField(pub u32, pub u32);
pub struct ElectricField(pub Coord);
impl AbstractMeasurement for ElectricField {
fn eval(&self, state: &SimSnapshot) -> String {
let e = state.get(self.0, self.1).e();
format!("E({}, {}): ({:.2e}, {:.2e}, {:.2e})", self.0, self.1, e.x(), e.y(), e.z())
let e = state.get(self.0).e();
format!("E{}: ({:.2e}, {:.2e}, {:.2e})", self.0, e.x(), e.y(), e.z())
}
}

8
src/render.rs
View File

@@ -125,7 +125,7 @@ impl<'a> RenderSteps<'a> {
let h = self.sim.height();
for y in 0..h {
for x in 0..w {
let cell = self.sim.get(x, y);
let cell = self.sim.get((x, y).into());
let value = measure(cell);
let scaled = if signed {
scale_signed_to_u8(value, typical)
@@ -170,7 +170,7 @@ impl<'a> RenderSteps<'a> {
let yend = (ystart + size).min(h);
for y in ystart..yend {
for x in xstart..xend {
field += measure(self.sim.get(x, y));
field += measure(self.sim.get((x, y).into()));
}
}
let xw = xend - xstart;
@@ -217,12 +217,12 @@ impl Renderer for NumericTermRenderer {
fn render(&mut self, state: &SimSnapshot, _measurements: &[Box<dyn AbstractMeasurement>]) {
for y in 0..state.height() {
for x in 0..state.width() {
let cell = state.get(x, y);
let cell = state.get((x, y).into());
print!(" {:>10.1e}", cell.ex());
}
print!("\n");
for x in 0..state.width() {
let cell = state.get(x, y);
let cell = state.get((x, y).into());
print!("{:>10.1e} {:>10.1e}", cell.ey(), cell.bz());
}
print!("\n");

42
src/sim.rs
View File

@@ -1,10 +1,12 @@
use crate::{flt::{Flt, Real}, consts};
use crate::coord::Coord;
use crate::mat::{self, GenericMaterial, Material};
use crate::vec3::Vec3;
use log::trace;
use ndarray::{Array2, Zip};
use ndarray::parallel::prelude::*;
use std::convert::From;
use std::iter::Sum;
pub type SimSnapshot = SimState<mat::Static>;
@@ -18,10 +20,10 @@ pub struct SimState<M=GenericMaterial> {
}
impl<M: Material + Default> SimState<M> {
pub fn new(width: u32, height: u32, feature_size: Flt) -> Self {
pub fn new(size: Coord, feature_size: Flt) -> Self {
Self {
cells: Array2::default((height as _, width as _)),
scratch: Array2::default((height as _, width as _)),
cells: Array2::default((size.y() as _, size.x() as _)),
scratch: Array2::default((size.y() as _, size.x() as _)),
feature_size: feature_size.into(),
..Default::default()
}
@@ -103,17 +105,19 @@ impl<M: Material + Default + Send + Sync> SimState<M> {
self.cells.view().into_par_iter().map(f).sum()
}
pub fn map_sum_enumerated<F: Fn(u32, u32, &Cell<M>) -> R + Sync + Send, R: Sum<R> + Send>(&self, f: F) -> R {
Zip::from(ndarray::indices_of(&self.cells)).and(&self.cells).into_par_iter().map(|((x, y), c)| f(x as _, y as _, c)).sum()
pub fn map_sum_enumerated<F: Fn(Coord, &Cell<M>) -> R + Sync + Send, R: Sum<R> + Send>(&self, f: F) -> R {
Zip::from(ndarray::indices_of(&self.cells)).and(&self.cells).into_par_iter().map(|((y, x), c)| {
f(Coord::new(x as _, y as _), c)
}).sum()
}
}
impl<M: Material> SimState<M> {
pub fn impulse_bz(&mut self, x: u32, y: u32, bz: Flt) {
self.get_mut(x, y).impulse_bz(bz);
pub fn impulse_bz(&mut self, c: Coord, bz: Flt) {
self.get_mut(c).impulse_bz(bz);
}
pub fn current(&self, x: u32, y: u32) -> Vec3 {
self.get(x, y).current_density() * (self.feature_size * self.feature_size)
pub fn current(&self, c: Coord) -> Vec3 {
self.get(c).current_density() * (self.feature_size * self.feature_size)
}
}
@@ -130,14 +134,14 @@ impl<M> SimState<M> {
self.feature_size.into()
}
pub fn impulse_ex(&mut self, x: u32, y: u32, ex: Flt) {
self.get_mut(x, y).state.e += Vec3::new(ex, 0.0, 0.0);
pub fn impulse_ex(&mut self, c: Coord, ex: Flt) {
self.get_mut(c).state.e += Vec3::new(ex, 0.0, 0.0);
}
pub fn impulse_ey(&mut self, x: u32, y: u32, ey: Flt) {
self.get_mut(x, y).state.e += Vec3::new(0.0, ey, 0.0);
pub fn impulse_ey(&mut self, c: Coord, ey: Flt) {
self.get_mut(c).state.e += Vec3::new(0.0, ey, 0.0);
}
pub fn impulse_ez(&mut self, x: u32, y: u32, ez: Flt) {
self.get_mut(x, y).state.e += Vec3::new(0.0, 0.0, ez);
pub fn impulse_ez(&mut self, c: Coord, ez: Flt) {
self.get_mut(c).state.e += Vec3::new(0.0, 0.0, ez);
}
pub fn width(&self) -> u32 {
@@ -146,11 +150,11 @@ impl<M> SimState<M> {
pub fn height(&self) -> u32 {
self.cells.shape()[0] as _
}
pub fn get(&self, x: u32, y: u32) -> &Cell<M> {
&self.cells[[y as _, x as _]]
pub fn get(&self, c: Coord) -> &Cell<M> {
&self.cells[c.row_major_idx()]
}
pub fn get_mut(&mut self, x: u32, y: u32) -> &mut Cell<M> {
&mut self.cells[[y as _, x as _]]
pub fn get_mut(&mut self, c: Coord) -> &mut Cell<M> {
&mut self.cells[c.row_major_idx()]
}
pub fn timestep(&self) -> Flt {