Add hamilton function

4 files changed, 71 insertions, 38 deletions

diff --git a/src/main.rs b/src/main.rs
index ab88302..773ea72 100644
--- a/src/main.rs
+++ b/src/main.rs
@@ -1,44 +1,42 @@
 pub mod plot;
 mod simulation;
+use num_complex::Complex;
 use plot::Graph;
+use simulation::time_evolution::*;
+use simulation::*;
+
+struct StaticHamiltonion;
+
+impl MatrixGen for StaticHamiltonion {
+    fn gen(_t: f64) -> ndarray::Array2<c64> {
+        let mut b = unsafe { ndarray::Array2::<c64>::uninitialized((4, 4)) };
+        expm::expm(&simulation::hamiltonians::X0(), &mut b);
+        b
+    }
+}
 
 fn main() {
-    println!("Hello, world!");
-    plot::Gnuplot::plot3d(
-        &[
-            [
-                0.14285714285714296,
-                0.5714285714285717,
-                0.0,
-                0.5714285714285717,
-                0.5714285714285717,
-                0.0,
-            ],
-            [
-                0.11910591478146304,
-                0.5652016093730936,
-                -0.01707887122942119,
-                0.6464319082591804,
-                0.498191752207936,
-                0.0002,
-            ],
-            [
-                0.21910591478146304,
-                0.7652016093730936,
-                -0.01707887122942119,
-                0.6464319082591804,
-                0.498191752207936,
-                0.0002,
-            ],
-            [
-                0.01910591478146304,
-                0.8652016093730936,
-                -0.01707887122942119,
-                0.6464319082591804,
-                0.498191752207936,
-                0.0002,
-            ],
-        ],
-        1,
-    )
+    let dt = 1.0;
+    let iterations = 100;
+
+    let state = simulation::twoxtwo_level::TwoXTwoLevel::new(
+        Complex::new(0.0, 0.0),
+        Complex::new(0.0, 0.0),
+        Complex::new(0.0, 0.0),
+        Complex::new(0.0, 0.0),
+    );
+    let trace = Vec::new();
+    for i in 0..iterations {
+        state.evolve::<StaticHamiltonion>(i as f64 * dt);
+        let fibre = state.fibrate();
+        trace.push([
+            fibre[0],
+            fibre[1],
+            fibre[2],
+            0.0,
+            0.0,
+            i as f64 / iterations as f64,
+        ]);
+    }
+    plot::Gnuplot::plot3d(&trace, 1)
 }
diff --git a/src/simulation/hamiltonians.rs b/src/simulation/hamiltonians.rs
new file mode 100644
index 0000000..03d8d1f
--- /dev/null
+++ b/src/simulation/hamiltonians.rs
@@ -0,0 +1,31 @@
+use super::c64;
+use ndarray::prelude::*;
+
+pub fn X0() -> Array2<c64> {
+    array![
+        [
+            c64::new(1.0, 0.0),
+            c64::new(0.0, 0.0),
+            c64::new(0.0, 0.0),
+            c64::new(0.0, 0.0),
+        ],
+        [
+            c64::new(0.0, 0.0),
+            c64::new(1.0, 0.0),
+            c64::new(0.0, 0.0),
+            c64::new(0.0, 0.0),
+        ],
+        [
+            c64::new(0.0, 0.0),
+            c64::new(0.0, 0.0),
+            c64::new(-1.0, 0.0),
+            c64::new(0.0, 0.0),
+        ],
+        [
+            c64::new(0.0, 0.0),
+            c64::new(0.0, 0.0),
+            c64::new(0.0, 0.0),
+            c64::new(-1.0, 0.0),
+        ]
+    ]
+}
diff --git a/src/simulation/mod.rs b/src/simulation/mod.rs
index 641444e..c72f275 100644
--- a/src/simulation/mod.rs
+++ b/src/simulation/mod.rs
@@ -1,3 +1,4 @@
+pub mod hamiltonians;
 pub mod s_5_fibration;
 pub mod time_evolution;
 pub mod two_level;
diff --git a/src/simulation/twoxtwo_level.rs b/src/simulation/twoxtwo_level.rs
index 8a34167..62425c8 100644
--- a/src/simulation/twoxtwo_level.rs
+++ b/src/simulation/twoxtwo_level.rs
@@ -25,11 +25,14 @@ impl State for TwoXTwoLevel {
         for val in conjugate.iter_mut() {
             *val = val.conj();
         }
+        let v = conjugate.dot(&X3().dot(&self.state))[(0, 0)];
 
         vec![
             conjugate.dot(&X0().dot(&self.state))[(0, 0)].re,
             conjugate.dot(&X1().dot(&self.state))[(0, 0)].re,
             conjugate.dot(&X2().dot(&self.state))[(0, 0)].re,
+            v.re,
+            v.im,
         ]
     }