classes.rs
6.6 kB · rust · 204 lines
1use crate::fills::cell_index;2use crate::tables::write_csv;3use mrlymath::bang::baseq::fill_from_corners;4use mrlymath::bang::counting;5use mrlymath::bang::factory::code_to_corners;6use mrlymath::bang::Code;7use mrlymath::rules::render;8use num_bigint::BigUint;9use std::collections::{BTreeMap, BTreeSet};10use std::path::Path;11use std::thread;1213const SIDES: [usize; 3] = [3, 5, 7];14const LEVELS: [usize; 3] = [1, 2, 3];15const SEQUENCE: usize = 9;16const TABLE: usize = 8;17const CLASS_SUM_LIMIT: usize = 6;18const PUBLISHED: [&str; 16] = [19 "2",20 "4",21 "12",22 "64",23 "700",24 "17424",25 "1053696",26 "160579584",27 "62856336636",28 "63812936890000",29 "168895157342195152",30 "1169048914836855865344",31 "21209591746609937928524800",32 "1010490883477487017627972550656",33 "126641164340871500483202065902080000",34 "41817338589698457759723104703370865147904",35];3637fn brute(design: Code, number: usize, dimension: usize, level: usize) -> u128 {38 let tile = render(39 |r| design >> cell_index(r, 2) & 1 == 1,40 number,41 dimension,42 2,43 )44 .expect("the tile renders");45 u128::from(tile.fractal(level).sum())46}4748fn closed(design: Code, number: usize, dimension: usize, level: usize) -> u128 {49 let corners = code_to_corners(design, dimension, 2).expect("the design fits its corners");50 fill_from_corners(&corners, number, dimension).pow(level as u32)51}5253fn profile(design: Code, dimension: usize) -> Vec<usize> {54 let mut out = vec![0usize; dimension + 1];55 for corner in 0..1usize << dimension {56 if design >> corner & 1 == 1 {57 out[corner.count_ones() as usize] += 1;58 }59 }60 out61}6263fn binomial(n: usize, k: usize) -> u128 {64 (0..k).fold(1u128, |acc, i| acc * (n - i) as u128 / (i as u128 + 1))65}6667pub fn a129824(dimension: usize) -> BigUint {68 (0..=dimension)69 .map(|k| BigUint::from(1 + binomial(dimension, k)))70 .product()71}7273fn weight_at_most_one(dimension: usize) -> Code {74 (0..1usize << dimension)75 .filter(|corner| corner.count_ones() <= 1)76 .map(|corner| 1u128 << corner)77 .sum()78}7980fn compare(dimension: usize) -> (usize, usize) {81 let designs = 1u128 << (1usize << dimension);82 let workers = thread::available_parallelism().map_or(1, |n| n.get());83 let chunk = (designs as usize).div_ceil(workers) as u128;84 let totals: Vec<(usize, usize)> = thread::scope(|scope| {85 let handles: Vec<_> = (0..workers as u128)86 .map(|w| {87 scope.spawn(move || {88 let mut checks = 0usize;89 let mut bad = 0usize;90 for design in w * chunk..((w + 1) * chunk).min(designs) {91 for &n in &SIDES {92 for &level in &LEVELS {93 checks += 1;94 if brute(design, n, dimension, level)95 != closed(design, n, dimension, level)96 {97 bad += 1;98 }99 }100 }101 }102 (checks, bad)103 })104 })105 .collect();106 handles107 .into_iter()108 .map(|handle| handle.join().expect("the worker finishes"))109 .collect()110 });111 totals.iter().fold((0, 0), |(c, b), &(x, y)| (c + x, b + y))112}113114fn distinct(dimension: usize) -> (usize, usize) {115 let mut seen: BTreeMap<Vec<u128>, BTreeSet<Vec<usize>>> = BTreeMap::new();116 for design in 0..1u128 << (1usize << dimension) {117 let sequence: Vec<u128> = (1..=SEQUENCE)118 .map(|n| closed(design, n, dimension, 1))119 .collect();120 seen.entry(sequence)121 .or_default()122 .insert(profile(design, dimension));123 }124 let collisions = seen.values().filter(|profiles| profiles.len() > 1).count();125 (seen.len(), collisions)126}127128pub fn report(path: &Path) {129 println!("fill classes at base 2: two generators, the closed form and the published terms");130 let sponge = weight_at_most_one(3);131 let carpet = weight_at_most_one(2);132 println!(133 "Menger sponge n 3: level 1 {} level 2 {}; carpet n 3: {} n 5: {}",134 brute(sponge, 3, 3, 1),135 brute(sponge, 3, 3, 2),136 brute(carpet, 3, 2, 1),137 brute(carpet, 5, 2, 1)138 );139 for dimension in [2usize, 3] {140 let (checks, bad) = compare(dimension);141 println!("D {dimension}: {checks} generator checks, {bad} mismatches");142 }143 for dimension in 1..=4 {144 let (count, collisions) = distinct(dimension);145 let closed = a129824(dimension);146 println!(147 "D {dimension}: distinct fill sequences {count}, closed form {closed}, match {}, profile collisions {collisions}",148 count.to_string() == closed.to_string()149 );150 }151 let agreed = PUBLISHED152 .iter()153 .enumerate()154 .filter(|(d, term)| a129824(*d).to_string() == **term)155 .count();156 let terms: Vec<String> = (0..PUBLISHED.len())157 .map(|d| a129824(d).to_string())158 .collect();159 println!(160 "A129824 closed form D 0..{}: {}",161 PUBLISHED.len() - 1,162 terms.join(", ")163 );164 println!("published terms matched: {agreed} of {}", PUBLISHED.len());165 println!(166 "{:<3}{:>12}{:>18}{:>16}{:>8}",167 "D", "designs", "A000616", "A129824", "ratios"168 );169 let header: Vec<String> = [170 "D",171 "total_designs",172 "shape_classes_A000616",173 "fractal_dim_classes_A129824",174 "limiting_ratio_classes",175 ]176 .iter()177 .map(|name| (*name).to_string())178 .collect();179 let mut records = Vec::new();180 for dimension in 1..=TABLE {181 let designs = BigUint::from(2u32).pow(1 << dimension);182 let shapes = (dimension <= CLASS_SUM_LIMIT)183 .then(|| counting::distinct_designs(dimension).expect("the class sums close"))184 .map(|value| value.to_string())185 .unwrap_or_default();186 let classes = a129824(dimension);187 let ratios = (1u32 << dimension) + 1;188 let shown = if designs.to_string().len() > 12 {189 format!("2^{}", 1 << dimension)190 } else {191 designs.to_string()192 };193 println!("{dimension:<3}{shown:>12}{shapes:>18}{classes:>16}{ratios:>8}");194 records.push(vec![195 dimension.to_string(),196 designs.to_string(),197 shapes,198 classes.to_string(),199 ratios.to_string(),200 ]);201 }202 println!("written counts.csv, A000616 by class sums to D = {CLASS_SUM_LIMIT}");203 write_csv(path, &header, &records);204}