|
|
| Line 1: |
Line 1: |
| == planegeometry.rs == | | == Billiard scales == |
| | To run this, you need Rust installed. |
| | # In any folder you want, create a directory <code>billiard-scales</code>. |
| | # Using the command line, run <code>cargo init</code> in the new directory <code>billiard-scales</code>. |
| | # Create or copy the following .rs files in <code>billiard-scales/src</code>. Copy the contents of the <code>Cargo.toml</code> provided into the new <code>Cargo.toml</code> |
| | # Run <code>cargo run --release</code>. |
| | === planegeometry.rs === |
| <syntaxhighlight lang="rs">// plane_geometry.rs | | <syntaxhighlight lang="rs">// plane_geometry.rs |
| // v0.0 by inthar | | // v0.0.0 by inthar |
| // Geometry of points and lines in the plane, implemented with rational numbers. The implementation is janky, though. | | // Geometry of points and lines in the plane, implemented with rational numbers. The implementation is janky, though. |
| // Can be used to enumerate ternary billiard scales. | | // Can be used to enumerate ternary billiard scales. |
| Line 637: |
Line 643: |
| } | | } |
| } | | } |
| | </syntaxhighlight> |
| | === lib.rs === |
| | <syntaxhighlight lang="rs"> |
| | // scaletheory.rs |
| | // Inthar's higher-rank scale theory library |
| | // v0.0.0 |
| | use std::collections::BTreeSet; |
| | use std::collections::HashSet; |
| | use num_rational::Rational64 as r64; |
| | mod plane_geometry; |
| | use plane_geometry::*; |
| | use std::cmp::{min, max}; |
| | |
| | type JIScaleProperty = fn(s: [r64]) -> bool; |
| | type CentsScaleProperty = fn(s: [f64]) -> bool; |
| | type AbstractScaleProperty = fn(s: &str) -> bool; |
| | type EquivalenceRelation<T> = fn(&T, &T) -> bool; |
|
| |
|
| #[cfg(test)]
| | pub fn are_rotationally_equivalent(scale1: &String, scale2: &String)-> bool{ |
| mod tests {
| | if scale1.len() != scale2.len() { |
| use num_rational::Rational64 as r64; | | return false; |
| use crate::plane_geometry::Slope;
| | } else { |
| use crate::plane_geometry::Point; | | let length = scale1.len(); |
| use crate::plane_geometry::Line;
| | for i in 0..length { |
| use crate::plane_geometry::intersection;
| | if slice_cyclic_string(scale1, i, length) == *scale2 { |
| use crate::plane_geometry::ConvexPolygon;
| | return true; |
| use crate::plane_geometry::PointLineConfiguration;
| | } |
| use crate::plane_geometry::projected_cube;
| | } |
|
| | false |
| // Tests for the struct `Slope`.
| |
| #[test]
| |
| fn slope_bad_slope() {
| |
| let bad_slope : Result<Slope, String> = Slope::new(0,0);
| |
| assert_eq!(bad_slope, | |
| Err("Attempted to create a Slope with both `numer` and `denom` equal to 0".to_string()));
| |
| } | | } |
|
| | } |
| #[test]
| | |
| fn slope_infinite_slope() {
| | // Assumes that equiv : T x T -> bool is an equivalence relation. |
| let infinite_slope : Result<Slope, String> = Slope::new(1,0);
| | pub fn equivalence_class_representatives<T>(set: Vec<T>, equiv: EquivalenceRelation<T>) -> Vec<T> |
| assert_eq!(infinite_slope, Ok(Slope::infinity()));
| | where T : Clone + PartialEq { |
| }
| | let mut result : Vec<T> = Vec::new(); |
|
| | for i in 0..set.len() { |
| #[test]
| | let mut found_class_equivalent = false; |
| fn slope_zero_slope() { | | 'a: for j in 0..result.len(){ // All vecs in result are nonempty. |
| let zero_slope : Result<Slope, String> = Slope::new(0,1);
| | if equiv(&result[j], &set[i]){ |
| assert_eq!(zero_slope, Ok(Slope::zero()));
| | found_class_equivalent = true; |
| } | | continue 'a; |
|
| | } |
| #[test]
| | } |
| fn slope_integer_slope() {
| | if found_class_equivalent == false { |
| let slope_four : Slope = Slope::integer(4i64); | | result.push(set[i].clone()); |
| assert_eq!(slope_four, Slope::raw(4,1)); | | } |
| }
| |
|
| |
| #[test]
| |
| fn slope_reduction() {
| |
| let slope_two : Result<Slope, String> = Slope::new(8,4);
| |
| assert_eq!(slope_two, Ok(Slope::raw(2,1)));
| |
| let slope_half : Result<Slope, String> = Slope::new(4,8);
| |
| assert_eq!(slope_half, Ok(Slope::raw(1,2))); | |
| }
| |
|
| |
| #[test]
| |
| fn slope_as_rational() {
| |
| let good_rational : Result<r64, String> = Slope::new(3,2).expect("Bad slope").as_rational();
| |
| assert_eq!(good_rational, Ok(r64::new(3,2)));
| |
| | |
| let bad_rational : Result<r64, String> = Slope::new(3,0).expect("Bad slope").as_rational();
| |
| assert_eq!(bad_rational,
| |
| Err("Attempted to convert a Slope with `denom` == 0 into a Rational64".to_string()));
| |
| } | | } |
| | result |
| | } |
| | |
| | fn gcd(m: i64, n: i64) -> i64 { |
| | let mut x = m; |
| | let mut y = n; |
| | while x != 0 && y != 0 { |
| | if x > y { |
| | x = x % y; |
| | } else { |
| | y = y % x; |
| | } |
| | } |
| | if y == 0 {x} else {y} |
| | } |
|
| |
|
| // Tests for the module `Point`.
| | fn gcd_u32(m: u32, n: u32) -> u32 { |
|
| | let mut x = m; |
| #[test]
| | let mut y = n; |
| fn point_create_point() {
| | while x != 0 && y != 0 { |
| let point : Point = Point::new(r64::new(2,1),r64::new(2,1));
| | if x > y { |
| assert_eq!(point, Point{ x: r64::new(2,1), y: r64::new(2,1)})
| | x = x % y; |
| }
| | } else { |
|
| | y = y % x; |
| #[test]
| | } |
| fn point_midpoint_vertical() {
| | } |
| let point1 : Point = Point::new(r64::new(2,1),r64::new(2,1));
| | if y == 0 {x} else {y} |
| let point2 : Point = Point::new(r64::new(2,1),r64::new(1,1));
| | } |
| assert_eq!(Point::midpoint(&point1, &point2), Point{ x: r64::new(2,1), y: r64::new(3,2)})
| |
| }
| |
|
| |
| #[test]
| |
| fn point_midpoint_nonvertical() {
| |
| let point1 : Point = Point::new(r64::new(1,1),r64::new(2,1));
| |
| let point2 : Point = Point::new(r64::new(2,1),r64::new(1,1));
| |
| assert_eq!(Point::midpoint(&point1, &point2), Point{ x: r64::new(3,2), y: r64::new(3,2)})
| |
| }
| |
|
| |
|
| #[test]
| | fn lcm(m: i64, n: i64) -> i64 { |
| fn point_midpoint_identical() {
| | m*n/gcd(m,n) |
| let point1 : Point = Point::new(r64::new(1,1),r64::new(1,1));
| | } |
| let point2 : Point = Point::new(r64::new(1,1),r64::new(1,1));
| |
| assert_eq!(Point::midpoint(&point1, &point2), Point{ x: r64::new(1,1), y: r64::new(1,1)})
| |
| }
| |
|
| |
| #[test]
| |
| fn point_collinear() {
| |
| let point1 : Point = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let point2 : Point = Point::new(r64::new(3,1),r64::new(4,1));
| |
| let point3 : Point = Point::new(r64::new(1,1),r64::new(4,1));
| |
| let point4 : Point = Point::new(r64::new(2,1),r64::new(5,2));
| |
|
| |
| assert!(Point::are_collinear(point1, point1, point1));
| |
| assert!(Point::are_collinear(point1, point1, point2));
| |
| assert!(Point::are_collinear(point1, point2, point1));
| |
| assert!(Point::are_collinear(point2, point1, point1));
| |
| assert!(Point::are_collinear(point1, point2, point4));
| |
| assert!(!Point::are_collinear(point1, point2, point3));
| |
| }
| |
|
| |
| #[test]
| |
| fn line_parallel_implies_equal_slope() {
| |
| let line1 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(1,1),r64::new(1,1))).unwrap();
| |
| let line2 = Line::from_points(Point::new(r64::new(0,1),r64::new(1,1)), Point::new(r64::new(1,1),r64::new(2,1))).unwrap();
| |
| assert_eq!(line1.slope(), line2.slope());
| |
| }
| |
|
| |
|
| // Tests for the module `Line`.
| | // Is the scale word s = s(x, y, z) abstractly mv3? |
| #[test]
| | pub fn maximum_variety(s: &str) -> usize { |
| fn line_equality() {
| | let mut result = 0; |
| let line1 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(1,1),r64::new(1,1))).unwrap();
| | let floor_half: usize = s.len()/2; |
| let line2 = Line::from_slope_and_point(Slope::integer(1), Point::new(r64::new(2,1),r64::new(2,1))); // equal to `line1`
| | for l in 1..(floor_half+1) { |
| assert_eq!(line1, line2);
| | let mut sizes: BTreeSet<(usize, usize, usize)> = BTreeSet::new(); |
| let line3 = Line::from_slope_and_point(Slope::integer(1), Point::new(r64::new(1,2),r64::new(1,2)));
| | for b in 0..(s.len()) { |
| assert_eq!(line3, line2);
| | let mut size: (usize, usize, usize) = (0, 0, 0); |
|
| | let sl = &slice_cyclic_string(s, b, l); |
| let line4 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(2,1),r64::new(1,1))).unwrap();
| | let chars_in_sl = sl.chars(); |
| assert_ne!(line1, line4);
| | for ch in chars_in_sl { |
| let line5 = Line::from_points(Point::new(r64::new(1,1),r64::new(2,1)), Point::new(r64::new(0,1),r64::new(1,1))).unwrap();
| | match ch { |
| assert_ne!(line1, line5);
| | 'x' => {size.0 += 1;}, |
| }
| | 'y' => {size.1 += 1;}, |
|
| | 'z' => {size.2 += 1;}, |
| #[test]
| | _ => {panic!()}, |
| fn line_intersection_parallel_equal() {
| | } |
| let line1 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(1,1),r64::new(1,1))).unwrap();
| | } |
| let line2 = Line::from_slope_and_point(Slope::integer(1), Point::new(r64::new(2,1),r64::new(2,1))); // equal to `line1` | | sizes.insert(size); |
| assert_eq!(intersection(line1, line2), None);
| | } |
| }
| | result = max(result, sizes.len()); // update result |
| | } |
| | result |
| | } |
|
| |
|
|
| | // Is the scale word s = s(x, y, z) pairwise well-formed (pwf)? |
| #[test]
| | pub fn is_pwf(s: &str) -> bool { |
| fn line_intersection_parallel_unequal() {
| | let s1 = s.replace("z", "y"); |
| let line1 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(1,1),r64::new(1,1))).unwrap();
| | let s2 = s.replace("z", "x"); |
| let line2 = Line::from_points(Point::new(r64::new(0,1),r64::new(1,2)), Point::new(r64::new(1,1),r64::new(3,2))).unwrap(); // not equal but parallel to `line1`
| | let s3 = s.replace("y", "x").replace("z", "y"); |
| assert_eq!(intersection(line1, line2), None);
| | return is_sv2(&s1)&&is_sv2(&s2)&&is_sv2(&s3); |
| | } |
|
| |
|
| let line3 = Line::from_points(Point::new(r64::new(0,1),r64::new(1,1)), Point::new(r64::new(1,1),r64::new(2,1))).unwrap(); // not equal but parallel to `line1`
| | // Is the scale word s = s(a,b,c) pairwise well-formed (pwf)? |
| assert_eq!(intersection(line1, line3), None);
| | pub fn is_pmos(s: &str) -> bool { |
| }
| | let s1 = s.replace("z", "y"); |
|
| | let s2 = s.replace("z", "x"); |
| #[test]
| | let s3 = s.replace("y", "x").replace("z", "y"); |
| fn line_intersection_nonparallel() {
| | maximum_variety(&s1) == 2 && maximum_variety(&s2) == 2 && maximum_variety(&s3) == 2 |
| let line1 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(1,1),r64::new(1,1))).unwrap();
| | } |
| let line2 = Line::from_points(Point::new(r64::new(0,1),r64::new(0,1)), Point::new(r64::new(1,1),r64::new(2,1))).unwrap();
| |
| assert_eq!(intersection(line1, line2), Some(Point::new(r64::new(0,1), r64::new(0,1))));
| |
| }
| |
|
| |
| #[test]
| |
| fn point_order_coordinate_vecs() {
| |
| let e1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let e2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
|
| |
|
| let mut points = vec![e1, -e1, e2, -e2, e1+e2, e1-e2, -e1-e2, -e1+e2];
| | // Is the scale word s = s(x, y, z), x > y > z > 0, monotone-mos? |
| points = Point::ordered_ccw(points);
| | pub fn is_mmos(scale: &str) -> bool { |
| assert!(points[0] == -e1-e2 && points[1] == -e2 && points[2] == e1-e2 && points[3] == e1 && points[4] == e1+e2 && points[5] == e2 && points[6] == -e1+e2 && points[7] == -e1 );
| | let s1 = scale.replace("y", "x").replace("z", "y"); // m = L |
| }
| | let s2 = scale.replace("z", "y"); // m = s |
|
| | let s3 = scale.replace("z", ""); // s = 0 |
| #[test]
| | maximum_variety(&s1) == 2 && maximum_variety(&s2) == 2 && maximum_variety(&s3) == 2 |
| fn convexpolygon_disallow_polygon_with_fewer_than_3_vertices() {
| | } |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let result0 = ConvexPolygon::new(vec![]);
| |
| assert_eq!(result0, None);
| |
| let result1 = ConvexPolygon::new(vec![p1]);
| |
| assert_eq!(result1, None);
| |
| let result2 = ConvexPolygon::new(vec![p1, p2]);
| |
| assert_eq!(result2, None);
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_vertices_of_triangle_are_sorted() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
|
| |
|
| let triangle = ConvexPolygon::new(vec![p1, p2, p3]);
| | // Is the scale s = s(x, y) a single period mos (aka well-formed, wf)? |
| assert!(triangle.is_some());
| | pub fn is_sv2(scale: &str) -> bool { |
| | if scale.len() <= 1{ |
| | return false; |
| | } |
| | let floor_half : usize = scale.len()/2; |
| | for slice_length in 1..(floor_half+1) { |
| | let mut sizes: BTreeSet<(usize, usize)> = BTreeSet::new(); |
| | for basepoint in 0..(scale.len()) { |
| | let mut size: (usize, usize) = (0, 0); |
| | let sl = &slice_cyclic_string(scale, basepoint, slice_length); |
| | let chars_in_sl = sl.chars(); |
| | for c in chars_in_sl { |
| | match c { |
| | 'x' => {size.0 += 1;}, |
| | 'y' => {size.1 += 1;}, |
| | _ => {panic!()}, |
| | } |
| | } |
| | sizes.insert(size); |
| | } |
| | if sizes.len() != 2 { |
| | return false; |
| | } |
| | } |
| | return true; |
| | } |
|
| |
|
| let vertices = triangle.unwrap().vertices; | | // Return the balance of `s` = `s`(x, y, z), defined as max { | |w|_{x_i} - |w'|_{x_i} | : x_i is a letter of `s` and k = len(w) = len(w') }. |
| assert_eq!(vertices, vec![p3, p1, p2]);
| | pub fn balance(scale: &str) -> usize { |
| | if scale.len() <= 1 { |
| | scale.len() |
| | } else { |
| | let mut result = 0; |
| | let floor_half: usize = scale.len()/2; |
| | for slice_length in 1..(floor_half+1) { |
| | let (mut min_x, mut max_x, mut min_y, mut max_y, mut min_z, mut max_z) = (usize::MAX, 0, usize::MAX, 0, usize::MAX, 0); // min and max for each length |
| | let mut sizes: BTreeSet<(usize, usize, usize)> = BTreeSet::new(); |
| | for basepoint in 0..(scale.len()) { |
| | let mut size: (usize, usize, usize) = (0, 0, 0); |
| | let sl = &slice_cyclic_string(scale, basepoint, slice_length); |
| | let chars_in_sl = sl.chars(); |
| | for ch in chars_in_sl { |
| | match ch { |
| | 'x' => {size.0 += 1;}, |
| | 'y' => {size.1 += 1;}, |
| | 'z' => {size.2 += 1;}, |
| | _ => {panic!()}, |
| | } |
| | } |
| | sizes.insert(size); |
| | } |
| | for vector in &sizes { |
| | min_x = min(min_x, vector.0); |
| | max_x = max(max_x, vector.0); |
| | min_y = min(min_y, vector.1); |
| | max_y = max(max_y, vector.1); |
| | min_z = min(min_z, vector.2); |
| | max_z = max(max_z, vector.2); |
| | } |
| | result = max(result, max(max(max_x.abs_diff(min_x), max_y.abs_diff(min_y)), max_z.abs_diff(min_z))); |
| | } |
| | result |
| } | | } |
|
| | } |
| #[test]
| |
| fn convexpolygon_vertices_of_square_are_sorted() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
|
| |
| let square = ConvexPolygon::new(vec![p1, p2, p3, p4]);
| |
| assert!(square.is_some());
| |
|
| |
|
| let vertices = square.unwrap().vertices; | | // Function for computing two properties, maximum variety and balance, while executing the subroutine to extract words only once. |
| assert_eq!(vertices, vec![p3, p1, p4, p2]);
| | pub fn max_variety_and_balance(s: &str) -> (usize, usize) { |
| | if s.len() <= 1 { |
| | (s.len(), 0) |
| | } else { |
| | let mut mv_result = 0; |
| | let mut balance_result = 0; |
| | let (mut min_x, mut max_x, mut min_y, mut max_y, mut min_z, mut max_z) = (usize::MAX, 0, usize::MAX, 0, usize::MAX, 0); |
| | let floor_half: usize = s.len()/2; |
| | for l in 1..(floor_half+1) { |
| | let mut sizes: BTreeSet<(usize, usize, usize)> = BTreeSet::new(); |
| | for b in 0..s.len() { |
| | let mut size: (usize, usize, usize) = (0, 0, 0); |
| | let sl = &slice_cyclic_string(s, b, l); |
| | let chars_in_sl = sl.chars(); |
| | for ch in chars_in_sl { |
| | match ch { |
| | 'x' => {size.0 += 1;}, |
| | 'y' => {size.1 += 1;}, |
| | 'z' => {size.2 += 1;}, |
| | _ => {panic!()}, |
| | } |
| | } |
| | sizes.insert(size); |
| | } |
| | for vector in &sizes { |
| | min_x = min(min_x, vector.0); |
| | max_x = max(max_x, vector.0); |
| | min_y = min(min_y, vector.1); |
| | max_y = max(max_y, vector.1); |
| | min_z = min(min_z, vector.2); |
| | max_z = max(max_z, vector.2); |
| | } |
| | mv_result = max(mv_result, sizes.len()); // update result |
| | balance_result = max(balance_result, max(max(max_x.abs_diff(min_x), max_y.abs_diff(min_y)), max_z.abs_diff(min_z))); |
| | } |
| | (mv_result, balance_result) |
| } | | } |
| | } |
|
| |
|
| #[test]
| |
| fn convexpolygon_vertices_of_5gon_are_sorted() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let p5 = Point::new(r64::new(1,2),r64::new(2,1));
| |
|
| |
| let pentagon = ConvexPolygon::new(vec![p1, p2, p3, p4, p5]);
| |
| assert!(pentagon.is_some());
| |
|
| |
|
| let vertices = pentagon.unwrap().vertices; | | // Return a Christoffel word with `a` x's and `b` y's. |
| assert_eq!(vertices, vec![p3, p1, p4, p5, p2]); | | // Algorithm from Bulgakova et al, 2023, "On balanced and abelian properties of circular words over a ternary alphabet". |
| | pub fn christoffel_word(a: u32, b: u32) -> String { |
| | let d = gcd_u32(a, b); |
| | if d == 1 { |
| | let mut result : String = String::from(""); |
| | let (mut current_x, mut current_y) = (0u32, 0u32); // Start from the (0,0) vector. |
| | while current_x < a || current_y < b { |
| | if (current_y+1)*a <= b * (current_x) { // if making the (0,1) step doesn't lead to going above the line |
| | current_y += 1; // append the b step and reflect that in the plane vector. |
| | result.push('y'); |
| | } else { |
| | current_x += 1; |
| | result.push('x'); |
| | } |
| | } |
| | result |
| | } else { |
| | std::iter::repeat(christoffel_word(a/d, b/d)).take(d.try_into().unwrap()).collect::<String>() |
| } | | } |
|
| | } |
| #[test]
| |
| fn convexpolygon_triangle_has_centroid_inside() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
|
| |
|
| let triangle = ConvexPolygon::new(vec![p1, p2, p3]);
| | /// Billiard scales with the given signature up to rotation. |
| let centroid = ConvexPolygon::centroid(&triangle.clone().unwrap()); | | pub fn billiard_scales(a: u32, b: u32, c: u32) -> Vec<String> { |
| assert!(triangle.clone().unwrap().has_point_inside(centroid));
| | let mut result = Vec::<String>::new(); |
| | let regions: Vec<ConvexPolygon> = plane_geometry::project_and_partition(a, b, c); |
| | let n = a + b + c; |
| | for region in regions { |
| | let mut word = String::new(); |
| | let mut next_point = region.centroid(); |
| | let first_point = region.centroid().clone(); |
| | for i in 0..n { |
| | let res = plane_geometry::advance(a, b, c, next_point).unwrap(); |
| | next_point = res.0; |
| | debug_assert!(i == n-1 || first_point != next_point); // There should not be a subperiod in the orbit of the "billiard ball". |
| | word.push_str(&(res.1)); |
| | } |
| | result.push(word); |
| | debug_assert_eq!(first_point, next_point); // The last point should come back to the first point. When i = n - 1, next_point gets set to point number n. |
| } | | } |
| | | result = equivalence_class_representatives::<String>(result, are_rotationally_equivalent); |
| #[test] | | result |
| fn convexpolygon_vertex_is_not_inside() {
| | } |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| | </syntaxhighlight> |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| | === main.rs === |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| | <syntaxhighlight lang="rust> |
| | use scale_theory; |
| | use std::env; |
| | use std::fs::File; |
| | use std::io::Write; |
| | use std::time::Instant; |
|
| |
|
| let triangle = ConvexPolygon::new(vec![p1, p2, p3]);
| | /// Gives the greatest common denominator of the two inputs, unless that's 2^63. |
| assert!(!triangle.unwrap().has_point_inside(p1)); | | /// 2^63 doesn't fit in an `i64`, so it returns -2^63, which does. |
| | pub fn gcd(u: i64, v: i64) -> i64 { |
| | // `wrapping_abs` gives a number's absolute value, unless that's 2^63. 2^63 |
| | // won't fit in `i64`, so it gives -2^63 instead. |
| | let mut v = v.wrapping_abs() as u64; |
| | if u == 0 { |
| | return v as i64; |
| | } |
| | let mut u = u.wrapping_abs() as u64; |
| | if v == 0 { |
| | return u as i64; |
| } | | } |
|
| |
| #[test]
| |
| fn convexpolygon_point_on_edge_is_not_inside() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
|
| |
|
| let triangle = ConvexPolygon::new(vec![p1, p2, p3]);
| | // `|` is bitwise OR. `trailing_zeros` quickly counts a binary number's |
| let p_edge = Point::midpoint(&p1, &p2);
| | // trailing zeros, giving its prime factorization's exponent on two. |
| assert!(!triangle.unwrap().has_point_inside(p_edge));
| | let gcd_exponent_on_two = (u | v).trailing_zeros(); |
| }
| |
| | |
| #[test]
| |
| fn convexpolygon_point_outside_is_not_inside() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
|
| |
|
| let triangle = ConvexPolygon::new(vec![p1, p2, p3]);
| | // `>>=` divides the left by two to the power of the right, storing that in |
| let p_outside = Point::new(r64::new(1,1),r64::new(1,1));
| | // the left variable. `u` divided by its prime factorization's power of two |
| assert!(!triangle.unwrap().has_point_inside(p_outside));
| | // turns it odd. |
| }
| | u >>= u.trailing_zeros(); |
| | v >>= v.trailing_zeros(); |
|
| |
|
| #[test] | | while u != v { |
| fn line_point_line_config_test_vertical_line() {
| | if u < v { |
| let p_left = Point::new(r64::new(-1,1),r64::new(0,1));
| | // Swap the variables' values with each other. |
| let p_on_line = Point::new(r64::new(0,1),r64::new(1,1));
| | core::mem::swap(&mut u, &mut v); |
| let p_right = Point::new(r64::new(1,1),r64::new(1,1)); | | } |
| let vertical_line_up = Line::from_slope_and_point(Slope::infinity(), Point::zero()); | | u -= v; |
| assert_eq!(vertical_line_up.point_line_config(p_left), PointLineConfiguration::Left); | | u >>= u.trailing_zeros(); |
| assert_eq!(vertical_line_up.point_line_config(p_on_line), PointLineConfiguration::OnTheLine);
| |
| assert_eq!(vertical_line_up.point_line_config(p_right), PointLineConfiguration::Right);
| |
| let vertical_line_down = Line::from_points(Point::zero(), Point::new(r64::new(0,1), r64::new(-1,1)));
| |
| assert_eq!(vertical_line_down.clone().unwrap().point_line_config(p_left), PointLineConfiguration::Right);
| |
| assert_eq!(vertical_line_down.clone().unwrap().point_line_config(p_on_line), PointLineConfiguration::OnTheLine);
| |
| assert_eq!(vertical_line_down.clone().unwrap().point_line_config(p_right), PointLineConfiguration::Left);
| |
| } | | } |
|
| |
|
| #[test] | | // `<<` multiplies the left by two to the power of the right. |
| fn line_point_line_config_test_pos_slope() {
| | (u << gcd_exponent_on_two) as i64 |
| let p_above = Point::new(r64::new(0,1),r64::new(1,1));
| | } |
| let p_on_line = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let p_below = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let line_up_right = Line::from_slope_and_point(Slope::integer(1), Point::zero());
| |
| assert_eq!(line_up_right.point_line_config(p_above), PointLineConfiguration::Left);
| |
| assert_eq!(line_up_right.point_line_config(p_on_line), PointLineConfiguration::OnTheLine);
| |
| assert_eq!(line_up_right.point_line_config(p_below), PointLineConfiguration::Right);
| |
| let line_down_left = Line::from_points(Point::zero(), Point::new(r64::new(-1,1), r64::new(-1,1)));
| |
| assert_eq!(line_down_left.clone().unwrap().point_line_config(p_above), PointLineConfiguration::Right);
| |
| assert_eq!(line_down_left.clone().unwrap().point_line_config(p_on_line), PointLineConfiguration::OnTheLine);
| |
| assert_eq!(line_down_left.clone().unwrap().point_line_config(p_below), PointLineConfiguration::Left);
| |
| }
| |
|
| |
| #[test] | |
| fn line_point_line_config_test_neg_slope() {
| |
| let p_below = Point::new(r64::new(-1,1),r64::new(0,1));
| |
| let p_on_line = Point::new(r64::new(-1,1),r64::new(1,1));
| |
| let p_above = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let line_down_right = Line::from_slope_and_point(Slope::integer(-1), Point::zero());
| |
| assert_eq!(line_down_right.point_line_config(p_below), PointLineConfiguration::Right);
| |
| assert_eq!(line_down_right.point_line_config(p_on_line), PointLineConfiguration::OnTheLine);
| |
| assert_eq!(line_down_right.point_line_config(p_above), PointLineConfiguration::Left);
| |
| let line_up_left = Line::from_points(Point::zero(), Point::new(r64::new(-1,1), r64::new(1,1)));
| |
| assert_eq!(line_up_left.clone().unwrap().point_line_config(p_below), PointLineConfiguration::Left);
| |
| assert_eq!(line_up_left.clone().unwrap().point_line_config(p_on_line), PointLineConfiguration::OnTheLine);
| |
| assert_eq!(line_up_left.clone().unwrap().point_line_config(p_above), PointLineConfiguration::Right);
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_if_all_vertices_are_to_one_side() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let square = ConvexPolygon::raw(vec![p3, p1, p4, p2]);
| |
|
| |
| let line_to_the_right = Line::from_slope_and_point(Slope::infinity(), Point::new(r64::new(2,1), r64::new(0,1)));
| |
| let line_to_the_left = Line::from_slope_and_point(Slope::infinity(), Point::new(r64::new(-1,1), r64::new(0,1)));
| |
| let line_above = Line::from_slope_and_point(Slope::zero(), Point::new(r64::new(0,1), r64::new(2,1)));
| |
| let line_below = Line::from_slope_and_point(Slope::zero(), Point::new(r64::new(0,1), r64::new(-1,1)));
| |
|
| |
| assert_eq!(square.subdivide(line_to_the_right), (Some(ConvexPolygon::raw(vec![p3, p1, p4, p2])), Option::<ConvexPolygon>::None) );
| |
| assert_eq!(square.subdivide(line_to_the_left), (Option::<ConvexPolygon>::None, Some(ConvexPolygon::raw(vec![p3, p1, p4, p2]))) );
| |
| assert_eq!(square.subdivide(line_above), (Option::<ConvexPolygon>::None, Some(ConvexPolygon::raw(vec![p3, p1, p4, p2]))) );
| |
| assert_eq!(square.subdivide(line_below), (Some(ConvexPolygon::raw(vec![p3, p1, p4, p2])), Option::<ConvexPolygon>::None) );
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_new_line_cuts_two_edges() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let square = ConvexPolygon::raw(vec![p3, p1, p4, p2]);
| |
|
| |
| let line = Line::from_slope_and_point(Slope::infinity(), Point::new(r64::new(1,2), r64::new(0,1)));
| |
| let (left_polygon, right_polygon) = square.subdivide(line);
| |
| assert_eq!(left_polygon, Some(ConvexPolygon::raw(vec![Point::new(r64::new(0,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,2),r64::new(0,1)),
| |
| Point::new(r64::new(1,2),r64::new(1,1)),
| |
| Point::new(r64::new(0,1),r64::new(1,1))
| |
| ])), "`left_polygon` is not correct");
| |
| assert_eq!(right_polygon, Some(ConvexPolygon::raw(vec![Point::new(r64::new(1,2),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(1,1)),
| |
| Point::new(r64::new(1,2),r64::new(1,1))
| |
| ])), "`right_polygon` is not correct");
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_new_line_grazes_one_vertex() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let square = ConvexPolygon::raw(vec![p3, p1, p4, p2]);
| |
|
| |
| let line = Line::from_slope_and_point(Slope::integer(1), Point::new(r64::new(0,1), r64::new(1,1)));
| |
| let (left_polygon, right_polygon) = square.subdivide(line);
| |
|
| |
| assert_eq!(left_polygon, None, "`left_polygon` is not correct");
| |
| assert_eq!(right_polygon, Some(ConvexPolygon::raw(vec![p3, p1, p4, p2])), "`right_polygon` is not correct");
| |
|
| |
|
| let line2 = Line::from_slope_and_point(Slope::integer(1), Point::new(r64::new(1,1), r64::new(0,1)));
| | fn main() -> std::io::Result<()> { |
| let (left_polygon, right_polygon) = square.subdivide(line2);
| | env::set_var("RUST_BACKTRACE", "1"); |
|
| | let mut f = File::create("billiard_scales_mediawiki")?; |
| assert_eq!(left_polygon, Some(ConvexPolygon::raw(vec![p3, p1, p4, p2])), "`left_polygon` is not correct");
| | let start_time = Instant::now(); |
| assert_eq!(right_polygon, None, "`right_polygon` is not correct");
| | for a in 1..=29 { |
| | for b in 1..=a { |
| | for c in 1..=b { |
| | if a + b + c <= 31 && gcd(gcd(a as i64, b as i64), c as i64) == 1 { |
| | writeln!(&f, "== {}x{}y{}z ==", a, b, c); |
| | let billiard_scales = scale_theory::billiard_scales(a, b, c); |
| | for scale in &billiard_scales{ |
| | writeln!(&f, "# {}\n#* maximum variety {}, balance {}", scale, scale_theory::maximum_variety(scale), scale_theory::balance(scale)); |
| | } |
| | } |
| | } |
| | } |
| } | | } |
| | | let elapsed_time = start_time.elapsed(); |
| #[test]
| | println!("Done in {:?}", elapsed_time); |
| fn convexpolygon_new_line_coincides_with_one_edge() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let square = ConvexPolygon::raw(vec![p3, p1, p4, p2]);
| |
|
| |
| let line = Line::from_slope_and_point(Slope::infinity(), Point::new(r64::new(0,1), r64::new(0,1)));
| |
| let (left_polygon, right_polygon) = square.subdivide(line);
| |
|
| |
| assert_eq!(left_polygon, None, "`left_polygon` is not correct");
| |
| assert_eq!(right_polygon, Some(ConvexPolygon::raw(vec![p3, p1, p4, p2])), "`right_polygon` is not correct");
| |
|
| |
|
| let line = Line::from_slope_and_point(Slope::integer(0), Point::new(r64::new(1,1), r64::new(0,1)));
| | Ok(()) |
| let (left_polygon, right_polygon) = square.subdivide(line);
| |
|
| |
| assert_eq!(left_polygon, Some(ConvexPolygon::raw(vec![p3, p1, p4, p2])), "`left_polygon` is not correct");
| |
| assert_eq!(right_polygon, None, "`right_polygon` is not correct");
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_new_line_cuts_one_edge() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let square = ConvexPolygon::raw(vec![p3, p1, p4, p2]);
| |
|
| |
| let line = Line::from_slope_and_point(Slope::integer(2), Point::new(r64::new(0,1), r64::new(0,1)));
| |
| let (left_polygon, right_polygon) = square.subdivide(line);
| |
|
| |
| assert_eq!(left_polygon, Some(ConvexPolygon::raw(vec![Point::new(r64::new(0,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,2),r64::new(1,1)),
| |
| Point::new(r64::new(0,1),r64::new(1,1)),
| |
| ])), "`left_polygon is not correct");
| |
| assert_eq!(right_polygon, Some(ConvexPolygon::raw(vec![Point::new(r64::new(0,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(1,1)),
| |
| Point::new(r64::new(1,2),r64::new(1,1))
| |
| ])), "`right_polygon` is not correct");
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_new_line_cuts_no_edges() {
| |
| let p1 = Point::new(r64::new(1,1),r64::new(0,1));
| |
| let p2 = Point::new(r64::new(0,1),r64::new(1,1));
| |
| let p3 = Point::new(r64::new(0,1),r64::new(0,1));
| |
| let p4 = Point::new(r64::new(1,1),r64::new(1,1));
| |
| let square = ConvexPolygon::raw(vec![p3, p1, p4, p2]);
| |
|
| |
| let line = Line::from_slope_and_point(Slope::integer(1), Point::new(r64::new(0,1), r64::new(0,1)));
| |
| let (left_polygon, right_polygon) = square.subdivide(line);
| |
|
| |
| assert_eq!(left_polygon, Some(ConvexPolygon::raw(vec![Point::new(r64::new(0,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(1,1)),
| |
| Point::new(r64::new(0,1),r64::new(1,1))
| |
| ])), "`left_polygon` is not correct");
| |
| assert_eq!(right_polygon, Some(ConvexPolygon::raw(vec![Point::new(r64::new(0,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(0,1)),
| |
| Point::new(r64::new(1,1),r64::new(1,1))
| |
| ])), "`right_polygon` is not correct");
| |
| }
| |
|
| |
| #[test]
| |
| fn convexpolygon_subdivide_test() { | |
| let hexagon = projected_cube(5, 2, 3);
| |
| assert_eq!(hexagon, ConvexPolygon::new(vec![Point::new(r64::new(1,1),r64::new(0,1)), Point::new(r64::new(0,1),r64::new(1,1)), Point::new(r64::new(1,1),r64::new(1,1)),
| |
| Point::new(r64::new(-5,3),r64::new(-2,3)), Point::new(r64::new(-2,3),r64::new(-2,3)), Point::new(r64::new(-5,3),r64::new(1,3))]) );
| |
| let line = Line::from_slope_and_point(Slope::integer(0), Point::new(r64::new(0,1), r64::new(-1,3)));
| |
| let divided = hexagon.unwrap().subdivide(line);
| |
| assert_eq!(divided.0, ConvexPolygon::new(vec![Point::new(r64::new(1,1),r64::new(0,1)), Point::new(r64::new(0,1),r64::new(1,1)), Point::new(r64::new(1,1),r64::new(1,1)),
| |
| Point::new(r64::new(-5,3),r64::new(-1,3)), Point::new(r64::new(1,6),r64::new(-1,3)), Point::new(r64::new(-5,3),r64::new(1,3))]) );
| |
| assert_eq!(divided.1, ConvexPolygon::new(vec![ Point::new(r64::new(-5,3),r64::new(-1,3)), Point::new(r64::new(1,6),r64::new(-1,3)),
| |
| Point::new(r64::new(-5,3),r64::new(-2,3)), Point::new(r64::new(-2,3),r64::new(-2,3)) ]) );
| |
| let line2 = Line::from_slope_and_point(Slope::integer(0), Point::new(r64::new(0,1), r64::new(0,1)));
| |
| let divided_again = (divided.0).unwrap().subdivide(line2);
| |
|
| |
| assert_eq!(divided_again.0, ConvexPolygon::new(vec![Point::new(r64::new(1,1),r64::new(0,1)), Point::new(r64::new(0,1),r64::new(1,1)), Point::new(r64::new(1,1),r64::new(1,1)),
| |
| Point::new(r64::new(-5,3),r64::new(0,1)), Point::new(r64::new(-5,3),r64::new(1,3))]) );
| |
| assert_eq!(divided_again.1, ConvexPolygon::new(vec![Point::new(r64::new(1,1),r64::new(0,1)), Point::new(r64::new(1,6),r64::new(-1,3)),
| |
| Point::new(r64::new(-5,3),r64::new(0,1)), Point::new(r64::new(-5,3),r64::new(-1,3))]) );
| |
| }
| |
| } | | } |
| | </syntaxhighlight> |
| | === Cargo.toml === |
| | <syntaxhighlight lang="toml"> |
| | [package] |
| | name = "billiard_scales" |
| | version = "0.0.0" |
| | edition = "2021" |
| | |
| | # See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html |
| | |
| | [dependencies] |
| | num-rational = "0.4" |
| | num-traits = "0.2.17" |
| </syntaxhighlight> | | </syntaxhighlight> |