hoomd_geometry/shape/
twelvetwelve.rs1use serde::{Deserialize, Serialize};
7
8use crate::IsPointInside;
9use hoomd_manifold::{Hyperbolic, Minkowski};
10use hoomd_vector::Metric;
11use std::f64::consts::PI;
12
13#[derive(Clone, Debug, PartialEq, Serialize, Deserialize)]
20pub struct TwelveTwelve {}
21
22impl IsPointInside<Hyperbolic<3>> for TwelveTwelve {
23 #[inline]
40 fn is_point_inside(&self, point: &Hyperbolic<3>) -> bool {
41 TwelveTwelve::distance_to_boundary(point) >= 0.0
42 }
43}
44
45impl TwelveTwelve {
46 #[inline]
80 #[must_use]
81 pub fn distance_to_boundary(point: &Hyperbolic<3>) -> f64 {
82 let theta =
83 (point.coordinates()[1].atan2(point.coordinates()[0])).rem_euclid(PI / 6.0) - PI / 12.0;
84 let boost = (point.coordinates()[2]).acosh();
85 let (b_sinh, b_cosh) = (boost.sinh(), point.coordinates()[2]);
86 let xi = Self::CUSP_TO_EDGE;
87 let (xi_sinh, xi_cosh) = (xi.sinh(), xi.cosh());
88 let edge_as_diameter: Hyperbolic<3> =
90 Hyperbolic::<3>::from_minkowski_coordinates(Minkowski::from([
91 xi_cosh * b_sinh * (theta.cos()) - xi_sinh * b_cosh,
92 b_sinh * (theta.sin()),
93 -xi_sinh * b_sinh * (theta.cos()) + xi_cosh * b_cosh,
94 ]));
95 let flipped = Hyperbolic::<3>::from_minkowski_coordinates(Minkowski::from([
96 -edge_as_diameter.coordinates()[0],
97 edge_as_diameter.coordinates()[1],
98 edge_as_diameter.coordinates()[2],
99 ]));
100 let sign = -(edge_as_diameter.coordinates()[0]).signum();
101 sign * (edge_as_diameter.distance(&flipped)) / 2.0
102 }
103 #[inline]
105 #[must_use]
106 pub fn gamma(eta: f64, theta: f64, point: &[f64; 3]) -> [f64; 3] {
107 let (eta_sinh_squared, two_eta_sinh, theta_sin, theta_cos) = (
108 (eta.sinh()).powi(2),
109 (2.0 * eta).sinh(),
110 theta.sin(),
111 theta.cos(),
112 );
113 [
114 (2.0 * (eta_sinh_squared) * ((theta_cos).powi(2)) + 1.0) * point[0]
115 + ((2.0 * theta).sin()) * (eta_sinh_squared) * point[1]
116 + (two_eta_sinh) * (theta_cos) * point[2],
117 ((2.0 * theta).sin()) * (eta_sinh_squared) * point[0]
118 + (2.0 * (eta_sinh_squared) * ((theta_sin).powi(2)) + 1.0) * point[1]
119 + (two_eta_sinh) * (theta_sin) * point[2],
120 (two_eta_sinh) * (theta_cos) * point[0]
121 + (two_eta_sinh) * (theta_sin) * point[1]
122 + ((2.0 * eta).cosh()) * point[2],
123 ]
124 }
125 #[inline]
128 #[must_use]
129 pub fn reorient(theta: f64, point: &[f64; 3]) -> f64 {
130 let (q_u, q_v) = (point[0] / (1.0 + point[2]), point[1] / (1.0 + point[2]));
131 let alpha = (1.0 + (PI / 6.0).cos()).sqrt();
132 let beta = (theta.cos()) * ((2.0 * ((PI / 6.0).cos())).sqrt());
133 let gamma = (theta.sin()) * ((2.0 * ((PI / 6.0).cos())).sqrt());
134 let p_x = alpha + beta * q_u + gamma * q_v;
135 let p_y = beta * q_v - gamma * q_u;
136 -2.0 * (p_y.atan2(p_x))
137 }
138 pub const TWELVETWELVE: f64 = 3.325_771_782_117_242;
140 pub const CUSP_TO_EDGE: f64 = 1.991_652_391_049_437;
142 pub const EDGE_LENGTH: f64 = 3.983_304_782_098_874;
144}
145
146#[cfg(test)]
147mod tests {
148 use super::*;
149 use approxim::assert_relative_eq;
150 use hoomd_manifold::{Hyperbolic, HyperbolicDisk, Minkowski};
151 use rand::{SeedableRng, distr::Distribution, rngs::StdRng};
152 use std::ops::Not;
153
154 #[test]
155 fn boundary_distance() {
156 let e = Hyperbolic::<3>::from_polar_coordinates(1.0, 0.1);
158 let e_edge_distance = TwelveTwelve::distance_to_boundary(&e);
159 let e_edge_distance_numeric = 1.028_489_242_583_61;
160 assert_relative_eq!(e_edge_distance, e_edge_distance_numeric, epsilon = 1e-12);
161
162 let f = Hyperbolic::<3>::from_polar_coordinates(0.6, 0.2 + PI / 6.0);
163 let f_edge_distance = TwelveTwelve::distance_to_boundary(&f);
164 let f_edge_distance_numeric = 1.393_774_804_611_314;
165 assert_relative_eq!(f_edge_distance, f_edge_distance_numeric, epsilon = 1e-12);
166 }
167
168 #[test]
169 fn inside_is_inside() {
170 let twelve_twelve = TwelveTwelve {};
171 let r = TwelveTwelve::CUSP_TO_EDGE;
172 let mut rng = StdRng::seed_from_u64(239);
173 let disk = HyperbolicDisk {
174 disk_radius: r.try_into().expect("hard-coded positive number"),
175 point: Hyperbolic::<3>::from_minkowski_coordinates(Minkowski::from([0.0, 0.0, 1.0])),
176 };
177 let random_point: Hyperbolic<3> = disk.sample(&mut rng);
178 assert!(twelve_twelve.is_point_inside(&random_point));
179
180 let point_1 = Hyperbolic::<3>::from_polar_coordinates(1.99, PI / 12.0);
181 assert!(twelve_twelve.is_point_inside(&point_1));
182
183 let point_2 = Hyperbolic::<3>::from_polar_coordinates(3.32, PI / 6.0);
184 assert!(twelve_twelve.is_point_inside(&point_2));
185 }
186
187 #[test]
188 fn outside_is_outside() {
189 let twelve_twelve = TwelveTwelve {};
190 let point_1 = Hyperbolic::<3>::from_polar_coordinates(2.0, PI / 12.0);
191 assert!((twelve_twelve.is_point_inside(&point_1)).not());
192
193 let point_2 = Hyperbolic::<3>::from_polar_coordinates(3.33, PI / 6.0);
194 assert!((twelve_twelve.is_point_inside(&point_2)).not());
195 }
196}