Skip to main content

outram_foam_basic_lib/primitives/
vector.rs

1// SPDX-License-Identifier: GPL-3.0-only
2// Copyright (C) 2026 OUTRAM PARK contributors
3// Derived from OpenFOAM (www.openfoam.com)
4// Copyright (C) 2004-2023 OpenFOAM Foundation
5// Copyright (C) 2016-2023 OpenCFD Ltd.
6//
7// This file is part of OUTRAM PARK.
8//
9// OUTRAM PARK is free software: you can redistribute it and/or modify it
10// under the terms of the GNU General Public License as published by the
11// Free Software Foundation, either version 3 of the License, or (at your
12// option) any later version.
13//
14// OUTRAM PARK is distributed in the hope that it will be useful, but
15// WITHOUT ANY WARRANTY; without even the implied warranty of
16// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
17// General Public License for more details.
18//
19// You should have received a copy of the GNU General Public License along
20// with OUTRAM PARK.  If not, see <https://www.gnu.org/licenses/>.
21
22use std::ops::{Add, AddAssign, Div, DivAssign, Mul, MulAssign, Neg, Sub, SubAssign};
23
24/// 3-component vector. Maps to `Foam::vector` (`Foam::Vector<scalar>`).
25/// Component layout: x, y, z.
26#[derive(Debug, Clone, Copy, PartialEq, Default)]
27pub struct Vector3 {
28    pub x: f64,
29    pub y: f64,
30    pub z: f64,
31}
32
33impl Vector3 {
34    pub const ZERO: Self = Self {
35        x: 0.0,
36        y: 0.0,
37        z: 0.0,
38    };
39    pub const ONE: Self = Self {
40        x: 1.0,
41        y: 1.0,
42        z: 1.0,
43    };
44    pub const X: Self = Self {
45        x: 1.0,
46        y: 0.0,
47        z: 0.0,
48    };
49    pub const Y: Self = Self {
50        x: 0.0,
51        y: 1.0,
52        z: 0.0,
53    };
54    pub const Z: Self = Self {
55        x: 0.0,
56        y: 0.0,
57        z: 1.0,
58    };
59
60    #[inline]
61    pub fn new(x: f64, y: f64, z: f64) -> Self {
62        Self { x, y, z }
63    }
64
65    /// Squared magnitude: |v|² = x² + y² + z²
66    #[inline]
67    pub fn mag_sqr(self) -> f64 {
68        self.x * self.x + self.y * self.y + self.z * self.z
69    }
70
71    /// Magnitude: |v|
72    #[inline]
73    pub fn mag(self) -> f64 {
74        self.mag_sqr().sqrt()
75    }
76
77    /// Squared distance to another vector
78    #[inline]
79    pub fn dist_sqr(self, other: Self) -> f64 {
80        (other - self).mag_sqr()
81    }
82
83    /// Distance to another vector
84    #[inline]
85    pub fn dist(self, other: Self) -> f64 {
86        self.dist_sqr(other).sqrt()
87    }
88
89    /// Dot (inner) product. C++ `operator&(Vector, Vector)`.
90    #[inline]
91    pub fn dot(self, other: Self) -> f64 {
92        self.x * other.x + self.y * other.y + self.z * other.z
93    }
94
95    /// Cross product. C++ `operator^(Vector, Vector)`.
96    #[inline]
97    pub fn cross(self, other: Self) -> Self {
98        Self {
99            x: self.y * other.z - self.z * other.y,
100            y: self.z * other.x - self.x * other.z,
101            z: self.x * other.y - self.y * other.x,
102        }
103    }
104
105    /// Normalise to unit vector; returns zero if `|v| < tol`.
106    pub fn normalise(self, tol: f64) -> Self {
107        let s = self.mag();
108        if s < tol {
109            Self::ZERO
110        } else {
111            self / s
112        }
113    }
114
115    /// Remove the component collinear with `unit_vec`: `self - (self·unit) * unit`.
116    #[inline]
117    pub fn remove_collinear(self, unit_vec: Self) -> Self {
118        self - self.dot(unit_vec) * unit_vec
119    }
120
121    /// Linear interpolation: `(1-t)*a + t*b`.
122    #[inline]
123    pub fn lerp(a: Self, b: Self, t: f64) -> Self {
124        let ot = 1.0 - t;
125        Self {
126            x: ot * a.x + t * b.x,
127            y: ot * a.y + t * b.y,
128            z: ot * a.z + t * b.z,
129        }
130    }
131}
132
133// --- Standard arithmetic operators ---
134
135impl Neg for Vector3 {
136    type Output = Self;
137    #[inline]
138    fn neg(self) -> Self {
139        Self {
140            x: -self.x,
141            y: -self.y,
142            z: -self.z,
143        }
144    }
145}
146
147impl Add for Vector3 {
148    type Output = Self;
149    #[inline]
150    fn add(self, r: Self) -> Self {
151        Self {
152            x: self.x + r.x,
153            y: self.y + r.y,
154            z: self.z + r.z,
155        }
156    }
157}
158
159impl Sub for Vector3 {
160    type Output = Self;
161    #[inline]
162    fn sub(self, r: Self) -> Self {
163        Self {
164            x: self.x - r.x,
165            y: self.y - r.y,
166            z: self.z - r.z,
167        }
168    }
169}
170
171impl Mul<f64> for Vector3 {
172    type Output = Self;
173    #[inline]
174    fn mul(self, s: f64) -> Self {
175        Self {
176            x: self.x * s,
177            y: self.y * s,
178            z: self.z * s,
179        }
180    }
181}
182
183impl Mul<Vector3> for f64 {
184    type Output = Vector3;
185    #[inline]
186    fn mul(self, v: Vector3) -> Vector3 {
187        v * self
188    }
189}
190
191impl Div<f64> for Vector3 {
192    type Output = Self;
193    #[inline]
194    fn div(self, s: f64) -> Self {
195        Self {
196            x: self.x / s,
197            y: self.y / s,
198            z: self.z / s,
199        }
200    }
201}
202
203impl AddAssign for Vector3 {
204    #[inline]
205    fn add_assign(&mut self, r: Self) {
206        *self = *self + r;
207    }
208}
209
210impl SubAssign for Vector3 {
211    #[inline]
212    fn sub_assign(&mut self, r: Self) {
213        *self = *self - r;
214    }
215}
216
217impl MulAssign<f64> for Vector3 {
218    #[inline]
219    fn mul_assign(&mut self, s: f64) {
220        *self = *self * s;
221    }
222}
223
224impl DivAssign<f64> for Vector3 {
225    #[inline]
226    fn div_assign(&mut self, s: f64) {
227        *self = *self / s;
228    }
229}
230
231// --- Free functions mirroring OpenFOAM globals ---
232
233#[inline]
234pub fn mag_sqr(v: Vector3) -> f64 {
235    v.mag_sqr()
236}
237
238#[inline]
239pub fn mag(v: Vector3) -> f64 {
240    v.mag()
241}
242
243/// Dot product. C++ `operator&`.
244#[inline]
245pub fn dot(a: Vector3, b: Vector3) -> f64 {
246    a.dot(b)
247}
248
249/// Cross product. C++ `operator^`.
250#[inline]
251pub fn cross(a: Vector3, b: Vector3) -> Vector3 {
252    a.cross(b)
253}
254
255#[inline]
256pub fn lerp(a: Vector3, b: Vector3, t: f64) -> Vector3 {
257    Vector3::lerp(a, b, t)
258}
259
260#[cfg(test)]
261mod tests {
262    use super::*;
263
264    #[test]
265    fn dot_and_cross() {
266        let x = Vector3::X;
267        let y = Vector3::Y;
268        let z = Vector3::Z;
269        assert_eq!(dot(x, x), 1.0);
270        assert_eq!(dot(x, y), 0.0);
271        assert_eq!(cross(x, y), z);
272        assert_eq!(cross(y, z), x);
273        assert_eq!(cross(z, x), y);
274    }
275
276    #[test]
277    fn cross_orthogonal_to_both_inputs() {
278        // For arbitrary vectors a, b: cross(a,b)·a == 0 and cross(a,b)·b == 0
279        let a = Vector3::new(1.0, 2.0, 3.0);
280        let b = Vector3::new(4.0, 5.0, 6.0);
281        let c = cross(a, b);
282        assert!(dot(c, a).abs() < 1e-14, "cross(a,b)·a = {}", dot(c, a));
283        assert!(dot(c, b).abs() < 1e-14, "cross(a,b)·b = {}", dot(c, b));
284    }
285
286    #[test]
287    fn mag_and_normalise() {
288        let v = Vector3::new(3.0, 4.0, 0.0);
289        assert_eq!(v.mag(), 5.0);
290        let n = v.normalise(1e-15);
291        assert!((n.mag() - 1.0).abs() < 1e-14);
292    }
293
294    #[test]
295    fn lerp() {
296        let a = Vector3::new(0.0, 0.0, 0.0);
297        let b = Vector3::new(2.0, 4.0, 6.0);
298        let m = Vector3::lerp(a, b, 0.5);
299        assert_eq!(m, Vector3::new(1.0, 2.0, 3.0));
300    }
301
302    #[test]
303    fn remove_collinear() {
304        let v = Vector3::new(3.0, 4.0, 0.0);
305        let n = v.normalise(1e-15);
306        let r = v.remove_collinear(n);
307        assert!(r.mag() < 1e-14);
308    }
309}