census.rs

1use crate::count::series;2use crate::design::Design;3use crate::fit::{fit, roots, show};4use std::collections::BTreeMap;56pub fn klein(code: u128, base: usize) -> u128 {7    let d = Design::full(code, base);8    (1..4).map(|g| d.image(g).full_code()).fold(code, u128::min)9}1011pub fn describe(s: &[u128]) -> String {12    let best = (0..s.len().min(8))13        .filter_map(|k| fit(&s[k..]).map(|f| (k, f)))14        .filter(|(_, f)| f.margin >= 2 && f.poly.iter().all(|&c| c.abs() < 1 << 40))15        .min_by_key(|(k, f)| (f.poly.len(), *k));16    match best.map(|(k, f)| (k, Some(f))).unwrap_or((0, None)) {17        (_, None) => "no integer recurrence in range".to_string(),18        (k, Some(f)) => {19            let (found, rest) = roots(&f.poly);20            let tail = if rest.len() > 1 {21                format!(", factor with no integer root {}", show(&rest))22            } else {23                String::new()24            };25            let list: Vec<String> = found.iter().map(|r| r.to_string()).collect();26            format!(27                "from level {k}, order {}, margin {}, roots [{}]{tail}",28                f.poly.len() - 1,29                f.margin,30                list.join(" ")31            )32        }33    }34}3536pub fn census(base: usize, top: usize, deep: usize, codes: &[u128]) {37    let mut classes: BTreeMap<u128, Vec<u128>> = BTreeMap::new();38    for &code in codes {39        classes.entry(klein(code, base)).or_default().push(code);40    }41    let mut groups: BTreeMap<Vec<u128>, Vec<u128>> = BTreeMap::new();42    for &rep in classes.keys() {43        let s = series(base, &Design::full(rep, base).tile, top);44        groups.entry(s).or_default().push(rep);45    }46    let zero = groups47        .iter()48        .filter(|(s, _)| s.iter().all(|&x| x == 0))49        .map(|(_, r)| r.len())50        .sum::<usize>();51    let silent = groups52        .iter()53        .filter(|(s, _)| s.iter().all(|&x| x == 0))54        .flat_map(|(_, r)| r.iter().map(|c| classes[c].len()))55        .sum::<usize>();56    println!("BASE {base}: {} codes, {} Klein classes, {} nonzero loop sequences, {} classes and {} codes never loop, levels 0..{top}", codes.len(), classes.len(), groups.len() - usize::from(zero > 0), zero, silent);57    let mut rows: Vec<(&Vec<u128>, &Vec<u128>)> = groups58        .iter()59        .filter(|(s, _)| s.iter().any(|&x| x > 0))60        .collect();61    rows.sort_by_key(|(s, _)| (*s).clone());62    let mut tally = [0usize; 3];63    let mut known = 0;64    for (s, reps) in rows {65        let shown: Vec<String> = s.iter().take(9).map(|x| x.to_string()).collect();66        let names: Vec<String> = reps67            .iter()68            .map(|&c| {69                let d = Design::full(c, base);70                match d.parity_code() {71                    Some(p) if base > 2 => format!("{c} (parity {p})"),72                    _ => c.to_string(),73                }74            })75            .collect();76        let said = describe(s);77        tally[usize::from(said.contains("no integer root"))78            + 2 * usize::from(said.starts_with("no"))] += 1;79        let size: usize = reps.iter().map(|r| classes[r].len()).sum();80        let found = oeis(s);81        known += usize::from(found.starts_with("OEIS A"));82        println!(83            "  {} | {} | codes {} | {size} codes | {found}",84            shown.join(" "),85            said,86            names.join(" ")87        );88    }89    println!("  fits of margin at least 2: {} with integer roots only, {} with a factor of no integer root, {} with no fit", tally[0], tally[1], tally[2]);90    println!(91        "  in the OEIS by the first seven terms from the first nonzero one: {known} of {}",92        tally.iter().sum::<usize>()93    );94    let open: Vec<u128> = groups95        .iter()96        .filter(|(s, _)| describe(s).starts_with("no"))97        .map(|(_, r)| r[0])98        .collect();99    if deep > top && !open.is_empty() {100        census_deep(base, deep, &open);101    }102}103104fn census_deep(base: usize, top: usize, codes: &[u128]) {105    println!("  the sequences with no fit, rerun to level {top}");106    for &code in codes {107        let d = Design::full(code, base);108        let s = series(base, &d.tile, top);109        let shown: Vec<String> = s.iter().map(|x| x.to_string()).collect();110        println!(111            "    {} | kept {} | {} | {}",112            d.name(),113            d.kept(),114            shown.join(" "),115            describe(&s)116        );117    }118}119120pub fn bases() {121    for (base, top) in [(4usize, 12usize), (5, 10)] {122        let mut codes: Vec<u128> = (1..16u128)123            .map(|p| Design::new(p, base, 2).full_code())124            .collect();125        let all = (1u128 << (base * base)) - 1;126        codes.extend((0..base * base).map(|c| all ^ (1 << c)));127        let mut reps: Vec<u128> = codes.iter().map(|&c| klein(c, base)).collect();128        reps.sort();129        reps.dedup();130        let mut distinct: BTreeMap<Vec<u128>, bool> = BTreeMap::new();131        println!("BASE {base}: the fifteen parity codes and every one-cell deletion, {} Klein classes, levels 0..{top}", reps.len());132        for rep in reps {133            let d = Design::full(rep, base);134            let s = series(base, &d.tile, top);135            if s.iter().all(|&x| x == 0) {136                println!("  {} | kept {} | no loop", d.name(), d.kept());137                continue;138            }139            let shown: Vec<String> = s.iter().take(8).map(|x| x.to_string()).collect();140            let found = oeis(&s);141            distinct.insert(s.clone(), found.starts_with("OEIS A"));142            println!(143                "  {} | kept {} | {} | {} | {found}",144                d.name(),145                d.kept(),146                shown.join(" "),147                describe(&s)148            );149        }150        let known = distinct.values().filter(|&&hit| hit).count();151        println!("  {} distinct nonzero loop sequences, in the OEIS by the first seven terms from the first nonzero one: {known} of {}", distinct.len(), distinct.len());152    }153}154155pub fn lengths(base: usize, top: usize, codes: &[u128]) {156    for &code in codes {157        let d = Design::full(code, base);158        println!(159            "NEW LOOPS BY LENGTH in block strands, void strands counted, {}",160            d.name()161        );162        let mut kept = crate::count::unit(true);163        let mut longest = 0;164        for level in 0..top {165            kept = crate::count::glue(base, &d.tile, &kept, level + 1 < top);166            let row: Vec<String> = kept167                .lengths168                .iter()169                .map(|(l, c)| format!("{l}:{c}"))170                .collect();171            longest = longest.max(kept.lengths.keys().last().copied().unwrap_or(0));172            println!("  gluing {level} to {}: {}", level + 1, row.join(" "));173        }174        println!(175            "  longest new loop over gluings 0 to {}: {longest} strands",176            top - 1177        );178    }179}180181fn fibonacci(n: i64) -> i128 {182    let (mut a, mut b) = (1i128, 0i128);183    for _ in 0..n + 1 {184        (a, b) = (b, a + b);185    }186    a187}188189fn gain_law(code: u128, number: usize, rule: usize, n: i64) -> Option<i128> {190    let p = |b: i128, e: i64| b.pow(e as u32);191    let first = |v: i128, f: i128| Some(if n == 0 { v } else { f });192    match (code, number, rule) {193        (7, 2, 2) => first(0, p(2, n) - 2),194        (11, 2, 2) => first(0, if n > 0 { p(2, n - 1) } else { 0 }),195        (9, 2, 2) => Some(1),196        (7, 3, 2) => Some(5 * p(3, n) - 7 * n as i128 - 5),197        (14, 3, 2) => first(1, if n > 0 { 2 * p(3, n - 1) } else { 0 }),198        (9, 3, 2) => Some(2 * p(3, n)),199        (6, 3, 2) => Some(p(2, n + 1)),200        (11, 3, 2) | (13, 3, 2) => Some(p(2, n + 1) - 2),201        (13, 3, 3) => first(0, fibonacci(2 * n - 3)),202        _ => None,203    }204}205206fn loop_law(code: u128, number: usize, rule: usize, n: i64) -> Option<i128> {207    let p = |b: i128, e: i64| b.pow(e as u32);208    let first = |f: &dyn Fn() -> i128| Some(if n == 0 { 0 } else { f() });209    match (code, number, rule) {210        (7, 2, 2) => first(&|| p(3, n - 1) - p(2, n) + 1),211        (11, 2, 2) => first(&|| p(3, n - 1) - p(2, n - 1)),212        (9, 2, 2) => Some(p(2, n) - 1),213        (7, 3, 2) => Some((p(8, n) - 1) / 7 - p(3, n) + n as i128 + 1),214        (14, 3, 2) => first(&|| 2 * p(5, n - 1) - p(3, n - 1)),215        (9, 3, 2) => Some(p(5, n) - p(3, n)),216        (6, 3, 2) => Some(p(4, n) - p(2, n)),217        (11, 3, 2) | (13, 3, 2) => Some((p(7, n) - 6 * p(2, n) + 5) / 15),218        _ => None,219    }220}221222pub fn gains() {223    println!(224        "NAMED DESIGNS: loops L by level and new loops J = L(level + 1) - kept L(level) per gluing"225    );226    let named = [227        (7u128, 2usize, 2usize),228        (11, 2, 2),229        (9, 2, 2),230        (7, 3, 2),231        (14, 3, 2),232        (9, 3, 2),233        (6, 3, 2),234        (11, 3, 2),235        (13, 3, 2),236        (13, 3, 3),237        (287, 3, 3),238    ];239    for (code, number, rule) in named {240        let d = Design::new(code, number, rule);241        let top = if number == 2 { 20 } else { 14 };242        let s = series(number, &d.tile, top);243        let k = d.kept() as i128;244        let j: Vec<i128> = (0..top)245            .map(|n| s[n + 1] as i128 - k * s[n] as i128)246            .collect();247        let gain: Vec<String> = j.iter().map(|x| x.to_string()).collect();248        let shown: Vec<String> = s.iter().take(12).map(|x| x.to_string()).collect();249        println!(250            "  {} | kept {k}\n    L {}\n    J {}\n    {} | {}",251            d.name(),252            shown.join(" "),253            gain.join(" "),254            describe(&s),255            oeis(&s)256        );257        if gain_law(code, number, rule, 0).is_some() {258            let held = (0..top)259                .filter(|&n| gain_law(code, number, rule, n as i64) == Some(j[n]))260                .count();261            println!(262                "    stated J law holds at {held} of {top} gluings 0 to {}",263                top - 1264            );265        }266        if loop_law(code, number, rule, 0).is_some() {267            let held = (0..=top)268                .filter(|&n| loop_law(code, number, rule, n as i64) == Some(s[n] as i128))269                .count();270            println!(271                "    stated L law holds at {held} of {} levels 0 to {top}",272                top + 1273            );274        }275    }276}277278static DUMP: std::sync::OnceLock<Option<String>> = std::sync::OnceLock::new();279280pub fn oeis(s: &[u128]) -> String {281    let dump = DUMP.get_or_init(|| {282        std::env::var("OEIS_STRIPPED")283            .ok()284            .and_then(|p| std::fs::read_to_string(p).ok())285    });286    let Some(dump) = dump else {287        return "OEIS not read".to_string();288    };289    let Some(k) = s.iter().position(|&x| x > 0) else {290        return String::new();291    };292    let window: Vec<String> = s[k..].iter().take(7).map(|x| x.to_string()).collect();293    let key = format!(",{},", window.join(","));294    let hits: Vec<&str> = dump295        .lines()296        .filter(|l| l.contains(&key))297        .filter_map(|l| l.split(' ').next())298        .take(3)299        .collect();300    if hits.is_empty() {301        "OEIS absent".to_string()302    } else {303        format!("OEIS {}", hits.join(" "))304    }305}