summaryrefslogtreecommitdiff
path: root/src/main.rs
diff options
context:
space:
mode:
Diffstat (limited to 'src/main.rs')
-rw-r--r--src/main.rs134
1 files changed, 96 insertions, 38 deletions
diff --git a/src/main.rs b/src/main.rs
index b6ce57c..988a6bd 100644
--- a/src/main.rs
+++ b/src/main.rs
@@ -3,12 +3,15 @@ use num::{
Signed,
};
+use berlekamp::Poly;
+
/// A simple conversion function from anything that can be converted
/// to a big integer to a big rational number.
fn conversion(n: impl ToBigInt) -> BigInt {
n.to_bigint().unwrap()
}
+#[allow(unused)]
fn print_poly(poly: &[BigInt]) {
struct Monomial(usize);
@@ -80,17 +83,39 @@ fn main() {
.parse::<usize>()
.unwrap();
- println!("poly is:");
- print_poly(&poly);
+ if std::env::args().len() >= 4 {
+ // Search for primes modulo which the polynomial is irreducible instead.
- println!("prime is {p}");
+ let primes = generate_primes(p);
- // let poly1 = [12, 12, 0, 1].map(conversion);
- // let poly2 = [-4, 1].map(conversion);
+ println!("Searching {} primes up to {p}...", primes.len());
- // let composition = berlekamp::composition(&poly1, &poly2);
+ for p in primes {
+ let factors = berlekamp::factors(&poly, p);
+
+ let total = factors.iter().fold(0usize, |sum, (deg, _)| sum + deg);
+
+ if total == 1 {
+ println!("The polynomial {} is irreducible mod {p}", Poly(&poly));
+
+ return;
+ }
+ }
+
+ println!("The polynomial is reducible modulo primes up to {p}.");
+
+ return;
+ }
- use berlekamp::Poly;
+ // println!("poly is:");
+ // print_poly(&poly);
+
+ // println!("prime is {p}");
+
+ // let poly1 = [-20, -21, 0, 1].map(conversion);
+ // let poly2 = [1, 3].map(conversion);
+
+ // let composition = berlekamp::composition(&poly1, &poly2);
// println!(
// "the composition of {} and {} is {}",
@@ -107,44 +132,36 @@ fn main() {
let factors = berlekamp::factors(&poly, p);
- match factors {
- Ok(factors) => {
- println!("the irreducible factors are:");
+ print!("{} ≡ ", Poly(&poly));
- for factor in factors {
- print_poly(&factor);
- }
- }
- Err(gcd) => {
- println!(
- "the gcd of the polynomial with its derivative is {}",
- Poly(&gcd)
- );
+ let mut total = 0usize;
- if gcd.len() == poly.len() {
- // the derivative is zero
- println!("the polynomial is a {p}-th power of some polynomial");
+ for (index, (mul, f)) in factors.into_iter().enumerate() {
+ total += mul;
- return;
+ match mul {
+ 0 => {
+ dbg!("zero multiplicity?");
}
-
- let (quotient, _) = berlekamp::poly_division(&poly, &gcd, conversion(p));
-
- println!(
- "the quotient of the polynomial by the above gcd is {}",
- Poly(&quotient)
- );
-
- let factors =
- berlekamp::factors(&quotient, p).expect("this quotient should be square-free");
-
- println!("the irreducible factors of the quotient, and hence of the polynomial, are:");
-
- for factor in factors {
- print_poly(&factor);
+ 1 => {
+ print!("{}({})", if index > 0 { " * " } else { "" }, Poly(&f));
+ }
+ _ => {
+ print!(
+ "{}({})^{}",
+ if index > 0 { " * " } else { "" },
+ Poly(&f),
+ mul
+ );
}
}
}
+
+ println!(" (mod {p})");
+
+ if total == 1 {
+ println!("the polynomial is irreducible");
+ }
}
#[allow(unused)]
@@ -315,3 +332,44 @@ fn print_matrix<const X: usize, const Y: usize, const T: usize>(matrix: &[BigInt
}
println!("]");
}
+
+#[allow(unused)]
+/// Return a vector of primes <= BOUND.
+///
+/// This is naïve and slow, so BOUND should not be too large.
+fn generate_primes(bound: usize) -> Vec<usize> {
+ if bound <= 1 {
+ return Vec::new();
+ }
+
+ let sqrt = (bound as f64).sqrt() as usize;
+
+ let mut result: Vec<usize> = Vec::with_capacity(bound);
+
+ let mut records: Vec<bool> = std::iter::repeat(true).take(bound).collect();
+
+ records[0] = false;
+ records[1] = false;
+
+ for i in 2..=(sqrt + 1) {
+ if !matches!(records.get(i), Some(true)) {
+ continue;
+ }
+
+ let mut multiple = i + i;
+
+ while let Some(rec) = records.get_mut(multiple) {
+ *rec = false;
+ multiple += i;
+ }
+ }
+
+ result.extend(
+ records
+ .iter()
+ .enumerate()
+ .filter_map(|(n, rec)| rec.then_some(n)),
+ );
+
+ result
+}