From a9fc0894e059bfc2626fd43b846ce43814320a0c Mon Sep 17 00:00:00 2001
From: Colin <colin@uninsane.org>
Date: Tue, 15 Dec 2020 23:54:33 -0800
Subject: [PATCH] MagneticLoop tracks H as well

---
 src/meas.rs | 52 ++++++++++++++++++++++++++++++++++++++++------------
 1 file changed, 40 insertions(+), 12 deletions(-)

diff --git a/src/meas.rs b/src/meas.rs
index b87ac13..c7f5625 100644
--- a/src/meas.rs
+++ b/src/meas.rs
@@ -95,16 +95,30 @@ impl std::iter::Sum for FieldSample {
 }
 
 #[derive(Default)]
-struct FieldSamplePair(FieldSample, FieldSample);
+struct FieldSamples<T>(T);
 
-impl std::iter::Sum for FieldSamplePair {
+impl std::iter::Sum for FieldSamples<[FieldSample; 2]> {
     fn sum<I>(iter: I) -> Self
         where I: Iterator<Item = Self>
     {
-        let mut s = FieldSamplePair::default();
+        let mut s = Self::default();
         for p in iter {
-            s.0 = FieldSample(s.0.0 + p.0.0, s.0.1 + p.0.1, s.0.2 + p.0.2);
-            s.1 = FieldSample(s.1.0 + p.1.0, s.1.1 + p.1.1, s.1.2 + p.1.2);
+            s.0[0] = FieldSample(s.0[0].0 + p.0[0].0, s.0[0].1 + p.0[0].1, s.0[0].2 + p.0[0].2);
+            s.0[1] = FieldSample(s.0[1].0 + p.0[1].0, s.0[1].1 + p.0[1].1, s.0[1].2 + p.0[1].2);
+        }
+        s
+    }
+}
+
+impl std::iter::Sum for FieldSamples<[FieldSample; 3]> {
+    fn sum<I>(iter: I) -> Self
+        where I: Iterator<Item = Self>
+    {
+        let mut s = Self::default();
+        for p in iter {
+            s.0[0] = FieldSample(s.0[0].0 + p.0[0].0, s.0[0].1 + p.0[0].1, s.0[0].2 + p.0[0].2);
+            s.0[1] = FieldSample(s.0[1].0 + p.0[1].0, s.0[1].1 + p.0[1].1, s.0[1].2 + p.0[1].2);
+            s.0[2] = FieldSample(s.0[2].0 + p.0[2].0, s.0[2].1 + p.0[2].1, s.0[2].2 + p.0[2].2);
         }
         s
     }
@@ -182,26 +196,40 @@ impl MagneticLoop {
 #[typetag::serde]
 impl AbstractMeasurement for MagneticLoop {
     fn eval(&self, state: &dyn GenericSim) -> String {
-        let FieldSamplePair(FieldSample(volume, directed_m, _m_vec), FieldSample(_, directed_b, _b_vec)) =
-            sum_over_region(state, &self.region, |coord, cell|
+        let FieldSamples([
+            FieldSample(volume, directed_m, _m_vec),
+            FieldSample(_, directed_b, _b_vec),
+            FieldSample(_, directed_h, _h_vec),
+        ]) = sum_over_region(state, &self.region, |coord, cell|
         {
             let normal = self.region.axis();
             let to_coord = *coord - *self.region.center();
             let tangent = normal.cross(to_coord).norm();
+
             let m = cell.mat().m();
             let directed_m = m.dot(tangent);
             let b = cell.b();
             let directed_b = b.dot(tangent);
-            FieldSamplePair(FieldSample(1, directed_m, m), FieldSample(1, directed_b, b))
+            let h = cell.h();
+            let directed_h = h.dot(tangent);
+            FieldSamples([
+                FieldSample(1, directed_m, m),
+                FieldSample(1, directed_b, b),
+                FieldSample(1, directed_h, h),
+            ])
         });
-        let mean_directed_m = directed_m / (volume as Flt);
-        let mean_directed_b = directed_b / (volume as Flt);
         let cross_section = self.region.cross_section() / (state.feature_size() * state.feature_size());
+        let mean_directed_m = directed_m / (volume as Flt);
         let cross_sectional_m = mean_directed_m * cross_section;
+        let mean_directed_b = directed_b / (volume as Flt);
         let cross_sectional_b = mean_directed_b * cross_section;
+        let mean_directed_h = directed_h / (volume as Flt);
+        let cross_sectional_h = mean_directed_h * cross_section;
         format!(
-            "M({}): {:.2e}; B({}): {:.2e}",
-            self.name, cross_sectional_m, self.name, cross_sectional_b
+            "M({}): {:.2e}; B({}): {:.2e}; H({}): {:.2e}",
+            self.name, cross_sectional_m,
+            self.name, cross_sectional_b,
+            self.name, cross_sectional_h,
         )
     }
 }