Tue, 22 Oct 2024 08:39:46 -0500
documentation
13 | 1 | /*! |
2 | Implementation of the surface of the 3D cube as a [`ManifoldPoint`]. | |
3 | */ | |
0 | 4 | |
12 | 5 | use serde_repr::*; |
6 | use serde::Serialize; | |
0 | 7 | use alg_tools::loc::Loc; |
8 | use alg_tools::norms::{Norm, L2}; | |
7 | 9 | use crate::manifold::{ManifoldPoint, EmbeddedManifoldPoint}; |
0 | 10 | |
13 | 11 | /// All the difference faces of a [`OnCube`]. |
12 | 12 | #[derive(Copy, Clone, Debug, Eq, PartialEq, Serialize_repr, Deserialize_repr)] |
13 | #[repr(u8)] | |
14 | pub enum Face {F1 = 1, F2 = 2, F3 = 3, F4 = 4, F5 = 5, F6 = 6} | |
0 | 15 | use Face::*; |
16 | ||
17 | pub type Point = Loc<f64, 2>; | |
18 | ||
19 | pub type AdjacentFaces = [Face; 4]; | |
20 | ||
12 | 21 | #[derive(Clone, Debug, Serialize)] |
0 | 22 | pub enum Path { |
23 | Direct { destination : Face }, | |
24 | Indirect { destination : Face, intermediate : Face }, | |
25 | } | |
26 | ||
27 | /// An iterator over paths on a cube, from a source face to a destination face. | |
28 | #[derive(Clone, Debug)] | |
29 | pub enum PathIter { | |
30 | Direct(Face), | |
31 | Indirect{ destination : Face, intermediate : AdjacentFaces, current : usize}, | |
32 | Exhausted, | |
33 | } | |
34 | ||
35 | impl std::iter::Iterator for PathIter { | |
36 | type Item = Path; | |
37 | ||
38 | fn next(&mut self) -> Option<Self::Item> { | |
39 | use PathIter::*; | |
40 | match self { | |
41 | &mut Exhausted => None, | |
42 | &mut Direct(destination) => { | |
43 | *self = Exhausted; | |
44 | Some(Path::Direct { destination }) | |
45 | }, | |
46 | &mut Indirect{destination, intermediate : ref i, ref mut current} => { | |
47 | if *current < i.len() { | |
48 | let intermediate = i[*current]; | |
49 | *current += 1; | |
50 | Some(Path::Indirect{ destination, intermediate }) | |
51 | } else { | |
52 | *self = Exhausted; | |
53 | None | |
54 | } | |
55 | } | |
56 | } | |
57 | } | |
58 | } | |
59 | ||
11 | 60 | impl std::fmt::Display for Face { |
61 | fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { | |
62 | let s = match *self { | |
63 | F1 => "F1", | |
64 | F2 => "F2", | |
65 | F3 => "F3", | |
66 | F4 => "F4", | |
67 | F5 => "F5", | |
68 | F6 => "F6", | |
69 | }; | |
70 | write!(f, "{}", s) | |
71 | } | |
72 | } | |
0 | 73 | |
74 | impl Face { | |
11 | 75 | /// Return an aray of all faces |
76 | pub fn all() -> [Face; 6] { | |
77 | [F1, F2, F3, F4, F5, F6] | |
78 | } | |
79 | ||
3 | 80 | /// Returns an array of the four faces adjacent to `self` in the |
81 | /// order [left, right, down, up] in the `self`-relative unfolding. | |
0 | 82 | pub fn adjacent_faces(&self) -> AdjacentFaces { |
83 | match *self { | |
3 | 84 | F1 => [F3, F2, F4, F5], |
85 | F2 => [F4, F5, F1, F6], | |
86 | F3 => [F5, F4, F1, F6], | |
87 | F4 => [F3, F2, F1, F6], | |
88 | F5 => [F2, F3, F1, F6], | |
89 | F6 => [F3, F2, F4, F5], | |
0 | 90 | } |
91 | } | |
92 | ||
93 | /// Returns the face opposing `self`. | |
94 | pub fn opposing_face(&self) -> Face { | |
95 | match *self { | |
96 | F1 => F6, | |
97 | F2 => F3, | |
98 | F3 => F2, | |
3 | 99 | F4 => F5, |
0 | 100 | F5 => F4, |
101 | F6 => F1, | |
102 | } | |
103 | } | |
104 | ||
105 | /// Converts a point on an adjacent face to the coordinate system of `self`. | |
106 | pub fn convert_adjacent(&self, adjacent : Face, p: &Point) -> Option<Point> { | |
107 | let Loc([x, y]) = *p; | |
108 | let mk = |x, y| Some(Loc([x, y])); | |
109 | match adjacent { | |
110 | F1 => match *self { | |
111 | F2 => mk(y, x - 1.0), | |
112 | F3 => mk(1.0 - y, -x), | |
113 | F4 => mk(x, -y), | |
114 | F5 => mk(1.0 - x, y - 1.0), | |
115 | F1 => mk(x, y), | |
116 | F6 => None, | |
117 | }, | |
118 | F2 => match *self { | |
119 | F1 => mk(y + 1.0, x), | |
120 | F4 => mk(x + 1.0, y), | |
121 | F5 => mk(x - 1.0, y), | |
122 | F6 => mk(2.0 - y, x), | |
123 | F2 => mk(x, y), | |
124 | F3 => None, | |
125 | }, | |
126 | F3 => match *self { | |
127 | F1 => mk(-y, 1.0 - x), | |
128 | F4 => mk(x - 1.0, y), | |
129 | F5 => mk(x + 1.0, y), | |
130 | F6 => mk(y - 1.0, 1.0 - x), | |
131 | F3 => mk(x, y), | |
132 | F2 => None, | |
133 | }, | |
134 | F4 => match *self { | |
135 | F1 => mk(x, -y), | |
136 | F2 => mk(x - 1.0, y), | |
137 | F3 => mk(x + 1.0, y), | |
138 | F6 => mk(x, y - 1.0), | |
139 | F4 => mk(x, y), | |
140 | F5 => None, | |
141 | }, | |
142 | F5 => match *self { | |
143 | F1 => mk(1.0 -x, y + 1.0), | |
144 | F2 => mk(x + 1.0, y), | |
145 | F3 => mk(x - 1.0, y), | |
146 | F6 => mk(1.0 -x, 2.0 - y), | |
147 | F5 => mk(x, y), | |
148 | F4 => None, | |
149 | }, | |
150 | F6 => match *self { | |
151 | F2 => mk(y, 2.0 - x), | |
152 | F3 => mk(1.0 - y, x + 1.0), | |
153 | F4 => mk(x, y + 1.0), | |
154 | F5 => mk(1.0 - x, 2.0 - y), | |
155 | F6 => mk(x, y), | |
156 | F1 => None, | |
157 | } | |
158 | } | |
159 | } | |
160 | ||
161 | /// Converts a point behind a path to the coordinate system of `self`. | |
162 | pub fn convert(&self, path : &Path, p: &Point) -> Point { | |
163 | use Path::*; | |
3 | 164 | //dbg!(*self, path); |
0 | 165 | match path { |
166 | &Direct{ destination : d} => self.convert_adjacent(d, p), | |
167 | &Indirect{ destination : d, intermediate : i } | |
3 | 168 | => {self.convert_adjacent(i, &i.convert_adjacent(d, p).unwrap())} |
0 | 169 | }.unwrap() |
170 | } | |
171 | ||
172 | ||
173 | /// Returns an iterator over all the paths from `self` to `other`. | |
174 | fn paths(&self, other : Face) -> PathIter { | |
3 | 175 | //dbg!(self, other); |
0 | 176 | if self.opposing_face() == other { |
177 | PathIter::Indirect { | |
178 | intermediate : self.adjacent_faces(), | |
179 | destination : other, | |
180 | current : 0 | |
181 | } | |
182 | } else { | |
183 | PathIter::Direct(other) | |
184 | } | |
185 | } | |
3 | 186 | |
13 | 187 | /// Indicates whether an unfolded point `p` is on this face, i.e., |
188 | /// has coordinates in [0,1]². | |
3 | 189 | pub fn is_in_face(&self, p: &Point) -> bool { |
190 | p.iter().map(|t| t.abs()).all(|t| 0.0 <= t && t <= 1.0) | |
191 | } | |
192 | ||
13 | 193 | /// Given an unfolded point `p`, possibly outside this face, finds |
194 | /// the edge, presented by an adjacent face, in whose direction it is. | |
195 | /// | |
196 | /// **TODO:** this does not correctly handle corners, i.e., when the point is not in | |
197 | /// the direction of an adjacent face. | |
3 | 198 | pub fn find_crossing(&self, p : &Point) -> Face { |
199 | let &Loc([x, y]) = p; | |
200 | use std::cmp::Ordering::*; | |
201 | let crossing = |t| match (0.0 <= t, t<=1.0) { | |
202 | (false, _) => Less, | |
203 | (_, false) => Greater, | |
204 | _ => Equal, | |
205 | }; | |
206 | ||
207 | // TODO: how to properly handle corners? Just throw an error? | |
208 | match (crossing(x), crossing(y)) { | |
209 | (Equal, Equal) => *self, | |
210 | (Less, _) => self.adjacent_faces()[0], | |
211 | (Greater, _) => self.adjacent_faces()[1], | |
212 | (Equal, Less) => self.adjacent_faces()[2], | |
213 | (Equal, Greater) => self.adjacent_faces()[3], | |
214 | } | |
215 | } | |
7 | 216 | |
217 | /// Get embedded 3D coordinates | |
218 | pub fn embedded_coords(&self, p : &Point) -> Loc<f64, 3> { | |
219 | let &Loc([x, y]) = p; | |
220 | Loc(match *self { | |
221 | F1 => [x, y, 0.0], | |
222 | F2 => [1.0, x, y], | |
223 | F3 => [0.0, 1.0-x, y], | |
224 | F4 => [x, 0.0, y], | |
225 | F5 => [1.0 - x, 1.0, y], | |
226 | F6 => [x, y, 1.0], | |
227 | }) | |
228 | } | |
0 | 229 | } |
230 | ||
12 | 231 | #[derive(Clone, Debug, PartialEq, Serialize)] |
0 | 232 | pub struct OnCube { |
233 | face : Face, | |
234 | point : Point, | |
235 | } | |
236 | ||
5 | 237 | impl OnCube { |
7 | 238 | /// Creates a new point on the cube, given a face and face-relative coordinates |
239 | /// in [0, 1]^2 | |
240 | pub fn new(face : Face, point : Point) -> Self { | |
241 | assert!(face.is_in_face(&point)); | |
242 | OnCube { face, point } | |
243 | } | |
244 | ||
5 | 245 | /// Calculates both the logarithmic map and distance to another point |
246 | fn log_dist(&self, other : &Self) -> (<Self as ManifoldPoint>::Tangent, f64) { | |
247 | let mut best_len = f64::INFINITY; | |
248 | let mut best_tan = Loc([0.0, 0.0]); | |
249 | for path in self.face.paths(other.face) { | |
250 | let tan = self.face.convert(&path, &other.point) - &self.point; | |
251 | let len = tan.norm(L2); | |
252 | if len < best_len { | |
253 | best_tan = tan; | |
254 | best_len = len; | |
255 | } | |
256 | } | |
257 | (best_tan, best_len) | |
258 | } | |
12 | 259 | |
13 | 260 | /// Returns the face of this point. |
12 | 261 | pub fn face(&self) -> Face { |
262 | self.face | |
263 | } | |
5 | 264 | } |
265 | ||
7 | 266 | |
267 | impl EmbeddedManifoldPoint for OnCube { | |
268 | type EmbeddedCoords = Loc<f64, 3>; | |
269 | ||
270 | /// Get embedded 3D coordinates | |
271 | fn embedded_coords(&self) -> Loc<f64, 3> { | |
272 | self.face.embedded_coords(&self.point) | |
273 | } | |
274 | } | |
275 | ||
0 | 276 | impl ManifoldPoint for OnCube { |
277 | type Tangent = Point; | |
278 | ||
8 | 279 | fn exp(self, tangent : &Self::Tangent) -> Self { |
3 | 280 | let mut face = self.face; |
281 | let mut point = self.point + tangent; | |
282 | loop { | |
283 | let next_face = face.find_crossing(&point); | |
284 | if next_face == face { | |
285 | break | |
286 | } | |
287 | point = next_face.convert_adjacent(face, &point).unwrap(); | |
288 | face = next_face; | |
289 | } | |
290 | OnCube { face, point } | |
0 | 291 | } |
292 | ||
293 | fn log(&self, other : &Self) -> Self::Tangent { | |
5 | 294 | self.log_dist(other).0 |
295 | } | |
296 | ||
297 | fn dist_to(&self, other : &Self) -> f64 { | |
298 | self.log_dist(other).1 | |
0 | 299 | } |
6
df9628092285
Add a zero function on manifolds
Tuomo Valkonen <tuomov@iki.fi>
parents:
5
diff
changeset
|
300 | |
df9628092285
Add a zero function on manifolds
Tuomo Valkonen <tuomov@iki.fi>
parents:
5
diff
changeset
|
301 | fn tangent_origin(&self) -> Self::Tangent { |
df9628092285
Add a zero function on manifolds
Tuomo Valkonen <tuomov@iki.fi>
parents:
5
diff
changeset
|
302 | Loc([0.0, 0.0]) |
df9628092285
Add a zero function on manifolds
Tuomo Valkonen <tuomov@iki.fi>
parents:
5
diff
changeset
|
303 | } |
0 | 304 | } |
305 | ||
306 | #[cfg(test)] | |
307 | mod tests { | |
308 | use super::*; | |
309 | ||
310 | #[test] | |
10 | 311 | fn convert_adjacent() { |
312 | let point = Loc([0.4, 0.6]); | |
313 | ||
314 | for f1 in [F1, F2, F3, F4, F5, F6] { | |
315 | for f2 in [F1, F2, F3, F4, F5, F6] { | |
316 | println!("{:?}-{:?}", f1, f2); | |
317 | match f1.convert_adjacent(f2, &point) { | |
318 | None => assert_eq!(f2.opposing_face(), f1), | |
319 | Some(q) => { | |
320 | match f2.convert_adjacent(f1, &q) { | |
321 | None => assert_eq!(f1.opposing_face(), f2), | |
322 | Some(p) => assert!((p-&point).norm(L2) < 1e-9), | |
323 | } | |
324 | } | |
325 | } | |
326 | } | |
327 | } | |
328 | } | |
329 | ||
330 | // This will fail, as different return path does not guarantee | |
331 | // that a point outside the face will be returned to its point of origin. | |
332 | // #[test] | |
333 | // fn convert_paths() { | |
334 | // let point = Loc([0.4, 0.6]); | |
335 | ||
336 | // for f1 in [F1, F2, F3, F4, F5, F6] { | |
337 | // for f2 in [F1, F2, F3, F4, F5, F6] { | |
338 | // for p1 in f2.paths(f1) { | |
339 | // for p2 in f1.paths(f2) { | |
340 | // println!("{:?}-{:?}; {:?} {:?}", f1, f2, p1, p2); | |
341 | // let v = &f2.convert(&p1, &point); | |
342 | // let q = f1.convert(&p2, v); | |
343 | // assert!((q-&point).norm(L2) < 1e-9, | |
344 | // "norm({}-{}) ≥ 1e-9 (dest {})", q, &point, &v); | |
345 | // } | |
346 | // } | |
347 | // } | |
348 | // } | |
349 | // } | |
350 | ||
351 | #[test] | |
352 | fn log_adjacent() { | |
0 | 353 | let p1 = OnCube{ face : F1, point : Loc([0.5, 0.5])}; |
354 | let p2 = OnCube{ face : F2, point : Loc([0.5, 0.5])}; | |
355 | ||
356 | assert_eq!(p1.log(&p2).norm(L2), 1.0); | |
357 | } | |
358 | ||
359 | #[test] | |
10 | 360 | fn log_opposing_equal() { |
0 | 361 | let p1 = OnCube{ face : F1, point : Loc([0.5, 0.5])}; |
362 | let p2 = OnCube{ face : F6, point : Loc([0.5, 0.5])}; | |
363 | ||
364 | assert_eq!(p1.log(&p2).norm(L2), 2.0); | |
365 | } | |
366 | ||
367 | #[test] | |
10 | 368 | fn log_opposing_unique_shortest() { |
0 | 369 | let p1 = OnCube{ face : F1, point : Loc([0.3, 0.25])}; |
370 | let p2 = OnCube{ face : F6, point : Loc([0.3, 0.25])}; | |
371 | ||
372 | assert_eq!(p1.log(&p2).norm(L2), 1.5); | |
373 | } | |
374 | } | |
375 |