better col-limit

2 files changed, 55 insertions, 35 deletions

diff --git a/README.md b/README.md
index 63b2fc3..0f0dcfe 100644
--- a/README.md
+++ b/README.md
@@ -66,6 +66,8 @@ $${33551 \choose 12816} = {33552 \choose 12815} \approx 6.0 \times 10^{9687}$$
 
 Again we run into a memory bottleneck — searching this far required 14 GB of memory.
 
-We might suspect that other sporadic solutions (if they exist) have small values of $l$, seeing
-as the two known sporadic solutions have $l=5,8$. Searching with $l\leq 30, n,m < 3\times 10^7$ (arguments `col-limit 30000000 30`)
-uses 13GB of memory and yields no more solutions; and the same is true of searching with $l\leq 10, n,m < 10^8$.
+We can sidestep the large memory requirements almost entirely if we assume that sporadic solutions have
+very small values of $k$ (this is motivated by the two known sporadic solutions having $k=2$). Then we just
+go through each entry in the triangle, only storing one row in memory at a time, and do a binary search to check if
+the entry is ${n\choose k}$ for some small $k$. Using this method we can verify quickly and with just a few megabytes of memory
+that no more solutions exist with ${n\choose k}<10^{64}, n,m < 3\times 10^6, k \leq 10$ (arguments `col-limit 10`).
diff --git a/src/main.rs b/src/main.rs
index f7b4053..67840ab 100644
--- a/src/main.rs
+++ b/src/main.rs
@@ -2,13 +2,13 @@
 use num_bigint::BigUint;
 use std::ffi::OsString;
 use std::process::ExitCode;
-use std::cmp::{min, max};
+use std::cmp::min;
 
 // bnum seems to be slightly faster than ruint,
 // and 3x faster than uint.
 type UInt = bnum::types::U256;
 
-fn is_choose_r(mut k: UInt, r: u32) -> bool {
+fn find_choose_r(mut k: UInt, r: u32) -> Option<UInt> {
 	// multiply k by r! to avoid future divisions
 	for i in 1..=r {
 		k *= UInt::from(i);
@@ -32,10 +32,14 @@ fn is_choose_r(mut k: UInt, r: u32) -> bool {
 		if mid_falling_r < k {
 			lo = mid + UInt::from(1u8);
 		} else {
-			return true;
+			return Some(mid);
 		}
 	}
-	false
+	None
+}
+
+fn is_choose_r(k: UInt, r: u32) -> bool {
+	find_choose_r(k, r).is_some()
 }
 
 fn superscript(number: &str) -> String {
@@ -170,6 +174,7 @@ fn search_row_limit(row_limit: u32) {
 	for row in 0..row_limit {
 		if row > 0 && row % 2 == 0 {
 			pascal_row[row as usize / 2] = pascal_row[row as usize / 2 - 1].wrapping_mul(2);
+			entries.push(PascalEntry::new(pascal_row[row as usize / 2], row, (row / 2) as u16));
 		}
 		for col in (2..=(std::cmp::max(3, row) - 1) / 2).rev() {
 			pascal_row[col as usize] =
@@ -185,27 +190,48 @@ fn search_row_limit(row_limit: u32) {
 	find_duplicates_in(&mut entries);
 }
 
-fn search_col_limit(row_limit: u64, col_limit: u16) {
-	println!("Searching up to column {col_limit} in rows up to {row_limit}");
-	let mut pascal_row = vec![0u64; row_limit as usize / 2 + 1];
-	pascal_row[0] = 1;
-	let mut entries: Vec<PascalEntry> = vec![];
-	for row in 0..row_limit {
-		if row > 0 && row % 2 == 0 {
-			pascal_row[row as usize / 2] = pascal_row[row as usize / 2 - 1].wrapping_mul(2);
+fn search_col_limit(col_limit: u16) {
+	let mut pascal_row = vec![UInt::from(1u8), UInt::from(4u8), UInt::from(6u8)];
+	let mut range = usize::MAX;
+	let mut cutoff = UInt::MAX >> (2 * col_limit + 1);
+	for i in 2..=col_limit {
+		cutoff /= UInt::from(i);
+	}
+	println!("entry cutoff > 10^{}", cutoff.ilog10());
+	for row in 5usize..usize::MAX {
+		if row % (1 << 15) == 0 {
+			println!("row = {row}");
 		}
-		for col in (2..=min(u64::from(col_limit), (max(3, row) - 1) / 2)).rev() {
-			pascal_row[col as usize] =
-				pascal_row[col as usize].wrapping_add(pascal_row[col as usize - 1]);
-			entries.push(PascalEntry::new(pascal_row[col as usize], row, col as u16));
+		for col in (2..min(row / 2 + 1, range)).rev() {
+			if 2 * col == row {
+				pascal_row.push(pascal_row[row / 2 - 1] * UInt::from(2u8));
+			} else if col == 2 {
+				pascal_row[col] += UInt::from(row - 1);
+			} else {
+				let prev = pascal_row[col - 1];
+				pascal_row[col] += prev;
+			}
+			if pascal_row[col] > cutoff {
+				// getting dicey
+				println!("abandoning column {col}");
+				range = col;
+				continue;
+			}
+			if col <= 4 {
+				// we already know the solutions with both columns ≤ 4
+				continue;
+			}
+			for col2 in 2..min(col, col_limit as usize + 1) {
+				if let Some(x) = find_choose_r(pascal_row[col], col2 as u32) {
+					println!("({row} choose {col}) = ({x} choose {col2})");
+				}
+			}
+		}
+		if range <= 2 {
+			println!("can't go any further without increasing size of UInt");
+			return;
 		}
-		pascal_row[1] = row;
 	}
-	println!(
-		"memory needed = {}MiB",
-		(entries.len() * size_of::<PascalEntry>()) >> 20
-	);
-	find_duplicates_in(&mut entries);
 }
 
 fn main() -> ExitCode {
@@ -249,23 +275,15 @@ fn main() -> ExitCode {
 			search_row_limit(row_limit);
 		}
 		Some("col-limit") => {
-			let row_limit: Option<u64> = match args.get(2) {
-				Some(s) => s.clone().into_string().ok().and_then(|x| x.parse().ok()),
-				None => Some(100_000),
-			};
-			let col_limit: Option<u16> = match args.get(3) {
+			let col_limit: Option<u16> = match args.get(2) {
 				Some(s) => s.clone().into_string().ok().and_then(|x| x.parse().ok()),
 				None => Some(30),
 			};
-			let Some(row_limit) = row_limit else {
-				eprintln!("row limit must be a nonnegative integer");
-				return ExitCode::FAILURE;
-			};
 			let Some(col_limit) = col_limit else {
 				eprintln!("column limit must be a nonnegative integer < 65536");
 				return ExitCode::FAILURE;
 			};
-			search_col_limit(row_limit, col_limit);
+			search_col_limit(col_limit);
 		}
 		_ => {
 			eprintln!("Bad command: {:?}", args[1]);