Added new fn is_zero_det_probabilistic

src/Exports.jl

@@ -377,6 +377,7 @@ export is_unit
 export is_unknown
 export is_unimodular
 export is_upper_triangular
+export is_zero_det_probabilistic
 export is_zero_entry
 export is_zero_row
 export isequal

src/flint/fmpq_mat.jl

@@ -419,6 +419,58 @@ function det(x::QQMatrix)
   return z
 end
 
+# Test whether the determinant of a rational matrix is zero:
+# -- false means definitely not zero
+# -- true means "probably" zero (but hard to quantify how probable)
+# We use random starting points for the primes so that it is
+# harder to supply an input which will fool the implementation.
+# For initial speed we use small 22-bit primes, but after 3 iters
+# we switch to larger 61-bit primes (which are faster amortized).
+
+@doc raw"""
+    is_zero_det_probabilistic(M::QQMatrix; modulus_bitsize::Int = 100)
+
+Same as `is_zero(det(M))` but faster, though may _very rarely_ erroneously
+return `true`.  The kwarg `modulus_bitsize` must be between 20 and 1000;
+larger values decrease the probability of a false positive (but require
+more computation time).  Under a "uniformity assumption" the probability
+of a false positive is about `2^(-modulus_bitsize)`.
+"""
+function is_zero_det_probabilistic(M::QQMatrix; modulus_bitsize::Int = 100)
+  @req  is_square(M)  "matrix must be square"
+  @req ((modulus_bitsize >= 20) && (modulus_bitsize <= 1000))  "modulus_bitsize must be between 20 and 1000 (but bigger than about 250 is usually senseless)"
+  # Dispose of two trivial cases:
+  (nrows(M) == 0) && return false
+  (nrows(M) == 1) && return is_zero(M)
+  # General method starts here:
+  small_prime_size = 21
+  large_prime_size = 60
+  p = 2^small_prime_size + rand(1:2^(small_prime_size-1))  # works well on my computer
+  log2_modulus = 0.0
+  iter_count = 0
+  while log2_modulus < modulus_bitsize
+    iter_count += 1
+    # After 3 iters, we switch to larger primes (to make progress faster)
+    if iter_count == 4
+      p = 2^large_prime_size + rand(1:2^30)
+    end
+    p = next_prime(p)  # ?? better to do  next_prime(p+rand(1:2^16)) ??
+    Fp = Native.GF(p; cached=false, check=false)
+    try
+      # an attempt to use map_entries! (from AbstractAlgebra/src/Matrix.jl) caused _det to error (no idea why)
+      Mp = change_base_ring(Fp, M)  # this will throw if p divides the denominator of some matrix entry
+      if !is_zero(_det(Mp))
+        return false
+      end
+      log2_modulus += log2(p)  # only approximate, but that's fine
+    catch exc
+      (exc isa ErrorException && exc.msg == "Unable to coerce") || rethrow()
+      # The prime was bad, so we just ignore it (& the exception)
+    end
+  end
+  return true
+end
+
 ###############################################################################
 #
 #   Gram-Schmidt orthogonalisation

src/flint/fmpz_mat.jl

@@ -684,6 +684,52 @@ function det_given_divisor(x::ZZMatrix, d::Integer, proved=true)
 end
 
 
+# Test whether the determinant of an integer matrix is zero:
+# -- false means definitely not zero
+# -- true means "probably" zero (but hard to quantify how probable)
+# We use random starting points for the primes so that it is
+# harder to supply an input which will fool the implementation.
+# For initial speed we use small 22-bit primes, but after 3 iters
+# we switch to larger 61-bit primes (which are faster amortized).
+
+@doc raw"""
+    is_zero_det_probabilistic(M::ZZMatrix; modulus_bitsize::Int = 100)
+
+Same as `is_zero(det(M))` but faster, though may _very rarely_ erroneously
+return `true`.  The kwarg `modulus_bitsize` must be between 20 and 1000;
+larger values decrease the probability of a false positive (but require
+more computation time).  Under a "uniformity assumption" the probability
+of a false positive is about `2^(-modulus_bitsize)`.
+"""
+function is_zero_det_probabilistic(M::ZZMatrix; modulus_bitsize::Int = 100)
+  @req  is_square(M)  "matrix must be square"
+  @req ((modulus_bitsize >= 20) && (modulus_bitsize <= 1000))  "modulus_bitsize must be between 20 and 1000 (but bigger than about 250 is usually senseless)"
+  # Dispose of two trivial cases:
+  (nrows(M) == 0) && return false
+  (nrows(M) == 1) && return is_zero(M)
+  # General method starts here:
+  small_prime_size = 21
+  large_prime_size = 60
+  p = 2^small_prime_size + rand(1:2^(small_prime_size-1))  # works well on my computer
+  log2_modulus = 0.0
+  Mp = zero_matrix(Native.GF(2), nrows(M), ncols(M)) # workspace, to avoid reallocating
+  while log2_modulus < modulus_bitsize
+    # After 3 iters, we switch to larger primes (to make progress faster)
+    if log2_modulus > 3*small_prime_size && log2_modulus < 3*small_prime_size + 4
+      p = 2^large_prime_size + rand(1:2^30)
+    end
+    p = next_prime(p)  # ?? better to do  next_prime(p+rand(1:2^16)) ??
+    log2_modulus += log2(p)  # only approximate, but that's fine
+    Fp = Native.GF(p; cached=false, check=false)
+    Mp = map_entries!(Fp, Mp, M)  # non-allocating change_base_ring
+    if !is_zero(_det(Mp))
+      return false
+    end
+  end
+  return true
+end
+
+
 @doc raw"""
     hadamard_bound2(M::ZZMatrix)
 

test/flint/fmpq_mat-test.jl

@@ -520,6 +520,12 @@ end
   @test det(A) == 27
 end
 
+@testset "QQMatrix.is_zero_det_probabilistic" begin
+  sz = 250
+  M = matrix(QQ, sz,sz, (x->1//x).(1:sz*sz))
+  @time !is_zero_det_probabilistic(M)
+end
+
 @testset "QQMatrix.hilbert" begin
   S = matrix_space(QQ, 4, 4)
 

test/flint/fmpz_mat-test.jl

@@ -511,6 +511,19 @@ end
   @test det_given_divisor(A, ZZRingElem(9)) == 27
 end
 
+@testset "ZZMatrix.is_zero_det_probabilistic" begin
+  M00 = zero_matrix(ZZ, 0,0)
+  @test !is_zero_det_probabilistic(M00)
+
+  sz = 500
+  V = zero_matrix(ZZ, sz,sz);
+  for i in 1:sz  for j in 1:sz  V[i,j] = ZZ(j)^(i-1);  end;  end
+  @test !is_zero_det_probabilistic(V)
+
+  V_reduced = mod(V, ZZ(499))
+  @test is_zero_det_probabilistic(V_reduced)
+end
+
 @testset "ZZMatrix.hadamard" begin
   S = matrix_space(ZZ, 4, 4)