diff --git a/examples/sr_latch.rs b/examples/sr_latch.rs index 279dba4..f266c45 100644 --- a/examples/sr_latch.rs +++ b/examples/sr_latch.rs @@ -23,7 +23,7 @@ fn main() { let m_to_um = |m: f32| (m * 1e6).round() as u32; let half_depth = depth * 0.5; - let base = "sr-latch-8-sr-longer-time-scale"; + let base = "sr-latch-9-higher-Hc"; let boundary_xy = 500e-6; let boundary_z = 300e-6; @@ -50,8 +50,16 @@ fn main() { let mut driver: Driver<_> = Driver::new_spirv(Meters::new(width, height, depth), feat_size); driver.set_steps_per_stim(1000); //driver.fill_region(&ferro1_region, mat::db::linear_iron()); - driver.fill_region(&ferro1_region, mat::MBFerromagnet::new(-0.3899, 0.3900, 310_000.0)); - driver.fill_region(&ferro2_region, mat::MBFerromagnet::new(-0.3899, 0.3900, 310_000.0)); + // Original, 3R1-LIKE ferromagnet (only a vague likeness), sr-latch-8: + // driver.fill_region(&ferro1_region, mat::MBFerromagnet::new(-0.3899, 0.3900, 310_000.0)); + // driver.fill_region(&ferro2_region, mat::MBFerromagnet::new(-0.3899, 0.3900, 310_000.0)); + // sr-latch-9; dead spot from B=[-0.03, 0.03]. This will help us see if the math is H-triggered + // or B-triggered + driver.fill_region(&ferro1_region, mat::MBFerromagnet::new(-0.3300, 0.3900, 310_000.0)); + driver.fill_region(&ferro2_region, mat::MBFerromagnet::new(-0.3300, 0.3900, 310_000.0)); + // mu_r=881.33, starting at H=25 to H=75. + // driver.fill_region(&ferro1_region, mat::MBFerromagnet::new(3.14159e-5, 0.055376, 44066.71)); + // driver.fill_region(&ferro2_region, mat::MBFerromagnet::new(3.14159e-5, 0.055376, 44066.71)); driver.fill_region(&set_region, mat::db::conductor(drive_conductivity)); driver.fill_region(&reset_region, mat::db::conductor(drive_conductivity)); driver.fill_region(&coupling_region, mat::db::conductor(drive_conductivity)); @@ -102,7 +110,7 @@ fn main() { add_drive_pulse(&reset_region, 55.0*current_duration, current_duration, peak_stim1); add_drive_pulse(&set_region, 55.0*current_duration, current_duration, peak_stim1); - let duration = 240.0*current_duration; + let duration = 120.0*current_duration; driver.add_measurement(meas::CurrentLoop::new("coupling", coupling_region.clone())); driver.add_measurement(meas::Current::new("coupling", coupling_region.clone())); @@ -132,7 +140,7 @@ fn main() { let _ = std::fs::create_dir_all(&prefix); driver.add_state_file(&*format!("{}/state.bc", prefix), 9600); - driver.add_serializer_renderer(&*format!("{}/frame-", prefix), 36000); + // driver.add_serializer_renderer(&*format!("{}/frame-", prefix), 36000); driver.add_csv_renderer(&*format!("{}/meas.csv", prefix), 200); driver.add_csv_renderer(&*format!("{}/meas-sparse.csv", prefix), 1600); diff --git a/src/sim/spirv_backend/src/mat.rs b/src/sim/spirv_backend/src/mat.rs index 186ca2d..517f4c0 100644 --- a/src/sim/spirv_backend/src/mat.rs +++ b/src/sim/spirv_backend/src/mat.rs @@ -55,6 +55,66 @@ impl MBPgram { } } +#[derive(Copy, Clone, Default)] +pub struct MHPgram { + /// X coordinate at which M is always zero. + h_intercept: f32, + /// relative mu value along the non-flat edges of the parallelogram. + mu_r: f32, + /// Vertical range of the graph + max_m: f32, +} + +impl MHPgram { + pub fn new(h_intercept: f32, mu_r: f32, max_m: f32) -> Self { + Self { h_intercept, mu_r, max_m } + } + + /// Return the new `M` + pub fn move_b(self, m: f32, target_b: f32) -> f32 { + const MU0_INV: f32 = 795774.715025073; + let target_mh = target_b*MU0_INV; + // The point may exist outside the parallelogram. + // The right slope is defined by: + // B(H)/mu0 = h_intercept + (h-h_intercept)*mu_r + // Left: + // B(H)/mu0 = -h_intercept + (h+h_intercept)*mu_r + let h_on_right = (target_mh-self.h_intercept)/self.mu_r + self.h_intercept; + let m_on_right = target_mh - h_on_right; + let h_on_left = (target_mh+self.h_intercept)/self.mu_r - self.h_intercept; + let m_on_left = target_mh - h_on_left; + + if m_on_right >= m { + if m_on_right <= self.max_m { + // rightward edge movement + m_on_right + } else { + // right of saturation + self.max_m + } + } else if m_on_left <= m { + if m_on_left >= -self.max_m { + // leftward edge movement + m_on_left + } else { + // left of saturation + -self.max_m + } + } else { + // interior movement + m + } + } + + pub fn move_b_vec(self, m: Vec3Std, target_b: Vec3Std) -> Vec3Std { + Vec3Std::new( + self.move_b(m.x(), target_b.x()), + self.move_b(m.y(), target_b.y()), + self.move_b(m.z(), target_b.z()), + ) + } +} + #[cfg(test)] mod test { use super::*; @@ -68,7 +128,7 @@ mod test { } #[test] - fn curve_interior() { + fn mb_curve_interior() { let curve = MBPgram::new(4.0, 6.0, 20.0f32); assert_eq_approx_tight(curve.move_b(0.0, 2.0), 0.0); @@ -82,7 +142,7 @@ mod test { } #[test] - fn curve_exterior() { + fn mb_curve_exterior() { let curve = MBPgram::new(4.0, 6.0, 20.0f32); assert_eq_approx_tight(curve.move_b(0.0, 6.0), 20.0); assert_eq_approx_tight(curve.move_b(0.0, 7.0), 20.0); @@ -100,7 +160,7 @@ mod test { } #[test] - fn curve_3r1() { + fn mb_curve_3r1() { // slope of 3r1 is about M=793210*B // This is almost identical to iron (795615)! let curve = MBPgram::new(-0.3899, 0.3900, 310000f32); @@ -130,4 +190,58 @@ mod test { // the magnetization has been completely erased: assert_eq_approx(curve.move_b(310000.0, -0.25), -198703.0, 1.0); } + + #[test] + fn mh_curve_edge_travel() { + const MU0: f32 = 1.2566370621219e-06; + let curve = MHPgram::new(50.0, 101.0, 500.0); + + assert_eq_approx(curve.move_b(0.0, 151.0*MU0), 100.0, 0.1); // rightward travel along edge + assert_eq_approx(curve.move_b(100.0, 252.0*MU0), 200.0, 0.1); + assert_eq_approx(curve.move_b(200.0, 555.0*MU0), 500.0, 0.1); + assert_eq_approx(curve.move_b(500.0, 500.0*MU0), 500.0, 0.1); // back (leftward) travel to H=0 + assert_eq_approx(curve.move_b(500.0, 455.0*MU0), 500.0, 0.1); // back (leftward) travel to H=-45 + assert_eq_approx(curve.move_b(500.0, 354.0*MU0), 400.0, 0.1); // back (leftward) travel to H=-46 + assert_eq_approx(curve.move_b(400.0, 253.0*MU0), 300.0, 0.1); // back (leftward) travel to H=-47 + assert_eq_approx(curve.move_b(300.0, -50.0*MU0), 0.0, 0.1); // back (leftward) travel to H=-50 + assert_eq_approx(curve.move_b(0.0, -151.0*MU0), -100.0, 0.1); // back (leftward) travel to H=-51 + assert_eq_approx(curve.move_b(-100.0, -555.0*MU0), -500.0, 0.1); // back (leftward) travel to H=-55 + assert_eq_approx(curve.move_b(-500.0, -456.0*MU0), -500.0, 0.1); // interior (rightward) travel to H=44 + assert_eq_approx(curve.move_b(-500.0, -354.0*MU0), -400.0, 0.1); // interior (rightward) travel to H=46 + assert_eq_approx(curve.move_b(-400.0, -253.0*MU0), -300.0, 0.1); // interior (rightward) travel to H=47 + assert_eq_approx(curve.move_b(-300.0, 50.0*MU0), 0.0, 0.1); // interior (rightward) travel to H=50 + } + + #[test] + fn mh_curve_interior() { + const MU0: f32 = 1.2566370621219e-06; + let curve = MHPgram::new(50.0, 101.0, 500.0); + // max M (right) happens at H=55; B = 555*MU0; + // min M (right) happens at H=45; B = -455*MU0; + // max M (left) happens at H=-45; B = 455*MU0; + // min M (left) happens at H=-55; B = -555*MU0; + + assert_eq_approx(curve.move_b(0.0, 40.0*MU0), 0.0, 0.1); + assert_eq_approx(curve.move_b(0.0, 50.0*MU0), 0.0, 0.1); + assert_eq_approx(curve.move_b(0.0, -40.0*MU0), 0.0, 0.1); + assert_eq_approx(curve.move_b(0.0, -50.0*MU0), 0.0, 0.1); + + assert_eq_approx(curve.move_b(-400.0, -400.0*MU0), -400.0, 0.1); // rightward travel from H=-54 to NOT H=-53.5 + assert_eq_approx(curve.move_b(-400.0, -355.0*MU0), -400.0, 0.1); // rightward travel from H=-54 to H=45 + assert_eq_approx(curve.move_b(-400.0, -354.0*MU0), -400.0, 0.1); // rightward travel from H=-54 to H=46 + assert_eq_approx(curve.move_b(-400.0, 5000.0*MU0), 500.0, 0.1); // rightward travel from H=-54 to H>>55 + } + + #[test] + fn mh_curve_exterior() { + const MU0: f32 = 1.2566370621219e-06; + let curve = MHPgram::new(50.0, 101.0, 500.0); + + assert_eq_approx(curve.move_b(500.0, 556.0*MU0), 500.0, 0.1); // exterior travel to H=56 + assert_eq_approx(curve.move_b(500.0, 5000.0*MU0), 500.0, 0.1); // exterior travel to H>>55 + + assert_eq_approx(curve.move_b(500.0, -5000.0*MU0), -500.0, 0.1); // exterior travel from H>>55 to H << -55 + assert_eq_approx(curve.move_b(-501.0, 454.0*MU0), 400.0, 0.1); // exterior travel from H<<-55 (rounding error) to H=54 + + } }