First stab at more functionality for Laurent polynomial rings (#2448)

Project.toml

@@ -23,7 +23,7 @@ UUIDs = "cf7118a7-6976-5b1a-9a39-7adc72f591a4"
 cohomCalg_jll = "5558cf25-a90e-53b0-b813-cadaa3ae7ade"
 
 [compat]
-AbstractAlgebra = "0.30.7"
+AbstractAlgebra = "0.30.9"
 AlgebraicSolving = "0.3.0"
 DocStringExtensions = "0.8, 0.9"
 GAP = "0.9.4"

experimental/Laurent/src/Laurent.jl

@@ -0,0 +1,213 @@
+module Laurent
+
+using ..Oscar
+
+include("Types.jl")
+
+import .AbstractAlgebra.Generic: LaurentMPolyWrapRing, LaurentMPolyWrap
+import .AbstractAlgebra: LaurentMPolyRing, LaurentMPolyRingElem
+import ..Oscar: MPolyAnyMap, ideal, image, preimage, base_ring, in, gens, hom, domain, codomain, +
+
+# Ideals and maps for multivariate Laurent polynomial rings
+# 
+# Tommy Hofmann:
+# If R = K[x1^+-,...,xn^+-] is a such a ring, then this is isomorphic as a K-algebra
+# to the affine algebra Raff = K[x1,...,xn,y1,...,yn]/(x1y1 - 1,...,xnyn - 1)
+# We piggyback on quotient rings for all serious computations. More precisely
+#
+# g = _polyringquo(R)
+#
+# constructs a map g : R -> Raff, which supports
+# - image(g, x)/g(x) for x an element or ideal of R or a map R -> ?
+# - preimage(g, x) for x an element or ideal of Raff
+#
+# In this way we can move everything to quotient rings
+
+################################################################################
+#
+#  Conversions to quotient rings
+#
+################################################################################
+
+function _polyringquo(R::LaurentMPolyWrapRing)
+  get_attribute!(R, :polyring) do
+    n = nvars(R)
+    C = base_ring(R)
+    Cx, x = polynomial_ring(C, append!(["x$i" for i in 1:n], ["x$i^-1" for i in 1:n]))
+    I = ideal(Cx, [x[i]*x[i + n] - 1 for i in 1:n])
+    Q, = quo(Cx, I)
+    return _LaurentMPolyBackend(R, Q)
+  end
+end
+
+function _evaluate_gens_cache(f::_LaurentMPolyBackend{D, C}) where {D, C}
+  if !isdefined(f, :_gens_cache)
+    f._gens_cache = gens(codomain(f))[1:nvars(domain(f))]
+  end
+  return f._gens_cache::Vector{elem_type(C)}
+end
+
+domain(f::_LaurentMPolyBackend) = f.R
+
+codomain(f::_LaurentMPolyBackend) = f.Q
+
+(f::_LaurentMPolyBackend)(p) = image(f, p)
+
+# conversion of elements
+# computes an exponent vector e (possible negative exponents) and a polynomial q
+# such that p = x^e * q
+# note that q is in the underlying multivariate polynomial ring
+# in particular <p> = <q> in the Laurent polynomial ring
+_split(p::LaurentMPolyWrap) = p.mindegs, p.mpoly
+
+function image(f::_LaurentMPolyBackend, p::LaurentMPolyWrap)
+  @assert parent(p) === domain(f)
+  Q = codomain(f)
+  n = nvars(domain(f))
+  _gens = _evaluate_gens_cache(f)
+  # I try to be a bit clever here
+  exps, poly = _split(p)
+  r = evaluate(poly, _gens)
+  v = zeros(Int, ngens(Q))
+  for i in 1:length(exps)
+    if exps[i] >= 0
+      v[i] = exps[i]
+    else
+      v[n + i] = -exps[i]
+    end
+  end
+  m = Q(monomial(base_ring(Q), v))
+  # TODO:
+  # One can remove the evaluate with some more tricks, since this is just
+  # K[x1,...xn] -> K[x1,...,xn,y1,...yn] -> K[x,y]/(xy - 1)
+  return m * evaluate(poly, _gens)
+end
+
+function preimage(f::_LaurentMPolyBackend, q::MPolyQuoRingElem)
+  @assert parent(q) === codomain(f)
+  return f.inv(q)
+end
+
+# conversion of ideals
+function image(f::_LaurentMPolyBackend{D, C, M}, I::LaurentMPolyIdeal) where {D, C, M}
+  @assert base_ring(I) === domain(f)
+  if isdefined(I, :data)
+    return I.data::ideal_type(C)
+  else
+    I.data = ideal(codomain(f), f.(gens(I)))
+    return I.data::ideal_type(C)
+  end
+end
+
+function preimage(f::_LaurentMPolyBackend, J::MPolyQuoIdeal)
+  @assert base_ring(J) === codomain(f)
+  I = ideal(domain(f), map(x -> preimage(f, x), gens(J)))
+  I.data = J
+  return I
+end
+
+# conversion for maps
+function image(f::_LaurentMPolyBackend, g::LaurentMPolyAnyMap)
+  Q = codomain(f)
+  _images = (g.image_of_gens)
+  __images = append!(copy(_images), inv.(_images))
+  return hom(Q, codomain(g), __images, check = false)
+end
+
+################################################################################
+#
+#  Maps from Laurent polynomial rings
+#
+################################################################################
+
+domain(f::LaurentMPolyAnyMap) = f.R
+
+codomain(f::LaurentMPolyAnyMap) = f.S
+
+function hom(R::LaurentMPolyRing, S::Ring, images::Vector; check::Bool = true)
+  @req length(images) == nvars(R) "Wrong number of images"
+  if check
+    @req all(is_unit, images) "Images of generators must be units"
+  end
+  _images = S.(images)
+  return LaurentMPolyAnyMap(R, S, _images)
+end
+
+function image(f::LaurentMPolyAnyMap{D, C}, g::LaurentMPolyRingElem) where {D, C}
+  @req parent(g) === domain(f) "Element not in domain of map"
+  return evaluate(g, f.image_of_gens::Vector{elem_type(C)})
+end
+
+function preimage(f::LaurentMPolyAnyMap, g)
+  @req parent(g) === codomain(f) "Element not in domain of map"
+  q = _polyringquo(domain(f))
+  qf = q(f)
+  return preimage(q, preimage(qf, g))
+end
+
+(f::LaurentMPolyAnyMap)(g::LaurentMPolyRingElem) = image(f, g)
+
+# preimage for ideals
+function preimage(f::MPolyAnyMap{X, <: LaurentMPolyRing, Y, Z}, I::LaurentMPolyIdeal) where {X, Y, Z}
+  R = domain(f)
+  S = codomain(f)
+  if coefficient_ring(R) === base_ring(S) && f.(gens(R)) == gens(S)
+    # We try to do something clever if f : K[x1,...xn] -> K[x1^+-,...xn^+-] is
+    # the natural inclusion
+    polygens = map(x -> _split(x)[2], gens(I))
+    # These are polynomials generating the same ideal as I
+    II = ideal(polygens) # this is an element of the internal polynomial ring underpinning
+                         # the Laurent polynomial ring R
+    _R = base_ring(II)
+    for g in gens(_R)
+      II = saturation(II, g*_R)
+    end
+    # We need to translate to an ideal of R
+    return ideal(R, map(x -> map_coefficients(identity, x, parent = R), gens(II)))
+  end
+  q = _polyringquo(codomain(f))
+  qI = q(I) # ideal in the quotient
+  # map from domain(f) to quotient
+  qf = hom(domain(f), codomain(q), [q(f(g)) for g in gens(domain(f))])
+  return preimage(qf, qI)
+end
+
+################################################################################
+#
+#  Ideals
+#
+################################################################################
+
+base_ring(I::LaurentMPolyIdeal{T}) where {T} = I.R::parent_type(T)
+
+gens(I::LaurentMPolyIdeal) = I.gens
+
+@enable_all_show_via_expressify LaurentMPolyIdeal
+
+function AbstractAlgebra.expressify(a::LaurentMPolyIdeal; context = nothing)
+  return Expr(:call, :ideal, [AbstractAlgebra.expressify(g, context = context) for g in collect(gens(a))]...)
+end
+
+function ideal(R::LaurentMPolyRing, x::Vector)
+  return LaurentMPolyIdeal(R, filter!(!iszero, R.(x)))
+end
+
+function in(x::LaurentMPolyRingElem, I::LaurentMPolyIdeal)
+  R = parent(x)
+  if parent(x) !== base_ring(I)
+    return false
+  end
+  f = _polyringquo(R)
+  return f(x) in f(I)
+end
+
+function +(I::LaurentMPolyIdeal, J::LaurentMPolyIdeal)
+  @req base_ring(I) === base_ring(J) "Rings must be equal"
+  R = base_ring(I)
+  f = _polyringquo(R)
+  IpJ = preimage(f, f(I) + f(J))
+  return IpJ
+end
+
+end # module
+

experimental/Laurent/src/Types.jl

@@ -0,0 +1,39 @@
+import .AbstractAlgebra.Generic: LaurentMPolyWrapRing, LaurentMPolyWrap
+import .AbstractAlgebra: LaurentMPolyRing, LaurentMPolyRingElem
+
+@attributes mutable struct LaurentMPolyAnyMap{D, C} <: Map{D, C, Hecke.Hecke.Map, LaurentMPolyAnyMap}
+  R::D
+  S::C
+  image_of_gens
+  data # map from _polyringquo(R) -> C
+
+  function LaurentMPolyAnyMap(R::D, S::C, image_of_gens) where {D, C}
+    @assert all(x -> parent(x) === S, image_of_gens)
+    return new{D, C}(R, S, image_of_gens)
+  end
+end
+
+mutable struct LaurentMPolyIdeal{T}
+  R
+  gens::Vector{T}
+  data # ideal of of _polyringquo(R)
+
+  function LaurentMPolyIdeal(R::LaurentMPolyRing, gens::Vector)
+    @assert all(x -> parent(x) === R, gens)
+    return new{eltype(gens)}(R, gens)
+  end
+end
+
+mutable struct _LaurentMPolyBackend{D, C, M}
+  R::D
+  Q::C
+  inv::M
+  _gens_cache
+
+  function _LaurentMPolyBackend(R::D, Q::C) where {D, C}
+    _inv = hom(Q, R, vcat(gens(R), inv.(gens(R))))
+    return new{D, C, typeof(_inv)}(R, Q, _inv)
+  end
+end
+
+

experimental/Laurent/test/runtests.jl

@@ -0,0 +1,34 @@
+@testset "Laurent" begin
+  for K in [QQ, GF(5)]
+    Kx, x = LaurentPolynomialRing(K, 2, "x")
+    I = ideal(Kx, [x[1]])
+    @test gens(I) == [x[1]]
+    @test one(Kx) in I
+    @test x[1]^-1 in I
+    @test base_ring(I) === Kx
+    _Kx, _x = LaurentPolynomialRing(K, 3, "xx")
+    @test !(_x[1] in I)
+
+    f = hom(Kx, K, [K(2), K(3)])
+    @test domain(f) === Kx
+    @test codomain(f) === K
+    @test f(x[1]^-1 + x[2]) == K(2)^-1 + K(3)
+
+    @test_throws ArgumentError hom(Kx, ZZ, [ZZ(2), ZZ(3)])
+
+    Ky, y = polynomial_ring(K, 2, "y")
+    f = hom(Ky, Kx, gens(Kx))
+    @test f(y[1]) == x[1]
+    @test preimage(f, I) == ideal(Ky, [one(Ky)])
+
+    f = hom(Ky, Kx, reverse(gens(Kx)))
+    @test f(y[1]) == x[2]
+    @test preimage(f, I) == ideal(Ky, [one(Ky)])
+
+    Q, = quo(Ky, ideal(Ky, [y[1] * y[2] - 1]))
+    f = hom(Kx, Q, gens(Q))
+    for q in gens(Q)
+      @test f(preimage(f, q)) == q
+    end
+  end
+end

src/Rings/MPolyQuo.jl

@@ -564,7 +564,13 @@ end
 
 *(a::MPolyQuoRingElem{S}, b::MPolyQuoRingElem{S}) where {S} = check_parent(a, b) && simplify(MPolyQuoRingElem(a.f*b.f, a.P))
 
-^(a::MPolyQuoRingElem, b::Base.Integer) = simplify(MPolyQuoRingElem(Base.power_by_squaring(a.f, b), a.P))
+function Base.:(^)(a::MPolyQuoRingElem, b::Base.Integer)
+  if b >= 0
+    simplify(MPolyQuoRingElem(Base.power_by_squaring(a.f, b), a.P))
+  else
+    return inv(a)^(-b)
+  end
+end
 
 *(a::MPolyQuoRingElem, b::QQFieldElem) = simplify(MPolyQuoRingElem(a.f * b, a.P))