diff --git a/src/uu/factor/build.rs b/src/uu/factor/build.rs
index acc4d9a3f..3a1b7df94 100644
--- a/src/uu/factor/build.rs
+++ b/src/uu/factor/build.rs
@@ -26,7 +26,7 @@ use std::path::Path;
 use std::u64::MAX as MAX_U64;
 
 #[cfg(test)]
-use numeric::is_prime;
+use miller_rabin::is_prime;
 
 #[cfg(test)]
 #[path = "src/numeric.rs"]
diff --git a/src/uu/factor/src/factor.rs b/src/uu/factor/src/factor.rs
index afaf4f8cd..5a458a143 100644
--- a/src/uu/factor/src/factor.rs
+++ b/src/uu/factor/src/factor.rs
@@ -23,6 +23,7 @@ use std::fmt;
 use std::io::{stdin, BufRead};
 use std::ops;
 
+mod miller_rabin;
 mod numeric;
 mod rho;
 mod table;
diff --git a/src/uu/factor/src/miller_rabin.rs b/src/uu/factor/src/miller_rabin.rs
new file mode 100644
index 000000000..c42aee500
--- /dev/null
+++ b/src/uu/factor/src/miller_rabin.rs
@@ -0,0 +1,52 @@
+use crate::numeric::*;
+
+fn witness(mut a: u64, exponent: u64, m: u64) -> bool {
+    if a == 0 {
+        return false;
+    }
+
+    let mul = if m < 1 << 63 {
+        sm_mul as fn(u64, u64, u64) -> u64
+    } else {
+        big_mul as fn(u64, u64, u64) -> u64
+    };
+
+    if pow(a, m - 1, m, mul) != 1 {
+        return true;
+    }
+    a = pow(a, exponent, m, mul);
+    if a == 1 {
+        return false;
+    }
+    loop {
+        if a == 1 {
+            return true;
+        }
+        if a == m - 1 {
+            return false;
+        }
+        a = mul(a, a, m);
+    }
+}
+
+// Deterministic Miller-Rabin primality-checking algorithm, adapted to extract
+// (some) dividers; it will fail to factor strong pseudoprimes.
+pub(crate) fn is_prime(num: u64) -> bool {
+    if num < 2 {
+        return false;
+    }
+    if num % 2 == 0 {
+        return num == 2;
+    }
+    let mut exponent = num - 1;
+    while exponent & 1 == 0 {
+        exponent >>= 1;
+    }
+
+    // These witnesses detect all composites up to at least 2^64.
+    // Discovered by Jim Sinclair, according to http://miller-rabin.appspot.com
+    let witnesses = [2, 325, 9_375, 28_178, 450_775, 9_780_504, 1_795_265_022];
+    !witnesses
+        .iter()
+        .any(|&wit| witness(wit % num, exponent, num))
+}
diff --git a/src/uu/factor/src/numeric.rs b/src/uu/factor/src/numeric.rs
index d9e6f9441..7dfba41fb 100644
--- a/src/uu/factor/src/numeric.rs
+++ b/src/uu/factor/src/numeric.rs
@@ -78,7 +78,7 @@ pub fn big_mul(mut a: u64, mut b: u64, m: u64) -> u64 {
 }
 
 // computes a.pow(b) % m
-fn pow(mut a: u64, mut b: u64, m: u64, mul: fn(u64, u64, u64) -> u64) -> u64 {
+pub(crate) fn pow(mut a: u64, mut b: u64, m: u64, mul: fn(u64, u64, u64) -> u64) -> u64 {
     let mut result = 1;
     while b > 0 {
         if b & 1 != 0 {
@@ -89,53 +89,3 @@ fn pow(mut a: u64, mut b: u64, m: u64, mul: fn(u64, u64, u64) -> u64) -> u64 {
     }
     result
 }
-
-fn witness(mut a: u64, exponent: u64, m: u64) -> bool {
-    if a == 0 {
-        return false;
-    }
-
-    let mul = if m < 1 << 63 {
-        sm_mul as fn(u64, u64, u64) -> u64
-    } else {
-        big_mul as fn(u64, u64, u64) -> u64
-    };
-
-    if pow(a, m - 1, m, mul) != 1 {
-        return true;
-    }
-    a = pow(a, exponent, m, mul);
-    if a == 1 {
-        return false;
-    }
-    loop {
-        if a == 1 {
-            return true;
-        }
-        if a == m - 1 {
-            return false;
-        }
-        a = mul(a, a, m);
-    }
-}
-
-// uses deterministic (i.e., fixed witness set) Miller-Rabin test
-pub fn is_prime(num: u64) -> bool {
-    if num < 2 {
-        return false;
-    }
-    if num % 2 == 0 {
-        return num == 2;
-    }
-    let mut exponent = num - 1;
-    while exponent & 1 == 0 {
-        exponent >>= 1;
-    }
-
-    // These witnesses detect all composites up to at least 2^64.
-    // Discovered by Jim Sinclair, according to http://miller-rabin.appspot.com
-    let witnesses = [2, 325, 9_375, 28_178, 450_775, 9_780_504, 1_795_265_022];
-    !witnesses
-        .iter()
-        .any(|&wit| witness(wit % num, exponent, num))
-}
diff --git a/src/uu/factor/src/rho.rs b/src/uu/factor/src/rho.rs
index 949497737..947eac478 100644
--- a/src/uu/factor/src/rho.rs
+++ b/src/uu/factor/src/rho.rs
@@ -1,4 +1,4 @@
-use crate::Factors;
+use crate::{miller_rabin, Factors};
 use numeric::*;
 use rand::distributions::{Distribution, Uniform};
 use rand::rngs::SmallRng;
@@ -41,8 +41,11 @@ fn find_divisor(num: u64) -> u64 {
 
 pub(crate) fn factor(num: u64) -> Factors {
     let mut factors = Factors::new();
-    if num == 1 { return factors; }
-    if is_prime(num) {
+    if num == 1 {
+        return factors;
+    }
+
+    if miller_rabin::is_prime(num) {
         factors.push(num);
         return factors;
     }