cross: split supporting types out of step.rs

this focuses the core "business logic" into a more narrow slice.
good for organization, and also to highlight where the complexity lies
(i.e. most valuable places to test).
This commit is contained in:
2022-08-24 17:45:28 -07:00
parent 1cfef7cac6
commit 9e49538ba7
2 changed files with 162 additions and 149 deletions

View File

@@ -0,0 +1,156 @@
use core::ops::{Index, IndexMut};
use crate::dim::DimSlice;
use crate::mat::Material;
use crate::real::Real;
use crate::vec::{Vec3, Vec3u};
mod support;
pub use support::{
SimMeta,
VolumeSampleNeg,
VolumeSamplePos,
};
pub struct StepEContext<'a, R, M> {
pub inv_feature_size: R,
pub time_step: R,
pub stim_e: Vec3<R>,
pub mat: &'a M,
/// Input field sampled near this location
pub in_h: VolumeSampleNeg<R>,
pub in_e: Vec3<R>,
}
impl<'a, R: Real, M: Material<R>> StepEContext<'a, R, M> {
pub fn step_flat_view<RM, RF, WF>(
meta: SimMeta<R>,
mat: &RM,
stim_e: &RF,
e: &mut WF,
h: &RF,
idx: Vec3u,
)
where
RM: Index<usize, Output=M> + ?Sized,
RF: Index<usize, Output=Vec3<R>> + ?Sized,
WF: Index<usize, Output=Vec3<R>> + IndexMut<usize> + ?Sized,
{
let dim = meta.dim();
let stim_e_matrix = DimSlice::new(dim, stim_e);
let mat_matrix = DimSlice::new(dim, mat);
let mut e_matrix = DimSlice::new(dim, e);
let h_matrix = DimSlice::new(dim, h);
let stim_e = stim_e_matrix[idx];
let mat = &mat_matrix[idx];
let in_e = e_matrix[idx];
let in_h = VolumeSampleNeg::from_indexable(&h_matrix, idx);
let update_state = StepEContext {
inv_feature_size: meta.inv_feature_size(),
time_step: meta.time_step(),
stim_e,
mat,
in_h,
in_e,
};
let new_e = update_state.step_e();
e_matrix[idx] = new_e;
}
pub fn step_e(self) -> Vec3<R> {
let twice_eps0 = R::twice_eps0();
let deltas = self.in_h.delta_h();
// \nabla x H
let nabla_h = deltas.nabla() * self.inv_feature_size;
// $\nabla x H = \epsilon_0 dE/dt + \sigma E$
// no-conductivity version:
// let delta_e = nabla_h * (self.time_step * EPS0_INV);
let sigma = self.mat.conductivity();
let e_prev = self.in_e;
let delta_e = (nabla_h - e_prev.elem_mul(sigma)).elem_div(
sigma*self.time_step + Vec3::uniform(twice_eps0)
)*(R::two()*self.time_step);
// println!("spirv-step_e delta_e: {:?}", delta_e);
e_prev + delta_e + self.stim_e
}
}
pub struct StepHContext<'a, R, M> {
pub inv_feature_size: R,
pub time_step: R,
pub stim_h: Vec3<R>,
pub mat: &'a M,
/// Input field sampled near this location
pub in_e: VolumeSamplePos<R>,
pub in_h: Vec3<R>,
pub in_m: Vec3<R>,
}
impl<'a, R: Real, M: Material<R>> StepHContext<'a, R, M> {
pub fn step_flat_view<RM, RF, WF>(
meta: SimMeta<R>,
mat: &RM,
stim_h: &RF,
e: &RF,
h: &mut WF,
m: &mut WF,
idx: Vec3u,
)
where
RM: Index<usize, Output=M> + ?Sized,
RF: Index<usize, Output=Vec3<R>> + ?Sized,
WF: Index<usize, Output=Vec3<R>> + IndexMut<usize> + ?Sized,
{
let dim = meta.dim();
let stim_h_matrix = DimSlice::new(dim, stim_h);
let mat_matrix = DimSlice::new(dim, mat);
let e_matrix = DimSlice::new(dim, e);
let mut h_matrix = DimSlice::new(dim, h);
let mut m_matrix = DimSlice::new(dim, m);
let stim_h = stim_h_matrix[idx];
let mat = &mat_matrix[idx];
let in_e = VolumeSamplePos::from_indexable(&e_matrix, dim, idx);
let in_h = h_matrix[idx];
let in_m = m_matrix[idx];
let update_state = StepHContext {
inv_feature_size: meta.inv_feature_size(),
time_step: meta.time_step(),
stim_h,
mat,
in_e,
in_h,
in_m,
};
let (new_h, new_m) = update_state.step_h();
h_matrix[idx] = new_h;
m_matrix[idx] = new_m;
}
pub fn step_h(self) -> (Vec3<R>, Vec3<R>) {
let mu0 = R::mu0();
let mu0_inv = R::mu0_inv();
let deltas = self.in_e.delta_e();
// println!("spirv-step_h delta_e_struct: {:?}", deltas);
// \nabla x E
let nabla_e = deltas.nabla() * self.inv_feature_size;
// println!("spirv-step_h nabla_e: {:?}", nabla_e);
let delta_b = nabla_e * (-self.time_step);
// Relation between these is: B = mu0*(H + M)
let old_h = self.in_h;
let old_m = self.in_m;
let old_b = (old_h + old_m) * mu0;
let new_b = old_b + delta_b + self.stim_h * mu0;
let mat = self.mat;
let new_m = mat.move_b_vec(old_m, new_b);
let new_h = new_b * mu0_inv - new_m;
// println!("spirv-step_h delta_h: {:?}", delta_h);
(new_h, new_m)
}
}

View File

@@ -1,14 +1,13 @@
use core::ops::{Index, IndexMut};
use core::ops::Index;
use crate::compound::Optional;
use crate::dim::DimSlice;
use crate::mat::Material;
use crate::real::Real;
use crate::vec::{Vec3, Vec3u};
#[cfg(feature = "serde")]
use serde::{Serialize, Deserialize};
#[cfg_attr(feature = "serde", derive(Deserialize, Serialize))]
#[cfg_attr(feature = "fmt", derive(Debug))]
#[derive(Copy, Clone, Default, PartialEq)]
@@ -102,7 +101,7 @@ fn prev_z(idx: Vec3u) -> Optional<Vec3u> {
impl<R: Real> VolumeSampleNeg<R> {
/// Calculate the delta in H values amongst this cell and its neighbors (left/up/out)
fn delta_h(self) -> FieldDeltas<R> {
pub fn delta_h(self) -> FieldDeltas<R> {
let mid = self.mid;
// let (dfy_dx, dfz_dx) = self.xm1.map(|xm1| {
// (mid.y() - xm1.y(), mid.z() - xm1.z())
@@ -189,7 +188,7 @@ fn next_z(dim: Vec3u, idx: Vec3u) -> Optional<Vec3u> {
impl<R: Real> VolumeSamplePos<R> {
/// Calculate the delta in E values amongst this cell and its neighbors (right/down/in)
fn delta_e(self) -> FieldDeltas<R> {
pub fn delta_e(self) -> FieldDeltas<R> {
let mid = self.mid;
// let (dfy_dx, dfz_dx) = self.xp1.map(|xp1| {
// (xp1.y() - mid.y(), xp1.z() - mid.z())
@@ -232,7 +231,7 @@ impl<R: Real> VolumeSamplePos<R> {
}
}
struct FieldDeltas<R> {
pub struct FieldDeltas<R> {
dfy_dx: R,
dfz_dx: R,
dfx_dy: R,
@@ -242,7 +241,7 @@ struct FieldDeltas<R> {
}
impl<R: Real> FieldDeltas<R> {
fn nabla(self) -> Vec3<R> {
pub fn nabla(self) -> Vec3<R> {
Vec3::new(
self.dfz_dy - self.dfy_dz,
self.dfx_dz - self.dfz_dx,
@@ -251,145 +250,3 @@ impl<R: Real> FieldDeltas<R> {
}
}
pub struct StepEContext<'a, R, M> {
pub inv_feature_size: R,
pub time_step: R,
pub stim_e: Vec3<R>,
pub mat: &'a M,
/// Input field sampled near this location
pub in_h: VolumeSampleNeg<R>,
pub in_e: Vec3<R>,
}
impl<'a, R: Real, M: Material<R>> StepEContext<'a, R, M> {
pub fn step_flat_view<RM, RF, WF>(
meta: SimMeta<R>,
mat: &RM,
stim_e: &RF,
e: &mut WF,
h: &RF,
idx: Vec3u,
)
where
RM: Index<usize, Output=M> + ?Sized,
RF: Index<usize, Output=Vec3<R>> + ?Sized,
WF: Index<usize, Output=Vec3<R>> + IndexMut<usize> + ?Sized,
{
let dim = meta.dim;
let stim_e_matrix = DimSlice::new(dim, stim_e);
let mat_matrix = DimSlice::new(dim, mat);
let mut e_matrix = DimSlice::new(dim, e);
let h_matrix = DimSlice::new(dim, h);
let stim_e = stim_e_matrix[idx];
let mat = &mat_matrix[idx];
let in_e = e_matrix[idx];
let in_h = VolumeSampleNeg::from_indexable(&h_matrix, idx);
let update_state = StepEContext {
inv_feature_size: meta.inv_feature_size,
time_step: meta.time_step,
stim_e,
mat,
in_h,
in_e,
};
let new_e = update_state.step_e();
e_matrix[idx] = new_e;
}
pub fn step_e(self) -> Vec3<R> {
let twice_eps0 = R::twice_eps0();
let deltas = self.in_h.delta_h();
// \nabla x H
let nabla_h = deltas.nabla() * self.inv_feature_size;
// $\nabla x H = \epsilon_0 dE/dt + \sigma E$
// no-conductivity version:
// let delta_e = nabla_h * (self.time_step * EPS0_INV);
let sigma = self.mat.conductivity();
let e_prev = self.in_e;
let delta_e = (nabla_h - e_prev.elem_mul(sigma)).elem_div(
sigma*self.time_step + Vec3::uniform(twice_eps0)
)*(R::two()*self.time_step);
// println!("spirv-step_e delta_e: {:?}", delta_e);
e_prev + delta_e + self.stim_e
}
}
pub struct StepHContext<'a, R, M> {
pub inv_feature_size: R,
pub time_step: R,
pub stim_h: Vec3<R>,
pub mat: &'a M,
/// Input field sampled near this location
pub in_e: VolumeSamplePos<R>,
pub in_h: Vec3<R>,
pub in_m: Vec3<R>,
}
impl<'a, R: Real, M: Material<R>> StepHContext<'a, R, M> {
pub fn step_flat_view<RM, RF, WF>(
meta: SimMeta<R>,
mat: &RM,
stim_h: &RF,
e: &RF,
h: &mut WF,
m: &mut WF,
idx: Vec3u,
)
where
RM: Index<usize, Output=M> + ?Sized,
RF: Index<usize, Output=Vec3<R>> + ?Sized,
WF: Index<usize, Output=Vec3<R>> + IndexMut<usize> + ?Sized,
{
let dim = meta.dim;
let stim_h_matrix = DimSlice::new(dim, stim_h);
let mat_matrix = DimSlice::new(dim, mat);
let e_matrix = DimSlice::new(dim, e);
let mut h_matrix = DimSlice::new(dim, h);
let mut m_matrix = DimSlice::new(dim, m);
let stim_h = stim_h_matrix[idx];
let mat = &mat_matrix[idx];
let in_e = VolumeSamplePos::from_indexable(&e_matrix, dim, idx);
let in_h = h_matrix[idx];
let in_m = m_matrix[idx];
let update_state = StepHContext {
inv_feature_size: meta.inv_feature_size,
time_step: meta.time_step,
stim_h,
mat,
in_e,
in_h,
in_m,
};
let (new_h, new_m) = update_state.step_h();
h_matrix[idx] = new_h;
m_matrix[idx] = new_m;
}
pub fn step_h(self) -> (Vec3<R>, Vec3<R>) {
let mu0 = R::mu0();
let mu0_inv = R::mu0_inv();
let deltas = self.in_e.delta_e();
// println!("spirv-step_h delta_e_struct: {:?}", deltas);
// \nabla x E
let nabla_e = deltas.nabla() * self.inv_feature_size;
// println!("spirv-step_h nabla_e: {:?}", nabla_e);
let delta_b = nabla_e * (-self.time_step);
// Relation between these is: B = mu0*(H + M)
let old_h = self.in_h;
let old_m = self.in_m;
let old_b = (old_h + old_m) * mu0;
let new_b = old_b + delta_b + self.stim_h * mu0;
let mat = self.mat;
let new_m = mat.move_b_vec(old_m, new_b);
let new_h = new_b * mu0_inv - new_m;
// println!("spirv-step_h delta_h: {:?}", delta_h);
(new_h, new_m)
}
}