Fix false-positive in check_primality on degenerate specialization

check_primality(polys, extra_relations) substitutes the non-leader
variables of each generator with random integers, then tests whether the
characteristic polynomial of a generic multiplication operator on the
resulting zero-dimensional quotient is irreducible. The random
specialization could land on a point where a generator's leading
coefficient (w.r.t. its leader) vanishes, dropping that generator's degree
in its leader. When every generator collapses to degree 1, the quotient
ring becomes 1-dimensional, its charpoly is a degree-1 (hence trivially
irreducible) polynomial, and the routine reports the ideal as prime even
when it is not.

This is a soundness bug: for the io_projections issue-#132 case the ideal
is not prime, yet ~6/500 RNG seeds drove check_primality(proj) to return
true. It surfaced as a master red on the julia-pre (1.13.0-rc1) Core lane,
whose changed global-RNG stream happened to hit such a draw at
test/bodies/io_projections.jl:55 (`@test !check_primality(proj)`).

Fix: only accept random specializations that preserve every generator's
degree in its leader, i.e. that stay on the open set where the leading
coefficients do not vanish (the set the ideal is implicitly saturated at,
per the docstring). Resample otherwise, with a bounded retry. The
extra_relations are evaluated at the same point, matching the previous
behavior.

Verified on Julia 1.13.0-rc1: false-positive rate over seeds 0..499 drops
from 6/500 to 0/500; io_projections body 10/10 (was 9/1); positive case
check_primality(proj, [projection_poly]) 50/50 true. No regression on Julia
1.12: io_projections 10/10, check_primality_zerodim 5/5.

Co-Authored-By: Chris Rackauckas <accounts@chrisrackauckas.com>
Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>

Author: ChrisRackauckas

src/primality_check.jl

@@ -61,10 +61,35 @@ function check_primality(
     ring = parent(leaders[1])
 
     Rspec, vspec = Nemo.polynomial_ring(Nemo.QQ, [var_to_str(l) for l in leaders])
-    eval_point = [v in keys(polys) ? v : ring(rand(1:100)) for v in gens(ring)]
-    all_polys = vcat(collect(values(polys)), extra_relations)
-    zerodim_ideal =
-        collect(map(p -> parent_ring_change(evaluate(p, eval_point), Rspec), all_polys))
+    # Generic degree of each generator in its leader. Specializing the other
+    # variables must preserve these degrees; the ideal is considered on the
+    # open set where the leading coefficients (w.r.t. the leaders) do not
+    # vanish, so a specialization that drops a degree corresponds to leaving
+    # that set and would collapse the quotient ring, yielding a spurious result.
+    leader_degrees = Dict(l => Nemo.degree(p, l) for (l, p) in polys)
+
+    local zerodim_ideal
+    attempts = 0
+    while true
+        attempts += 1
+        eval_point = [v in keys(polys) ? v : ring(rand(1:100)) for v in gens(ring)]
+        specialized = Dict(l => evaluate(polys[l], eval_point) for l in leaders)
+        if all(Nemo.degree(specialized[l], l) == leader_degrees[l] for l in leaders)
+            all_polys = vcat(
+                [specialized[l] for l in leaders],
+                [evaluate(p, eval_point) for p in extra_relations],
+            )
+            zerodim_ideal =
+                collect(map(p -> parent_ring_change(p, Rspec), all_polys))
+            break
+        end
+        if attempts >= 100
+            error(
+                "Failed to find a generic specialization for the primality check " *
+                    "after $attempts attempts (the leading coefficients keep vanishing).",
+            )
+        end
+    end
 
     return check_primality_zerodim(zerodim_ideal)
 end