Xenharmonic Wiki

User:Inthar/Code: Difference between revisions

← Older editNewer edit →
VisualWikitext
Inthar (talk | contribs)
autopatrolled, emailconfirmed
15,668 edits
Inthar (talk | contribs)
autopatrolled, emailconfirmed
15,668 edits
No edit summary
Line 619: Line 619:
/// Advance a projected point, `pt`, by one letter according to what integer coordinate it next hits as it traverses in the given velocity,
/// Advance a projected point, `pt`, by one letter according to what integer coordinate it next hits as it traverses in the given velocity,
/// and return the next point and the corresponding letter, one of "x", "y", and "z".
/// and return the next point and the corresponding letter, one of "x", "y", and "z".
/// Assumes gcd(a, b, c) == 1.
pub fn advance(a: u32, b: u32, c: u32, pt: Point) -> Result<(Point, String), String> {
pub fn advance(a: u32, b: u32, c: u32, pt: Point) -> Result<(Point, String), String> {
     let projected_cube = projected_cube(a, b, c).unwrap();
     let projected_cube = projected_cube(a, b, c).unwrap();
Line 642: Line 643:
         }
         }
     }
     }
}
</syntaxhighlight>
=== lib.rs ===
<syntaxhighlight lang="rs">
use std::collections::BTreeSet;
use num_rational::Rational64 as r64;
mod plane_geometry;
use plane_geometry::*;
use std::cmp::{min, max};
type EquivalenceRelation<T> = fn(&T, &T) -> bool;
/// Treating `scale` as a cyclic string (that is, "scale[i] == scale[i % scale.len()]"),
/// take a slice of length `slice_length` from `basepoint`; assumes `slice_length` <= `scale`.len()
fn slice_cyclic_string(scale: &str, basepoint: usize, slice_length: usize) -> String {
let arr = scale;
if basepoint + slice_length < scale.len() {
return String::from(&arr[basepoint..(basepoint+slice_length)]);
} else {
let mut ret = String::new();
ret.push_str(&arr[basepoint..]);
ret.push_str(&arr[..(basepoint + slice_length - scale.len())]);
return ret;
}
}
pub fn are_rotationally_equivalent(scale1: &String, scale2: &String)-> bool{
    if scale1.len() != scale2.len() {
        return false;
    } else {
        let length = scale1.len();
        for i in 0..length {
            if slice_cyclic_string(scale1, i, length) == *scale2 {
                return true;
            }
        }
        false
    }
}
/// Assuming that "`equiv` : T x T -> bool" is an equivalence relation,
/// return the equivalence class representatives `set` under `equiv`.
pub fn equivalence_class_representatives<T>(set: Vec<T>, equiv: EquivalenceRelation<T>) -> Vec<T>
where T : Clone + PartialEq {
    let mut result : Vec<T> = Vec::new();
    for i in 0..set.len() {
        let mut found_class_equivalent = false;
        'a: for j in 0..result.len(){ // All vecs in result are nonempty.
            if equiv(&result[j], &set[i]){
                found_class_equivalent = true;
                continue 'a;
            }
        }
        if found_class_equivalent == false {
            result.push(set[i].clone());
        }
    }
    result
}
/// Gives the greatest common denominator of the two inputs, unless that's 2^63.
/// 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;
    }
    // `|` is bitwise OR. `trailing_zeros` quickly counts a binary number's
    // trailing zeros, giving its prime factorization's exponent on two.
    let gcd_exponent_on_two = (u | v).trailing_zeros();
    // `>>=` divides the left by two to the power of the right, storing that in
    // the left variable. `u` divided by its prime factorization's power of two
    // turns it odd.
    u >>= u.trailing_zeros();
    v >>= v.trailing_zeros();
    while u != v {
        if u < v {
            // Swap the variables' values with each other.
            core::mem::swap(&mut u, &mut v);
        }
        u -= v;
        u >>= u.trailing_zeros();
    }
    // `<<` multiplies the left by two to the power of the right.
    (u << gcd_exponent_on_two) as i64
}
fn lcm(m: i64, n: i64) -> i64 {
m*n/gcd(m,n)
}
/// Computes the maximum variety of a ternary scale word.
pub fn maximum_variety(scale: &str) -> usize {
    let mut result = 0;
let floor_half: usize = scale.len()/2;
for slice_length in 1..(floor_half+1) {
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);
}
        result = max(result, sizes.len()); // update result
}
result
}
/// 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') }.
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
    }
}
/// Function for computing two properties, maximum variety and balance, while executing the subroutine to extract words only once at a time for both properties.
pub fn max_variety_and_balance(scale: &str) -> (usize, usize) {
    if scale.len() <= 1 {
        (scale.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 = scale.len()/2;
    for slice_length in 1..(floor_half+1) {
    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);
            }
            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)
    }
}
/// Return a Christoffel word with `a` x's and `b` y's.
/// 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(a as i64, b as i64) as u32;
    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>()
    }
}
/// Billiard scales with the given signature up to rotation.
pub fn billiard_scales(a: u32, b: u32, c: u32) -> Vec<String> {
    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);
    result
}
}
</syntaxhighlight>
</syntaxhighlight>
Line 896: Line 650:
use billiard_scales;
use billiard_scales;
use std::env;
use std::env;
use std::io;
use std::fs::File;
use std::fs::File;
use std::io::Write;
use std::io::Write;
use std::time::Instant;


fn main() -> std::io::Result<()> {
fn main() -> std::io::Result<()> {
     env::set_var("RUST_BACKTRACE", "1");
     env::set_var("RUST_BACKTRACE", "1");
    let mut input = String::new();
     let mut f = File::create("billiard_scales_mediawiki")?;
     let mut f = File::create("billiard_scales_mediawiki")?;
     let start_time = Instant::now();
    println!("This script will generate the list of all ternary billiard scales of the form ax by cz.");
    for n in 3..=31 {
    println!("Enter a > 0");
         writeln!(&f, "== {} notes ==", n);
    std::io::stdin().read_line(&mut input).unwrap();
         for a in 1..=(n - 2) {
     let mut result_a: Result<i64, String> = match input.trim().parse::<i64>(){
             for b in 1..=(n - a - 1) {
        Ok(n) => if n > 0 {Ok(n)} else {Err("input must be a positive integer".to_string())}
                let c = n - a - b;
        Err(e) => Err(e.to_string())
                if a >= b && b >= c && billiard_scales::gcd(billiard_scales::gcd(a as i64, b as i64), c as i64) == 1 {
    };
                    writeln!(&f, "=== {}x{}y{}z ===", a, b, c);
    while let Err(ref e) = result_a {
                    let billiard_scales = billiard_scales::billiard_scales(a, b, c);
         input.clear();
                    for scale in &billiard_scales{
        println!("Input must be a positive integer. Please try again.");
                        writeln!(&f, "# {}\n#* maximum variety {}, balance {}", scale, billiard_scales::maximum_variety(scale), billiard_scales::balance(scale));
        std::io::stdin().read_line(&mut input).unwrap();
                    }
        result_a = match input.trim().parse::<i64>(){
                }
            Ok(n) => if n > 0 {Ok(n)} else {Err("input must be a positive integer".to_string())}
            }
            Err(e) => Err(e.to_string())
        }
        };
    }
    let a: u32 = result_a.expect("Unknown error") as u32;
    input.clear();
    println!("Enter b > 0");
    std::io::stdin().read_line(&mut input).unwrap();
    let mut result_b: Result<i64, String> = match input.trim().parse::<i64>(){
        Ok(n) => if n > 0 {Ok(n)} else {Err("input must be a positive integer".to_string())}
        Err(e) => Err(e.to_string())
    };
    while let Err(ref e) = result_b {
        input.clear();
         println!("Input must be a positive integer. Please try again.");
        std::io::stdin().read_line(&mut input).unwrap();
        result_b = match input.trim().parse::<i64>(){
            Ok(n) => if n > 0 {Ok(n)} else {Err("input must be a positive integer".to_string())}
             Err(e) => Err(e.to_string())
        };
    }
    let b: u32 = result_b.expect("Unknown error") as u32;
    input.clear();
    println!("Enter c > 0");
    std::io::stdin().read_line(&mut input).unwrap();
    let mut result_c: Result<i64, String> = match input.trim().parse::<i64>(){
        Ok(n) => if n > 0 {Ok(n)} else {Err("input must be a positive integer".to_string())}
        Err(e) => Err(e.to_string())
    };
    while let Err(ref e) = result_c {
        input.clear();
        println!("Input must be a positive integer. Please try again.");
        std::io::stdin().read_line(&mut input).unwrap();
        result_c = match input.trim().parse::<i64>(){
            Ok(n) => if n > 0 {Ok(n)} else {Err("input must be a positive integer".to_string())}
            Err(e) => Err(e.to_string())
        };
    }
    let c: u32 = result_c.expect("Unknown error") as u32;
    input.clear();
    let d = billiard_scales::gcd(billiard_scales::gcd(a as i64, b as i64), c as i64);
    if d > 1 {
        println!("Note: a, b, and c are not coprime. You will get scale words with {} periods.", d);
    }
    writeln!(&f, "Billiard scales of signature {}x{}y{}z", a, b, c);
    let billiard_scales = billiard_scales::billiard_scales((a as usize).try_into().unwrap(), (b as usize).try_into().unwrap(), (c as usize).try_into().unwrap());
    for scale in &billiard_scales{
        writeln!(&f, "# {}\n#* maximum variety {}, balance {}", scale, billiard_scales::maximum_variety(scale), billiard_scales::balance(scale));
     }
     }
    let elapsed_time = start_time.elapsed();
     println!("Done. Output has been saved in `billiard_scales_mediawiki` in the current directory.");
     println!("Done in {:?}", elapsed_time);


     Ok(())
     Ok(())
Retrieved from "https://en.xen.wiki/w/User:Inthar/Code"