From 2bca8fc10b57949357d2d2c4cf2fb271155bad89 Mon Sep 17 00:00:00 2001 From: Tommy Hofmann Date: Wed, 14 Jun 2023 08:39:10 +0200 Subject: [PATCH] First stab at more functionality for Laurent polynomial rings (#2448) --- Project.toml | 2 +- experimental/Laurent/src/Laurent.jl | 213 ++++++++++++++++++++++++++ experimental/Laurent/src/Types.jl | 39 +++++ experimental/Laurent/test/runtests.jl | 34 ++++ src/Rings/MPolyQuo.jl | 8 +- 5 files changed, 294 insertions(+), 2 deletions(-) create mode 100644 experimental/Laurent/src/Laurent.jl create mode 100644 experimental/Laurent/src/Types.jl create mode 100644 experimental/Laurent/test/runtests.jl diff --git a/Project.toml b/Project.toml index d2261e5c46b0..31db1716ac74 100644 --- a/Project.toml +++ b/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" diff --git a/experimental/Laurent/src/Laurent.jl b/experimental/Laurent/src/Laurent.jl new file mode 100644 index 000000000000..c94e4a350aad --- /dev/null +++ b/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

= 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 + diff --git a/experimental/Laurent/src/Types.jl b/experimental/Laurent/src/Types.jl new file mode 100644 index 000000000000..55e611baa062 --- /dev/null +++ b/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 + + diff --git a/experimental/Laurent/test/runtests.jl b/experimental/Laurent/test/runtests.jl new file mode 100644 index 000000000000..bdbb275a7f49 --- /dev/null +++ b/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 diff --git a/src/Rings/MPolyQuo.jl b/src/Rings/MPolyQuo.jl index fbb10d43e752..2f09d3838a69 100644 --- a/src/Rings/MPolyQuo.jl +++ b/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))