reefer 0.3.0

Optimizing proc-macro for geometric algebra
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155

// Direct port of example_timespace_lorentz.html using Reefer's spacetime algebra.

// Speed of light in kilometres per second (matching the ganja.js example units)
const LIGHT_SPEED_KM_PER_S: f64 = 299_792.458;

#[allow(unused_macros)]
#[reefer::algebra(f64, 1, 3)]
mod spacetime {
    // Create a timespace Algebra with 1,3 metric.
    basis!(t = P0); // time-like
    basis!(e1 = N0); // x
    basis!(e2 = N1); // y
    basis!(e3 = N2); // z

    #[derive(Debug, Clone, Copy)]
    shape!(Scalar { 1 });
    #[derive(Debug, Clone, Copy)]
    shape!(Event { t, e1, e2, e3 });

    // Convenience helpers to mirror the JS example
    impl Event {
        pub fn new(t: f64, e1: f64, e2: f64, e3: f64) -> Self {
            Self { t, e1, e2, e3 }
        }

        pub fn components(&self) -> (f64, f64, f64, f64) {
            (self.t, self.e1, self.e2, self.e3)
        }

        #[allow(dead_code)]
        pub fn time(&self) -> f64 {
            self.t
        }
    }
}

impl core::fmt::Display for spacetime::Event {
    fn fmt(&self, f: &mut core::fmt::Formatter<'_>) -> core::fmt::Result {
        let (t, x, y, z) = self.components();
        write!(
            f,
            "Event {{ t: {:.9} μs, x: {:.9} km, y: {:.9} km, z: {:.9} km }}",
            t * 1.0e6,
            x * LIGHT_SPEED_KM_PER_S,
            y * LIGHT_SPEED_KM_PER_S,
            z * LIGHT_SPEED_KM_PER_S,
        )
    }
}

fn micros(x: f64) -> f64 {
    // The spacetime unit we use is lightseconds (for both time and space)
    x * 1.0e-6 // microseconds to lightseconds
}

fn km(x: f64) -> f64 {
    // kilometres to lightseconds for each axis
    x / LIGHT_SPEED_KM_PER_S
}

fn lorentz_gamma(beta: f64) -> f64 {
    1.0 / (1.0 - beta * beta).sqrt()
}

fn analytic_boost(event: &spacetime::Event, beta: f64) -> spacetime::Event {
    let gamma = lorentz_gamma(beta);
    let (t, x, y, z) = event.components();
    let t_prime = gamma * (t - beta * x);
    let x_prime = gamma * (x - beta * t);
    spacetime::Event::new(t_prime, x_prime, y, z)
}

fn interval_squared(event: &spacetime::Event) -> f64 {
    let (t, x, y, z) = event.components();
    t * t - x * x - y * y - z * z
}

fn assert_event_close(
    actual: &spacetime::Event,
    expected: &spacetime::Event,
    label: &str,
    tolerance: f64,
) {
    let (at, ax, ay, az) = actual.components();
    let (et, ex, ey, ez) = expected.components();
    assert!(
        (at - et).abs() <= tolerance,
        "time mismatch for {label}: expected {et}, got {at}"
    );
    assert!(
        (ax - ex).abs() <= tolerance,
        "x mismatch for {label}: expected {ex}, got {ax}"
    );
    assert!(
        (ay - ey).abs() <= tolerance,
        "y mismatch for {label}: expected {ey}, got {ay}"
    );
    assert!(
        (az - ez).abs() <= tolerance,
        "z mismatch for {label}: expected {ez}, got {az}"
    );
}

#[test]
fn test_lorentz_transform() {
    // Define two events in the earth's reference frame.
    let strike_tree = spacetime::Event::new(micros(10.0), km(0.0), 0.0, 0.0);
    let strike_pole = spacetime::Event::new(micros(10.0), km(20.0), 0.0, 0.0);
    let reference_event = spacetime::Event::new(1.0, 1.0, 0.0, 0.0);

    // Calculate the Lorentz transform to go to the rocket's frame via a rotor.
    let rocket_speed: f64 = 0.5; // Half the speed of light.
    let phi = rocket_speed.atanh();
    let a = -0.5 * phi;

    // have to explicitly pass everything into the expr! macro closure as a geometric object for now
    // TODO: still need to figure out what's allowed with the hygiene rules, tricky
    let boost_expr =
        spacetime::expr!(|ev: spacetime::Event, a: spacetime::Scalar| (a * e1t).exp().sandwich(ev));
    let boost = |ev: spacetime::Event| boost_expr(ev, spacetime::Scalar { _1: a });

    // Transform the events to the rocket's frame.
    let strike_tree_r = boost(strike_tree);
    let strike_pole_r = boost(strike_pole);
    let reference_event_r = boost(reference_event);

    // Analytic expectations using the standard Lorentz transform.
    let expected_tree = analytic_boost(&strike_tree, rocket_speed);
    let expected_pole = analytic_boost(&strike_pole, rocket_speed);
    let expected_reference = analytic_boost(&reference_event, rocket_speed);

    let tol = 1.0e-12;
    assert_event_close(&strike_tree_r, &expected_tree, "strike_tree", tol);
    assert_event_close(&strike_pole_r, &expected_pole, "strike_pole", tol);
    assert_event_close(
        &reference_event_r,
        &expected_reference,
        "reference_event",
        tol,
    );

    // Spacetime interval should remain invariant under the boost.
    for (label, original, boosted) in [
        ("strike_tree", strike_tree, strike_tree_r),
        ("strike_pole", strike_pole, strike_pole_r),
        ("reference_event", reference_event, reference_event_r),
    ] {
        let original_interval = interval_squared(&original);
        let boosted_interval = interval_squared(&boosted);
        assert!(
            (original_interval - boosted_interval).abs() <= tol,
            "Minkowski interval mismatch for {label}: expected {original_interval}, got {boosted_interval}"
        );
    }
}