[deps] cargo update

- don't clone primes for every thread, instead put it behind an Arc
  (atomic reference counted)
- print the pure calculation time at the end of the program
- don't show thread debug prints

Cargo.lock

@@ -3,184 +3,23 @@
 version = 4
 
 [[package]]
-name = "anstream"
+name = "either"
-version = "1.0.0"
+version = "1.15.0"
 source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "824a212faf96e9acacdbd09febd34438f8f711fb84e09a8916013cd7815ca28d"
+checksum = "48c757948c5ede0e46177b7add2e67155f70e33c07fea8284df6576da70b3719"
-dependencies = [
- "anstyle",
- "anstyle-parse",
- "anstyle-query",
- "anstyle-wincon",
- "colorchoice",
- "is_terminal_polyfill",
- "utf8parse",
-]
-
-[[package]]
-name = "anstyle"
-version = "1.0.14"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "940b3a0ca603d1eade50a4846a2afffd5ef57a9feac2c0e2ec2e14f9ead76000"
-
-[[package]]
-name = "anstyle-parse"
-version = "1.0.0"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "52ce7f38b242319f7cabaa6813055467063ecdc9d355bbb4ce0c68908cd8130e"
-dependencies = [
- "utf8parse",
-]
-
-[[package]]
-name = "anstyle-query"
-version = "1.1.5"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "40c48f72fd53cd289104fc64099abca73db4166ad86ea0b4341abe65af83dadc"
-dependencies = [
- "windows-sys",
-]
-
-[[package]]
-name = "anstyle-wincon"
-version = "3.0.11"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "291e6a250ff86cd4a820112fb8898808a366d8f9f58ce16d1f538353ad55747d"
-dependencies = [
- "anstyle",
- "once_cell_polyfill",
- "windows-sys",
-]
-
-[[package]]
-name = "clap"
-version = "4.6.1"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "1ddb117e43bbf7dacf0a4190fef4d345b9bad68dfc649cb349e7d17d28428e51"
-dependencies = [
- "clap_builder",
- "clap_derive",
-]
-
-[[package]]
-name = "clap_builder"
-version = "4.6.0"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "714a53001bf66416adb0e2ef5ac857140e7dc3a0c48fb28b2f10762fc4b5069f"
-dependencies = [
- "anstream",
- "anstyle",
- "clap_lex",
- "strsim",
-]
-
-[[package]]
-name = "clap_derive"
-version = "4.6.1"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "f2ce8604710f6733aa641a2b3731eaa1e8b3d9973d5e3565da11800813f997a9"
-dependencies = [
- "heck",
- "proc-macro2",
- "quote",
- "syn",
-]
-
-[[package]]
-name = "clap_lex"
-version = "1.1.0"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "c8d4a3bb8b1e0c1050499d1815f5ab16d04f0959b233085fb31653fbfc9d98f9"
-
-[[package]]
-name = "colorchoice"
-version = "1.0.5"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "1d07550c9036bf2ae0c684c4297d503f838287c83c53686d05370d0e139ae570"
 
 [[package]]
 name = "hcn"
-version = "2.0.0"
+version = "1.0.0-multithread"
 dependencies = [
- "clap",
+ "itertools",
 ]
 
 [[package]]
-name = "heck"
+name = "itertools"
-version = "0.5.0"
+version = "0.14.0"
 source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "2304e00983f87ffb38b55b444b5e3b60a884b5d30c0fca7d82fe33449bbe55ea"
+checksum = "2b192c782037fadd9cfa75548310488aabdbf3d2da73885b31bd0abd03351285"
-
-[[package]]
-name = "is_terminal_polyfill"
-version = "1.70.2"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "a6cb138bb79a146c1bd460005623e142ef0181e3d0219cb493e02f7d08a35695"
-
-[[package]]
-name = "once_cell_polyfill"
-version = "1.70.2"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "384b8ab6d37215f3c5301a95a4accb5d64aa607f1fcb26a11b5303878451b4fe"
-
-[[package]]
-name = "proc-macro2"
-version = "1.0.106"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "8fd00f0bb2e90d81d1044c2b32617f68fcb9fa3bb7640c23e9c748e53fb30934"
 dependencies = [
- "unicode-ident",
+ "either",
-]
-
-[[package]]
-name = "quote"
-version = "1.0.45"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "41f2619966050689382d2b44f664f4bc593e129785a36d6ee376ddf37259b924"
-dependencies = [
- "proc-macro2",
-]
-
-[[package]]
-name = "strsim"
-version = "0.11.1"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "7da8b5736845d9f2fcb837ea5d9e2628564b3b043a70948a3f0b778838c5fb4f"
-
-[[package]]
-name = "syn"
-version = "2.0.117"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "e665b8803e7b1d2a727f4023456bbbbe74da67099c585258af0ad9c5013b9b99"
-dependencies = [
- "proc-macro2",
- "quote",
- "unicode-ident",
-]
-
-[[package]]
-name = "unicode-ident"
-version = "1.0.24"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "e6e4313cd5fcd3dad5cafa179702e2b244f760991f45397d14d4ebf38247da75"
-
-[[package]]
-name = "utf8parse"
-version = "0.2.2"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "06abde3611657adf66d383f00b093d7faecc7fa57071cce2578660c9f1010821"
-
-[[package]]
-name = "windows-link"
-version = "0.2.1"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "f0805222e57f7521d6a62e36fa9163bc891acd422f971defe97d64e70d0a4fe5"
-
-[[package]]
-name = "windows-sys"
-version = "0.61.2"
-source = "registry+https://github.com/rust-lang/crates.io-index"
-checksum = "ae137229bcbd6cdf0f7b80a31df61766145077ddf49416a728b02cb3921ff3fc"
-dependencies = [
- "windows-link",
 ]

Cargo.toml

@@ -1,10 +1,10 @@
 [package]
 name = "hcn"
-version = "2.0.0"
+version = "1.0.0-multithread"
-edition = "2024"
+edition = "2021"
 
 [dependencies]
-clap = { version = "4", features = ["derive"] }
+itertools = "0.14"
 
 [lints.rust]
 unsafe_code = "forbid"

README.md

@@ -1,72 +0,0 @@
-# hcn
-
-Find highly composite numbers up to a search limit.
-
-```bash
-cargo run --release -- [SEARCH_LIMIT]
-```
-
-If no limit is given, the default is `1_000_000_000`.
-Use `--help` to see the generated command-line help and `--version` to print
-the package version.
-
-## Structure
-
-- `src/main.rs`
-  - Command-line interface: parses the optional `SEARCH_LIMIT` argument with
-    `clap` and rejects non-positive limits before the search starts.
-  - Search setup: defines the default limit and the fixed prime list used for
-    candidate generation.
-  - Candidate generation: recursively builds only exponent sequences that can
-    produce highly composite record candidates.
-  - Record filtering: sorts candidates and prints each new divisor-count record.
-
-## Algorithm
-
-A number with prime factorization
-
-```text
-n = 2^a * 3^b * 5^c * ...
-```
-
-has
-
-```text
-(a + 1) * (b + 1) * (c + 1) * ...
-```
-
-divisors.
-
-The exact primes do not affect the divisor count; only the exponents do.
-For a fixed exponent list, the smallest possible number is made by putting the
-largest exponent on the smallest prime:
-
-```text
-a >= b >= c >= ...
-```
-
-If a larger exponent appears on a larger prime, swapping those two exponents
-keeps the divisor count unchanged but makes the number smaller.
-
-So a highly composite number cannot have that shape. If it did, the smaller
-swapped number would already have the same number of divisors, so the original
-number would not be the first record.
-
-That means record candidates are only numbers like:
-
-```text
-2^a * 3^b * 5^c * ...  where a >= b >= c >= ...
-```
-
-All other numbers are redundant for finding new records.
-
-The program:
-
-1. Recursively generates numbers from the first primes.
-2. Only tries exponent sequences where each exponent is no larger than the
-   previous one.
-3. Stops a branch as soon as the number exceeds the search limit.
-4. Computes the divisor count directly from the exponents.
-5. Sorts candidates by number and prints each new divisor-count record.
-
-This avoids factoring every number in the range.

src/main.rs

@@ -1,278 +1,126 @@
-use std::{num::NonZeroU64, time::Instant};
+mod sieve;
 
-use clap::Parser;
+use std::{
+    collections::HashMap,
+    env,
+    sync::Arc,
+    thread::{self, available_parallelism},
+    time::Instant,
+};
 
-const DEFAULT_SEARCH_LIMIT: u64 = 1_000_000_000;
+use itertools::Itertools as _;
-const FIRST_PRIMES: &[u64] = &[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31];
+use sieve::get_primes;
 
-#[derive(Debug, Parser)]
+const MAX_SIEVED_PRIMES: usize = 100_000_000;
-#[command(version, about = "Find highly composite numbers up to a search limit")]
+
-struct Cli {
+fn count_prime_factor_exponent(prime: u64, mut nr: u64) -> u64 {
-    #[arg(value_name = "SEARCH_LIMIT", default_value_t = NonZeroU64::new(DEFAULT_SEARCH_LIMIT).expect("default search limit is non-zero"))]
+    let mut count = 0;
-    search_limit: NonZeroU64,
+
+    while nr % prime == 0 {
+        nr /= prime;
+        count += 1;
+    }
+
+    count
 }
 
-#[derive(Clone, Copy, Debug, Eq, PartialEq)]
+fn prime_factors(nr: u64, primes: &[u64]) -> u64 {
-struct Candidate {
+    let mut num_teilers = 1;
-    nr: u64,
-    divisors: u64,
-}
 
-fn max_exponent(limit: u64) -> u64 {
+    // safety: nr <= 2.pow(53) must be fulfilled,
-    let mut exponent = 0;
+    //         otherwise we have a precision loss
-    let mut nr = 1_u64;
+    assert!(nr <= 2u64.pow(53));
+    #[allow(
+        clippy::cast_sign_loss,
+        clippy::cast_possible_truncation,
+        clippy::cast_precision_loss
+    )]
+    let nr_sqrt = (nr as f64).sqrt().ceil() as u64;
 
-    while let Some(next_nr) = nr.checked_mul(2) {
+    for prime in primes {
-        if next_nr > limit {
+        if *prime > nr_sqrt {
             break;
         }
 
-        nr = next_nr;
+        let count = count_prime_factor_exponent(*prime, nr);
-        exponent += 1;
+        if count > 0 {
+            num_teilers *= count + 1;
+        }
     }
 
-    exponent
+    num_teilers
 }
 
-fn collect_hcn_candidates(
+fn the_thread(mut start: u64, end: u64, primes: &[u64]) -> HashMap<u64, u64> {
-    limit: u64,
+    let mut max_teilers = 0;
-    prime_idx: usize,
+    let mut max_teilers_hashmap = HashMap::new();
-    max_exponent: u64,
-    nr: u64,
-    divisors: u64,
-    candidates: &mut Vec<Candidate>,
-) {
-    let Some(&prime) = FIRST_PRIMES.get(prime_idx) else {
-        return;
-    };
 
-    let mut next_nr = nr;
+    let mut step_size = 10;
 
-    for exponent in 1..=max_exponent {
+    if start == 0 {
-        let Some(multiplied_nr) = next_nr.checked_mul(prime) else {
+        start = 2;
-            break;
+        step_size = 2;
-        };
-
-        if multiplied_nr > limit {
-            break;
     }
 
-        next_nr = multiplied_nr;
+    for nr in (start..=end).step_by(step_size) {
-        let next_divisors = divisors * (exponent + 1);
+        let teilers = prime_factors(nr, primes);
-
+        if teilers > max_teilers {
-        candidates.push(Candidate {
+            max_teilers = teilers;
-            nr: next_nr,
+            max_teilers_hashmap.insert(nr, teilers);
-            divisors: next_divisors,
-        });
-
-        collect_hcn_candidates(
-            limit,
-            prime_idx + 1,
-            exponent,
-            next_nr,
-            next_divisors,
-            candidates,
-        );
         }
+    }
+
+    max_teilers_hashmap
 }
 
-fn has_enough_primes(limit: u64) -> bool {
+fn calculate_chunk_bounds(i: usize, num_threads: usize, max_nr: usize) -> (usize, usize) {
-    let mut primorial = 1_u64;
+    let chunk_size = max_nr / num_threads;
-
+    let start = i * chunk_size;
-    for &prime in FIRST_PRIMES {
+    let end = ((i + 1) * chunk_size).min(max_nr);
-        let Some(next_primorial) = primorial.checked_mul(prime) else {
+    (start, end)
-            return true;
-        };
-
-        if next_primorial > limit {
-            return true;
-        }
-
-        primorial = next_primorial;
-    }
-
-    false
-}
-
-fn hcn_records(limit: u64) -> Vec<Candidate> {
-    assert!(
-        has_enough_primes(limit),
-        "not enough primes for search limit"
-    );
-
-    let mut candidates = Vec::new();
-    collect_hcn_candidates(limit, 0, max_exponent(limit), 1, 1, &mut candidates);
-    candidates.sort_unstable_by_key(|candidate| candidate.nr);
-
-    let mut records = Vec::new();
-    let mut max_divisors = 0;
-
-    for candidate in candidates {
-        if candidate.divisors > max_divisors {
-            records.push(candidate);
-            max_divisors = candidate.divisors;
-        }
-    }
-
-    records
 }
 
 fn main() {
-    let search_limit = Cli::parse().search_limit.get();
+    let max_nr = env::args()
+        .nth(1)
+        .expect("Usage: hcn <max_nr>")
+        .parse::<usize>()
+        .expect("Invalid max_nr");
 
     let start = Instant::now();
-    let records = hcn_records(search_limit);
+    println!("Precalculating primes...");
+    let primes = get_primes(MAX_SIEVED_PRIMES);
+    println!("{} primes. Took {:?}", primes.len(), start.elapsed());
 
-    for Candidate { nr, divisors } in records {
+    let primes = Arc::new(primes);
-        println!("{nr}: {divisors} ({:?} since start)", start.elapsed());
+
-    }
+    #[allow(clippy::unwrap_used)]
-}
+    let num_threads = available_parallelism().unwrap().get();
-
+    let mut threads = Vec::with_capacity(num_threads);
-#[cfg(test)]
+
-mod tests {
+    let now = Instant::now();
-    use clap::Parser;
+    for i in 0..num_threads {
-
+        let (start, end) = calculate_chunk_bounds(i, num_threads, max_nr);
-    use super::{Candidate, Cli, DEFAULT_SEARCH_LIMIT, hcn_records, max_exponent};
+        let primes = primes.clone();
-
+        threads.push(thread::spawn(move || {
-    fn count_divisors(mut nr: u64) -> u64 {
+            the_thread(start as u64, end as u64, &primes)
-        let mut divisor_count = 1;
+        }));
-        let mut divisor = 2;
-
-        while divisor * divisor <= nr {
-            let mut exponent = 0;
-
-            while nr.is_multiple_of(divisor) {
-                nr /= divisor;
-                exponent += 1;
-            }
-
-            if exponent > 0 {
-                divisor_count *= exponent + 1;
-            }
-
-            divisor += if divisor == 2 { 1 } else { 2 };
-        }
-
-        if nr > 1 {
-            divisor_count *= 2;
-        }
-
-        divisor_count
-    }
-
-    fn brute_force_hcn_records(limit: u64) -> Vec<Candidate> {
-        let mut records = Vec::new();
-        let mut max_divisors = 0;
-
-        for nr in 2..=limit {
-            let divisors = count_divisors(nr);
-
-            if divisors > max_divisors {
-                records.push(Candidate { nr, divisors });
-                max_divisors = divisors;
-            }
-        }
-
-        records
-    }
-
-    #[test]
-    fn max_exponent_counts_powers_of_two_up_to_limit() {
-        assert_eq!(max_exponent(1), 0);
-        assert_eq!(max_exponent(2), 1);
-        assert_eq!(max_exponent(3), 1);
-        assert_eq!(max_exponent(16), 4);
-    }
-
-    #[test]
-    fn finds_highly_composite_records_up_to_one_thousand() {
-        let records = hcn_records(1_000);
-
-        assert_eq!(
-            records,
-            vec![
-                Candidate { nr: 2, divisors: 2 },
-                Candidate { nr: 4, divisors: 3 },
-                Candidate { nr: 6, divisors: 4 },
-                Candidate {
-                    nr: 12,
-                    divisors: 6,
-                },
-                Candidate {
-                    nr: 24,
-                    divisors: 8,
-                },
-                Candidate {
-                    nr: 36,
-                    divisors: 9,
-                },
-                Candidate {
-                    nr: 48,
-                    divisors: 10,
-                },
-                Candidate {
-                    nr: 60,
-                    divisors: 12,
-                },
-                Candidate {
-                    nr: 120,
-                    divisors: 16,
-                },
-                Candidate {
-                    nr: 180,
-                    divisors: 18,
-                },
-                Candidate {
-                    nr: 240,
-                    divisors: 20,
-                },
-                Candidate {
-                    nr: 360,
-                    divisors: 24,
-                },
-                Candidate {
-                    nr: 720,
-                    divisors: 30,
-                },
-                Candidate {
-                    nr: 840,
-                    divisors: 32,
-                },
-            ],
-        );
-    }
-
-    #[test]
-    fn matches_brute_force_records_up_to_ten_thousand() {
-        assert_eq!(hcn_records(10_000), brute_force_hcn_records(10_000));
-    }
-
-    #[test]
-    fn parses_default_search_limit_without_arg() {
-        assert_eq!(
-            Cli::try_parse_from(["hcn"])
-                .expect("CLI should parse without an explicit search limit")
-                .search_limit
-                .get(),
-            DEFAULT_SEARCH_LIMIT
-        );
-    }
-
-    #[test]
-    fn parses_search_limit_arg() {
-        assert_eq!(
-            Cli::try_parse_from(["hcn", "42000"])
-                .expect("CLI should parse an explicit search limit")
-                .search_limit
-                .get(),
-            42_000
-        );
-    }
-
-    #[test]
-    fn rejects_invalid_search_limit_arg() {
-        assert!(Cli::try_parse_from(["hcn", "nope"]).is_err());
-        assert!(Cli::try_parse_from(["hcn", "0"]).is_err());
-        assert!(Cli::try_parse_from(["hcn", "100", "200"]).is_err());
     }
+
+    let mut results = HashMap::new();
+    for thread in threads {
+        let thread_result = thread.join().expect("Thread failed to join (panicked?)");
+        results.extend(thread_result);
+    }
+
+    let mut max_teilers = 0;
+
+    for (nr, teilers) in results.into_iter().sorted() {
+        if teilers > max_teilers {
+            println!("{nr}: {teilers}");
+            max_teilers = teilers;
+        }
+    }
+
+    println!("Took {:?}", now.elapsed());
 }

src/sieve.rs

@@ -0,0 +1,50 @@
+pub fn get_primes(max_nr: usize) -> Vec<u64> {
+    let mut sieve = vec![true; max_nr + 1];
+
+    let mut primes = Vec::with_capacity(1_000_000);
+
+    sieve[0] = false;
+    sieve[1] = false;
+
+    let mut nr = 2;
+    primes.push(2);
+
+    loop {
+        check_off_multiples_of_nr(nr, &mut sieve);
+
+        if find_next_prime(&mut nr, &mut sieve).is_none() {
+            break;
+        }
+
+        primes.push(nr);
+    }
+
+    primes
+}
+
+fn find_next_prime(nr: &mut u64, sieve: &mut [bool]) -> Option<u64> {
+    let start_idx = usize::try_from(*nr).ok()? + 1;
+
+    if let Some((idx, _)) = sieve
+        .iter()
+        .enumerate()
+        .skip(start_idx)
+        .find(|(_, &is_prime)| is_prime)
+    {
+        *nr = idx as u64;
+        return Some(*nr);
+    }
+
+    None
+}
+
+#[allow(clippy::unwrap_used)]
+fn check_off_multiples_of_nr(nr: u64, sieve: &mut [bool]) {
+    let nr = usize::try_from(nr).unwrap();
+
+    sieve
+        .iter_mut()
+        .skip(nr.checked_mul(2).unwrap())
+        .step_by(nr)
+        .for_each(|nr| *nr = false);
+}