cmtool_core/coordinates/
vec3.rs

1// SPDX-License-Identifier: GPL-3.0-or-later
2
3use crate::coordinates::{CartesianCoordinates, Coords3};
4
5#[derive(Clone, Copy, Debug)]
6pub struct CartesianVec3(pub Coords3);
7
8#[derive(Clone, Copy, Debug)]
9pub struct CylindricalVec3(pub Coords3, pub f64);
10
11pub trait Vec3 {
12    fn sub(&self, other: &Self) -> Self;
13    fn add(&self, other: &Self) -> Self;
14    fn cross(&self, other: &Self) -> Self;
15    fn dot(&self, other: &Self) -> f64;
16    fn norm(&self) -> f64;
17    fn normalized(&self) -> Self;
18    fn scale(&self, lambda: f64) -> Self;
19}
20
21impl CartesianVec3 {
22    pub fn from_point(a: CartesianCoordinates, b: CartesianCoordinates) -> Self {
23        Self([b.0[0] - a.0[0], b.0[1] - a.0[1], b.0[2] - a.0[2]])
24    }
25
26    pub fn from_point_origin(a: CartesianCoordinates) -> Self {
27        Self(a.0)
28    }
29}
30
31impl Vec3 for CartesianVec3 {
32    fn add(&self, other: &Self) -> Self {
33        Self([
34            self.0[0] + other.0[0],
35            self.0[1] + other.0[1],
36            self.0[2] + other.0[2],
37        ])
38    }
39
40    fn sub(&self, other: &Self) -> Self {
41        Self([
42            self.0[0] - other.0[0],
43            self.0[1] - other.0[1],
44            self.0[2] - other.0[2],
45        ])
46    }
47
48    fn scale(&self, lambda: f64) -> Self {
49        Self(self.0.map(|v| v * lambda))
50    }
51
52    fn cross(&self, other: &Self) -> Self {
53        Self([
54            self.0[1] * other.0[2] - self.0[2] * other.0[1],
55            self.0[2] * other.0[0] - self.0[0] * other.0[2],
56            self.0[0] * other.0[1] - self.0[1] * other.0[0],
57        ])
58    }
59
60    fn dot(&self, other: &Self) -> f64 {
61        self.0[0] * other.0[0] + self.0[1] * other.0[1] + self.0[2] * other.0[2]
62    }
63
64    fn norm(&self) -> f64 {
65        f64::sqrt(self.0[0] * self.0[0] + self.0[1] * self.0[1] + self.0[2] * self.0[2])
66    }
67
68    fn normalized(&self) -> Self {
69        let _norm = self.norm();
70        self.scale(1. / _norm)
71    }
72}
73
74//https://en.wikipedia.org/wiki/Vector_fields_in_cylindrical_and_spherical_coordinates
75impl CylindricalVec3 {
76    pub fn from_cartesian_vec(value: CartesianVec3, base_theta: f64) -> Self {
77        let vx = value.0[0];
78        let vy = value.0[1];
79        let vz = value.0[2];
80
81        let vr = vx * base_theta.cos() + vy * base_theta.sin();
82        let vtheta = -vx * base_theta.sin() + vy * base_theta.cos();
83
84        CylindricalVec3([vr, vtheta, vz], base_theta)
85    }
86
87    pub fn to_cartesian_vec(&self) -> CartesianVec3 {
88        let vr = self.0[0];
89        let vtheta = self.0[1];
90        let vz = self.0[2];
91        let base_theta = self.1;
92
93        let vx = vr * base_theta.cos() - vtheta * base_theta.sin();
94        let vy = vr * base_theta.sin() + vtheta * base_theta.cos();
95
96        CartesianVec3([vx, vy, vz])
97    }
98}
99
100impl CartesianVec3 {
101    pub fn from_cylindrical_vec(value: CylindricalVec3) -> Self {
102        let vr = value.0[0];
103        let vtheta = value.0[1];
104        let vz = value.0[2];
105        let base_theta = value.1;
106
107        let vx = vr * base_theta.cos() - vtheta * base_theta.sin();
108        let vy = vr * base_theta.sin() + vtheta * base_theta.cos();
109
110        CartesianVec3([vx, vy, vz])
111    }
112
113    pub fn to_cylindrical_vec(&self, base_theta: f64) -> CylindricalVec3 {
114        let vx = self.0[0];
115        let vy = self.0[1];
116        let vz = self.0[2];
117
118        let vr = vx * base_theta.cos() + vy * base_theta.sin();
119        let vtheta = -vx * base_theta.sin() + vy * base_theta.cos();
120
121        CylindricalVec3([vr, vtheta, vz], base_theta)
122    }
123}
124
125#[cfg(test)]
126mod tests {
127    use super::*;
128    #[test]
129    fn test_cartesian_cylindrical_round_trip() {
130        let original = CartesianVec3([3.0, 4.0, 5.0]);
131        let cyl = original.to_cylindrical_vec(0.);
132        let converted = cyl.to_cartesian_vec();
133
134        for i in 0..3 {
135            assert!((converted.0[i] - original.0[i]).abs() < 1e-10);
136        }
137    }
138
139    #[test]
140    fn test_cartesian_cylindrical_1() {
141        let base_theta = 0.0;
142
143        let cart = CartesianVec3([3.0, 4.0, 5.0]);
144        let cyl = cart.to_cylindrical_vec(base_theta);
145
146        let vr_expected = 3.0;
147        let vtheta_expected = 4.0;
148        let vz_expected = 5.0;
149
150        assert!((cyl.0[0] - vr_expected).abs() < 1e-10);
151        assert!((cyl.0[1] - vtheta_expected).abs() < 1e-10);
152        assert!((cyl.0[2] - vz_expected).abs() < 1e-10);
153        assert!((cyl.1 - base_theta).abs() < 1e-10);
154
155        // Conversion inverse
156        let cart_converted = cyl.to_cartesian_vec();
157
158        for i in 0..3 {
159            assert!((cart_converted.0[i] - cart.0[i]).abs() < 1e-10);
160        }
161    }
162
163    #[test]
164    fn test_cartesian_to_cylindrical() {
165        let base_theta = std::f64::consts::FRAC_PI_4;
166
167        let cart = CartesianVec3([1.0, 0.0, 2.0]);
168        let cyl = cart.to_cylindrical_vec(base_theta);
169
170        let vr_expected = 1.0 * base_theta.cos() + 0.0 * base_theta.sin();
171        let vtheta_expected = -base_theta.sin() + 0.0 * base_theta.cos();
172        let vz_expected = 2.0;
173
174        assert!((cyl.0[0] - vr_expected).abs() < 1e-10);
175        assert!((cyl.0[1] - vtheta_expected).abs() < 1e-10);
176        assert!((cyl.0[2] - vz_expected).abs() < 1e-10);
177
178        // Puis conversion inverse
179        let cart_converted = cyl.to_cartesian_vec();
180
181        for i in 0..3 {
182            assert!((cart_converted.0[i] - cart.0[i]).abs() < 1e-10);
183        }
184    }
185
186    #[test]
187    fn test_cylindrical_to_cartesian_base_theta_0() {
188        let r = 2.0_f64.sqrt();
189        let theta = std::f64::consts::FRAC_PI_4;
190        let z = 3.0;
191        let base_theta = 0.0;
192
193        let cyl = CylindricalVec3([r, theta, z], base_theta);
194        let cart = cyl.to_cartesian_vec();
195        // vx = vr * cos(0) - vtheta * sin(0) = vr
196        // vy = vr * sin(0) + vtheta * cos(0) = vtheta
197        let x_expected = r;
198        let y_expected = theta;
199        let z_expected = 3.0;
200
201        assert!((cart.0[0] - x_expected).abs() < 1e-10);
202        assert!((cart.0[1] - y_expected).abs() < 1e-10);
203        assert!((cart.0[2] - z_expected).abs() < 1e-10);
204    }
205    #[test]
206    fn test_sub() {
207        let v1 = CartesianVec3([3.0, 4.0, 5.0]);
208        let v2 = CartesianVec3([1.0, 1.0, 1.0]);
209        let result = v1.sub(&v2);
210        assert_eq!(result.0, [2.0, 3.0, 4.0]);
211    }
212
213    #[test]
214    fn test_scale() {
215        let v = CartesianVec3([1.0, -2.0, 3.0]);
216        let result = v.scale(2.0);
217        assert_eq!(result.0, [2.0, -4.0, 6.0]);
218    }
219
220    #[test]
221    fn test_dot() {
222        let v1 = CartesianVec3([1.0, 2.0, 3.0]);
223        let v2 = CartesianVec3([4.0, -5.0, 6.0]);
224        let result = v1.dot(&v2);
225        assert_eq!(result, 12.0); // 1*4 + 2*(-5) + 3*6 = 4 -10 +18 = 12
226    }
227
228    #[test]
229    fn test_cross() {
230        let v1 = CartesianVec3([1.0, 0.0, 0.0]);
231        let v2 = CartesianVec3([0.0, 1.0, 0.0]);
232        let result = v1.cross(&v2);
233        assert_eq!(result.0, [0.0, 0.0, 1.0]);
234    }
235
236    #[test]
237    fn test_norm() {
238        let v = CartesianVec3([3.0, 4.0, 0.0]);
239        let result = v.norm();
240        assert!((result - 5.0).abs() < 1e-10);
241    }
242
243    #[allow(clippy::needless_range_loop)]
244    #[test]
245    fn test_normalized() {
246        let v = CartesianVec3([0.0, 3.0, 4.0]);
247        let normed = v.normalized();
248        let expected = [0.0, 0.6, 0.8];
249        for i in 0..3 {
250            assert!((normed.0[i] - expected[i]).abs() < 1e-10);
251        }
252    }
253}
cmtool_core/coordinates/vec3.rs

cmtool_core/coordinates/
vec3.rs
