Clean up HeckeMiscPoly.jl (#1853)

src/HeckeMoreStuff.jl

@@ -262,26 +262,6 @@ function dot(a::Vector{<:NumFieldElem}, b::Vector{ZZRingElem})
   return d
 end
 
-function (R::zzModPolyRing)(g::QQPolyRingElem)
-  return fmpq_poly_to_nmod_poly(R, g)
-end
-
-function (R::fpPolyRing)(g::QQPolyRingElem)
-  return fmpq_poly_to_gfp_poly(R, g)
-end
-
-function (R::ZZModPolyRing)(g::QQPolyRingElem)
-  return fmpq_poly_to_fmpz_mod_poly(R, g)
-end
-
-function (R::FpPolyRing)(g::QQPolyRingElem)
-  return fmpq_poly_to_gfp_fmpz_poly(R, g)
-end
-
-function (R::FqPolyRing)(g::QQPolyRingElem)
-  return fmpq_poly_to_fq_default_poly(R, g)
-end
-
 function bits(x::ArbFieldElem)
   return ccall((:arb_bits, libflint), Int, (Ref{ArbFieldElem},), x)
 end

src/Nemo.jl

@@ -453,12 +453,12 @@ include("embedding/embedding.jl")
 include("Rings.jl")
 
 include("matrix.jl")
+include("poly.jl")
 
 include("Infinity.jl")
 
 include("HeckeMiscFiniteField.jl")
 include("HeckeMiscMatrix.jl")
-include("HeckeMiscPoly.jl")
 include("HeckeMoreStuff.jl")
 
 # More functionality for Julia types

src/flint/fmpq_poly.jl

@@ -837,6 +837,116 @@ promote_rule(::Type{QQPolyRingElem}, ::Type{Rational{T}}) where T <: Union{Int,
 
 (f::QQPolyRingElem)(a::Rational) = evaluate(f, QQFieldElem(a))
 
+###############################################################################
+#
+#   Conversion
+#
+###############################################################################
+
+function fmpq_poly_to_nmod_poly_raw!(r::zzModPolyRingElem, a::QQPolyRingElem)
+  ccall((:fmpq_poly_get_nmod_poly, libflint), Nothing,
+        (Ref{zzModPolyRingElem}, Ref{QQPolyRingElem}), r, a)
+  return r
+end
+
+function (R::zzModPolyRing)(g::QQPolyRingElem)
+  r = R()
+  fmpq_poly_to_nmod_poly_raw!(r, g)
+  return r
+end
+
+function fmpq_poly_to_gfp_poly_raw!(r::fpPolyRingElem, a::QQPolyRingElem)
+  ccall((:fmpq_poly_get_nmod_poly, libflint), Nothing,
+        (Ref{fpPolyRingElem}, Ref{QQPolyRingElem}), r, a)
+  return r
+end
+
+function (R::fpPolyRing)(g::QQPolyRingElem)
+  r = R()
+  fmpq_poly_to_gfp_poly_raw!(r, g)
+  return r
+end
+
+function fmpq_poly_to_fq_default_poly_raw!(r::FqPolyRingElem, a::QQPolyRingElem, t1::ZZPolyRingElem=ZZPolyRingElem(), t2::ZZRingElem=ZZRingElem())
+  ccall((:fmpq_poly_get_numerator, libflint), Nothing,
+        (Ref{ZZPolyRingElem}, Ref{QQPolyRingElem}), t1, a)
+  ccall((:fq_default_poly_set_fmpz_poly, libflint), Nothing,
+        (Ref{FqPolyRingElem}, Ref{ZZPolyRingElem}, Ref{FqField}),
+        r, t1, base_ring(parent(r)))
+  ccall((:fmpq_poly_get_denominator, libflint), Nothing,
+        (Ref{ZZRingElem}, Ref{QQPolyRingElem}), t2, a)
+  if !isone(t2)
+    ccall((:fq_default_poly_scalar_div_fq_default, libflint), Nothing,
+          (Ref{FqPolyRingElem}, Ref{FqPolyRingElem}, Ref{FqFieldElem}, Ref{FqField}),
+          r, r, coefficient_ring(r)(t2), coefficient_ring(r))
+  end
+  return r
+end
+
+function (R::FqPolyRing)(g::QQPolyRingElem)
+  r = R()
+  fmpq_poly_to_fq_default_poly_raw!(r, g)
+  return r
+end
+
+function fmpq_poly_to_fmpz_mod_poly_raw!(r::ZZModPolyRingElem, a::QQPolyRingElem, t1::ZZPolyRingElem=ZZPolyRingElem(), t2::ZZRingElem=ZZRingElem())
+  ccall((:fmpq_poly_get_numerator, libflint), Nothing,
+        (Ref{ZZPolyRingElem}, Ref{QQPolyRingElem}), t1, a)
+  ccall((:fmpz_mod_poly_set_fmpz_poly, libflint), Nothing,
+        (Ref{ZZModPolyRingElem}, Ref{ZZPolyRingElem}, Ref{fmpz_mod_ctx_struct}), r, t1, base_ring(parent(r)).ninv)
+  ccall((:fmpq_poly_get_denominator, libflint), Nothing,
+        (Ref{ZZRingElem}, Ref{QQPolyRingElem}), t2, a)
+  if !isone(t2)
+    res = ccall((:fmpz_invmod, libflint), Cint,
+                (Ref{ZZRingElem}, Ref{ZZRingElem}, Ref{ZZRingElem}),
+                t2, t2, modulus(base_ring(r)))
+    @assert res != 0
+    ccall((:fmpz_mod_poly_scalar_mul_fmpz, libflint), Nothing,
+          (Ref{ZZModPolyRingElem}, Ref{ZZModPolyRingElem}, Ref{ZZRingElem}, Ref{fmpz_mod_ctx_struct}),
+          r, r, t2, base_ring(parent(r)).ninv)
+  end
+  return r
+end
+
+function (R::ZZModPolyRing)(g::QQPolyRingElem)
+  r = R()
+  fmpq_poly_to_fmpz_mod_poly_raw!(r, g)
+  return r
+end
+
+function fmpq_poly_to_gfp_fmpz_poly_raw!(r::FpPolyRingElem, a::QQPolyRingElem, t1::ZZPolyRingElem=ZZPolyRingElem(), t2::ZZRingElem=ZZRingElem())
+  ccall((:fmpq_poly_get_numerator, libflint), Nothing,
+        (Ref{ZZPolyRingElem}, Ref{QQPolyRingElem}), t1, a)
+  ccall((:fmpz_mod_poly_set_fmpz_poly, libflint), Nothing,
+        (Ref{FpPolyRingElem}, Ref{ZZPolyRingElem}, Ref{fmpz_mod_ctx_struct}), r, t1, base_ring(parent(r)).ninv)
+  ccall((:fmpq_poly_get_denominator, libflint), Nothing,
+        (Ref{ZZRingElem}, Ref{QQPolyRingElem}), t2, a)
+  if !isone(t2)
+    res = ccall((:fmpz_invmod, libflint), Cint,
+                (Ref{ZZRingElem}, Ref{ZZRingElem}, Ref{ZZRingElem}),
+                t2, t2, modulus(base_ring(r)))
+    @assert res != 0
+    ccall((:fmpz_mod_poly_scalar_mul_fmpz, libflint), Nothing,
+          (Ref{FpPolyRingElem}, Ref{FpPolyRingElem}, Ref{ZZRingElem}, Ref{fmpz_mod_ctx_struct}),
+          r, r, t2, base_ring(parent(r)).ninv)
+  end
+  return r
+end
+
+function (R::FpPolyRing)(g::QQPolyRingElem)
+  r = R()
+  fmpq_poly_to_gfp_fmpz_poly_raw!(r, g)
+  return r
+end
+
+function (a::ZZPolyRing)(b::QQPolyRingElem)
+  (!isone(denominator(b))) && error("Denominator has to be 1")
+  z = a()
+  ccall((:fmpq_poly_get_numerator, libflint), Nothing,
+        (Ref{ZZPolyRingElem}, Ref{QQPolyRingElem}), z, b)
+  return z
+end
+
 ###############################################################################
 #
 #   Parent object call overloads
@@ -900,3 +1010,23 @@ function (a::QQPolyRing)(b::ZZPolyRingElem)
   z.parent = a
   return z
 end
+
+###############################################################################
+#
+#   Roots mod p
+#
+###############################################################################
+
+#TODO:
+# expand systematically for all finite fields
+# and for ZZRingElem/QQFieldElem poly
+
+# TODO: Can't we use the generic AbstractAlgebra `roots(::Field, ::PolyRingElem)`?
+# It doesn't multiply by the denominator, but if the denominator is 0 in R, then
+# then this function builds the zero polynomial and gives the root 0?
+function roots(R::T, f::QQPolyRingElem) where {T<:Union{fqPolyRepField, fpField}}
+  Rt, t = polynomial_ring(R, :t, cached=false)
+  fp = polynomial_ring(ZZ, cached=false)[1](f * denominator(f))
+  fpp = Rt(fp)
+  return roots(fpp)
+end

src/flint/fmpz.jl

@@ -3328,3 +3328,21 @@ end
 end
 
 using .BitsMod
+
+###############################################################################
+#
+#   Resultant
+#
+###############################################################################
+
+# From the flint docs:
+# Computes the resultant of f and g divided by d [...]. It is assumed that the
+# resultant is exactly divisible by d and the result res has at most nb bits.
+# This bypasses the computation of general bounds.
+function resultant(f::ZZPolyRingElem, g::ZZPolyRingElem, d::ZZRingElem, nb::Int)
+  z = ZZRingElem()
+  ccall((:fmpz_poly_resultant_modular_div, libflint), Nothing,
+        (Ref{ZZRingElem}, Ref{ZZPolyRingElem}, Ref{ZZPolyRingElem}, Ref{ZZRingElem}, Int),
+        z, f, g, d, nb)
+  return z
+end

src/flint/fmpz_mod.jl

@@ -32,6 +32,9 @@ function Base.hash(a::ZZModRingElem, h::UInt)
 end
 
 lift(a::ZZModRingElem) = data(a)
+lift(::ZZRing, a::ZZModRingElem) = lift(a)
+(::ZZRing)(a::ZZModRingElem) = lift(a)
+ZZRingElem(a::ZZModRingElem) = lift(a)
 
 function zero(R::ZZModRing)
   return ZZModRingElem(ZZRingElem(0), R)

src/flint/fmpz_mod_poly.jl

@@ -1012,3 +1012,23 @@ function (R::ZZModPolyRing)(f::ZZModPolyRingElem)
   parent(f) != R && error("Unable to coerce polynomial")
   return f
 end
+
+################################################################################
+#
+#  Rand
+#
+################################################################################
+
+@doc raw"""
+    rand(Rt::PolyRing{T}, n::Int) where T <: ResElem{ZZRingElem} -> PolyRingElem{T}
+
+Return a random polynomial of degree $n$.
+"""
+function Base.rand(Rt::PolyRing{T}, n::Int) where {T<:ResElem{ZZRingElem}}
+  f = Rt()
+  R = base_ring(Rt)
+  for i = 0:n
+    setcoeff!(f, i, rand(R))
+  end
+  return f
+end

src/flint/fmpz_poly.jl

@@ -66,6 +66,10 @@ function height(a::ZZPolyRingElem)
   return z
 end
 
+normalise(f::ZZPolyRingElem, ::Int) = degree(f) + 1
+
+set_length!(f::ZZPolyRingElem, ::Int) = nothing
+
 ###############################################################################
 #
 #   Similar and zero
@@ -896,6 +900,36 @@ promote_rule(::Type{ZZPolyRingElem}, ::Type{T}) where {T <: Integer} = ZZPolyRin
 
 promote_rule(::Type{ZZPolyRingElem}, ::Type{ZZRingElem}) = ZZPolyRingElem
 
+###############################################################################
+#
+#   Conversion
+#
+###############################################################################
+
+function fmpz_poly_to_nmod_poly_raw!(r::zzModPolyRingElem, a::ZZPolyRingElem)
+  ccall((:fmpz_poly_get_nmod_poly, libflint), Nothing,
+        (Ref{zzModPolyRingElem}, Ref{ZZPolyRingElem}), r, a)
+  return r
+end
+
+function (Rx::zzModPolyRing)(f::ZZPolyRingElem)
+  r = Rx()
+  fmpz_poly_to_nmod_poly_raw!(r, f)
+  return r
+end
+
+function fmpz_poly_to_gfp_poly_raw!(r::fpPolyRingElem, a::ZZPolyRingElem)
+  ccall((:fmpz_poly_get_nmod_poly, libflint), Nothing,
+        (Ref{fpPolyRingElem}, Ref{ZZPolyRingElem}), r, a)
+  return r
+end
+
+function (Rx::fpPolyRing)(f::ZZPolyRingElem)
+  r = Rx()
+  fmpz_poly_to_gfp_poly_raw!(r, f)
+  return r
+end
+
 ###############################################################################
 #
 #   Parent object call overloads
@@ -935,3 +969,57 @@ end
 (a::ZZPolyRing)(b::Vector{T}) where {T <: Integer} = a(map(ZZRingElem, b))
 
 (a::ZZPolyRing)(b::ZZPolyRingElem) = b
+
+###############################################################################
+#
+#  Sturm sequence
+#
+###############################################################################
+
+function _divide_by_content(f::ZZPolyRingElem)
+  p = primpart(f)
+  if sign(leading_coefficient(f)) == sign(leading_coefficient(p))
+    return p
+  else
+    return -p
+  end
+end
+
+function sturm_sequence(f::ZZPolyRingElem)
+  g = f
+  h = _divide_by_content(derivative(g))
+  seq = ZZPolyRingElem[g, h]
+  while true
+    r = _divide_by_content(pseudorem(g, h))
+    # r has the same sign as pseudorem(g, h)
+    # To get a pseudo remainder sequence for the Sturm sequence,
+    # we need r to be the pseudo remainder of |lc(b)|^(a - b + 1),
+    # so we need some adjustment. See
+    # https://en.wikipedia.org/wiki/Polynomial_greatest_common_divisor#Sturm_sequence_with_pseudo-remainders
+    if leading_coefficient(h) < 0 && isodd(degree(g) - degree(h) + 1)
+      r = -r
+    end
+    if r != 0
+      push!(seq, -r)
+      g, h = h, -r
+    else
+      break
+    end
+  end
+  return seq
+end
+
+###############################################################################
+#
+#   Mulhigh
+#
+###############################################################################
+
+function mulhigh_n(a::ZZPolyRingElem, b::ZZPolyRingElem, n::Int)
+  c = parent(a)()
+  #careful: as part of the interface, the coeffs 0 - (n-1) are random garbage
+  ccall((:fmpz_poly_mulhigh_n, libflint), Nothing,
+        (Ref{ZZPolyRingElem}, Ref{ZZPolyRingElem}, Ref{ZZPolyRingElem}, Cint),
+        c, a, b, n)
+  return c
+end

src/flint/fq_default.jl

@@ -956,3 +956,27 @@ function preimage(f::Generic.EuclideanRingResidueMap{ZZRing, FqField}, x)
   parent(x) !== codomain(f) && error("Not an element of the codomain")
   return lift(ZZ, x)
 end
+
+###############################################################################
+#
+#   Power detection
+#
+###############################################################################
+
+function is_power(a::Union{fpFieldElem, FpFieldElem, fqPolyRepFieldElem, FqPolyRepFieldElem, FqFieldElem}, m::Int)
+  if iszero(a)
+    return true, a
+  end
+  s = order(parent(a))
+  if gcd(s - 1, m) == 1
+    return true, a^invmod(ZZ(m), s - 1)
+  end
+  St, t = polynomial_ring(parent(a), :t, cached=false)
+  f = t^m - a
+  rt = roots(f)
+  if length(rt) > 0
+    return true, rt[1]
+  else
+    return false, a
+  end
+end