From abbe2ac324296f1df3efc1562c767937aaedaec0 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Lars=20G=C3=B6ttgens?= Date: Fri, 11 Oct 2024 11:09:48 +0200 Subject: [PATCH] Unsafe funcs for FreeAssAlgebra, Poly, MPoly (#1847) * Better conformance tests for FreeAssAlgebra * Add `zero!`, `one!`, `neg!` for FreeAssAlgElem * Add `neg!`, `one!` for Poly * Add `zero!`, `one!`, `neg!` for MPoly * Add `add!`, `sub!` for FreeAssAlgElem * Fix `neg!` for `UnivPoly` --- src/generic/FreeAssociativeAlgebra.jl | 142 +++++++++++++++++++- src/generic/MPoly.jl | 42 +++++- src/generic/Poly.jl | 23 ++++ src/generic/UnivPoly.jl | 2 +- src/generic/imports.jl | 1 + test/Rings-conformance-tests.jl | 25 ++++ test/generic/FreeAssociativeAlgebra-test.jl | 11 +- test/generic/MPoly-test.jl | 15 +++ 8 files changed, 249 insertions(+), 12 deletions(-) diff --git a/src/generic/FreeAssociativeAlgebra.jl b/src/generic/FreeAssociativeAlgebra.jl index 53a7498cb3..b82756a4a7 100644 --- a/src/generic/FreeAssociativeAlgebra.jl +++ b/src/generic/FreeAssociativeAlgebra.jl @@ -71,7 +71,7 @@ function isone(a::FreeAssociativeAlgebraElem{T}) where T if length(a) < 1 return isone(zero(base_ring(a))) else - return a.length == 1 && isone(a.coeffs[1]) && isempty(a.exps[1]) + return length(a) == 1 && isone(a.coeffs[1]) && isempty(a.exps[1]) end end @@ -646,6 +646,146 @@ function divexact( return combine_like_terms!(FreeAssociativeAlgebraElem{T}(R, zcoeffs, copy(a.exps), n)) end +############################################################################### +# +# Unsafe arithmetic functions +# +############################################################################### + +function zero!(a::FreeAssociativeAlgebraElem{T}) where T <: RingElement + a.length = 0 + return a +end + +function one!(a::FreeAssociativeAlgebraElem{T}) where T <: RingElement + a.length = 1 + fit!(a, 1) + a.coeffs[1] = one(base_ring(parent(a))) + a.exps[1] = Int[] + return a +end + +function neg!(a::FreeAssociativeAlgebraElem{T}) where T <: RingElement + for i in 1:length(a) + a.coeffs[i] = neg!(a.coeffs[i]) + end + return a +end + +function neg!(z::FreeAssociativeAlgebraElem{T}, a::FreeAssociativeAlgebraElem{T}) where T <: RingElement + if z === a + return neg!(a) + end + z.length = length(a) + fit!(z, length(a)) + for i in 1:length(a) + if isassigned(z.coeffs, i) + z.coeffs[i] = neg!(z.coeffs[i], a.coeffs[i]) + else + z.coeffs[i] = -a.coeffs[i] + end + # mutating z.exps[i] is not allowed since it could be aliased + z.exps[i] = a.exps[i] + end + return z +end + +function add!(a::FreeAssociativeAlgebraElem{T}, b::FreeAssociativeAlgebraElem{T}) where T <: RingElement + iszero(b) && return a + return add!(zero(a), a, b) +end + +function add!(z::FreeAssociativeAlgebraElem{T}, a::FreeAssociativeAlgebraElem{T}, b::FreeAssociativeAlgebraElem{T}) where T <: RingElement + if z === a + return add!(z, b) + elseif z === b + return add!(z, a) + end + z.coeffs = empty!(z.coeffs) + z.exps = empty!(z.exps) + i = j = 1 + while i <= a.length && j <= b.length + c = word_cmp(a.exps[i], b.exps[j]) + if c < 0 + push!(z.coeffs, b.coeffs[j]) + push!(z.exps, b.exps[j]) + j += 1 + elseif c > 0 + push!(z.coeffs, a.coeffs[i]) + push!(z.exps, a.exps[i]) + i += 1 + else + s = a.coeffs[i] + b.coeffs[j] + if !iszero(s) + push!(z.coeffs, s) + push!(z.exps, a.exps[i]) + end + i += 1 + j += 1 + end + end + while i <= a.length + push!(z.coeffs, a.coeffs[i]) + push!(z.exps, a.exps[i]) + i += 1 + end + while j <= b.length + push!(z.coeffs, b.coeffs[j]) + push!(z.exps, b.exps[j]) + j += 1 + end + z.length = length(z.coeffs) + return z +end + +function sub!(a::FreeAssociativeAlgebraElem{T}, b::FreeAssociativeAlgebraElem{T}) where T <: RingElement + iszero(b) && return a + return sub!(zero(a), a, b) +end + +function sub!(z::FreeAssociativeAlgebraElem{T}, a::FreeAssociativeAlgebraElem{T}, b::FreeAssociativeAlgebraElem{T}) where T <: RingElement + if z === a + return sub!(z, b) + elseif z === b + return sub!(zero(a), a, b) + end + z.coeffs = empty!(z.coeffs) + z.exps = empty!(z.exps) + i = j = 1 + while i <= a.length && j <= b.length + c = word_cmp(a.exps[i], b.exps[j]) + if c < 0 + push!(z.coeffs, -b.coeffs[j]) + push!(z.exps, b.exps[j]) + j += 1 + elseif c > 0 + push!(z.coeffs, a.coeffs[i]) + push!(z.exps, a.exps[i]) + i += 1 + else + s = a.coeffs[i] - b.coeffs[j] + if !iszero(s) + push!(z.coeffs, s) + push!(z.exps, a.exps[i]) + end + i += 1 + j += 1 + end + end + while i <= a.length + push!(z.coeffs, a.coeffs[i]) + push!(z.exps, a.exps[i]) + i += 1 + end + while j <= b.length + push!(z.coeffs, -b.coeffs[j]) + push!(z.exps, b.exps[j]) + j += 1 + end + z.length = length(z.coeffs) + return z +end + ################################################################################ # diff --git a/src/generic/MPoly.jl b/src/generic/MPoly.jl index 0a4d93e5b6..96861478bd 100644 --- a/src/generic/MPoly.jl +++ b/src/generic/MPoly.jl @@ -3886,6 +3886,43 @@ end # ############################################################################### +function zero!(a::MPoly{T}) where {T <: RingElement} + a.length = 0 + return a +end + +function one!(a::MPoly{T}) where {T <: RingElement} + a.length = 1 + fit!(a, 1) + a.coeffs[1] = one(base_ring(a)) + a.exps = zero(a.exps) + return a +end + +function neg!(a::MPoly{T}) where {T <: RingElement} + for i in 1:length(a) + a.coeffs[i] = neg!(a.coeffs[i]) + end + return a +end + +function neg!(z::MPoly{T}, a::MPoly{T}) where {T <: RingElement} + if z === a + return neg!(a) + end + z.length = length(a) + fit!(z, length(a)) + for i in 1:length(a) + if isassigned(z.coeffs, i) + z.coeffs[i] = neg!(z.coeffs[i], a.coeffs[i]) + else + z.coeffs[i] = -a.coeffs[i] + end + end + z.exps[:,1:length(a)] .= a.exps[:,1:length(a)] + return z +end + function add!(a::MPoly{T}, b::MPoly{T}, c::MPoly{T}) where {T <: RingElement} t = b + c a.coeffs = t.coeffs @@ -3933,11 +3970,6 @@ function fit!(a::MPoly{T}, n::Int) where {T <: RingElement} end return nothing end -# -function zero!(a::MPoly{T}) where {T <: RingElement} - a.length = 0 - return a -end @doc raw""" setcoeff!(a::MPoly{T}, i::Int, c::T) where T <: RingElement diff --git a/src/generic/Poly.jl b/src/generic/Poly.jl index af747c5981..cb75383418 100644 --- a/src/generic/Poly.jl +++ b/src/generic/Poly.jl @@ -198,6 +198,29 @@ function zero!(c::Poly{T}) where T <: RingElement return c end +function one!(c::Poly{T}) where T <: RingElement + fit!(c, 1) + c = set_length!(c, 1) + c.coeffs[1] = one(base_ring(c)) + return c +end + +function neg!(a::Poly{T}) where T <: RingElement + for i in 1:length(a) + a.coeffs[i] = neg!(a.coeffs[i]) + end + return a +end + +function neg!(z::Poly{T}, a::Poly{T}) where T <: RingElement + fit!(z, length(a)) + z = set_length!(z, length(a)) + for i in 1:length(a) + z.coeffs[i] = neg!(z.coeffs[i], a.coeffs[i]) + end + return z +end + function mul!(c::Poly{T}, a::Poly{T}, b::Poly{T}) where T <: RingElement lena = length(a) lenb = length(b) diff --git a/src/generic/UnivPoly.jl b/src/generic/UnivPoly.jl index 34601f8ffe..8f022a7e91 100644 --- a/src/generic/UnivPoly.jl +++ b/src/generic/UnivPoly.jl @@ -1015,7 +1015,7 @@ function one!(a::UnivPoly{T}) where {T <: RingElement} end function neg!(z::UnivPoly{T}, a::UnivPoly{T}) where {T <: RingElement} - z.p = neg!(a.p) + z.p = neg!(z.p, a.p) return z end diff --git a/src/generic/imports.jl b/src/generic/imports.jl index 9b036ce001..2ea33d8985 100644 --- a/src/generic/imports.jl +++ b/src/generic/imports.jl @@ -173,6 +173,7 @@ import ..AbstractAlgebra: mul! import ..AbstractAlgebra: mul_classical import ..AbstractAlgebra: mul_karatsuba import ..AbstractAlgebra: mullow +import ..AbstractAlgebra: neg! import ..AbstractAlgebra: number_of_columns import ..AbstractAlgebra: number_of_generators import ..AbstractAlgebra: number_of_rows diff --git a/test/Rings-conformance-tests.jl b/test/Rings-conformance-tests.jl index e53cfff977..442fbe8394 100644 --- a/test/Rings-conformance-tests.jl +++ b/test/Rings-conformance-tests.jl @@ -42,6 +42,18 @@ function test_elem(Rx::AbstractAlgebra.PolyRing) return Rx(elem_type(R)[test_elem(R) for i in 1:rand(0:6)]) end +function test_elem(Rx::AbstractAlgebra.MPolyRing) + R = base_ring(Rx) + num_gens = ngens(Rx) + iszero(num_gens) && return Rx(test_elem(R)) + len_bound = 8 + exp_bound = rand(1:5) + len = rand(0:len_bound) + coeffs = [test_elem(R) for _ in 1:len] + exps = [[rand(0:exp_bound) for _ in 1:num_gens] for _ in 1:len] + return Rx(coeffs, exps) +end + function test_elem(S::Union{AbstractAlgebra.MatSpace, AbstractAlgebra.MatRing}) R = base_ring(S) @@ -72,6 +84,19 @@ function test_elem(Rx::AbstractAlgebra.SeriesRing) end end +function test_elem(S::AbstractAlgebra.FreeAssociativeAlgebra) + f = S() + g = gens(S) + R = base_ring(S) + isempty(g) && return S(test_elem(R)) + len_bound = 8 + exp_bound = 6 + for i in 1:rand(0:len_bound) + f += test_elem(R) * prod(rand(g) for _ in 1:rand(0:exp_bound); init = S(1)) + end + return f +end + # helper function equality(a::T, b::T) where T <: AbstractAlgebra.NCRingElement diff --git a/test/generic/FreeAssociativeAlgebra-test.jl b/test/generic/FreeAssociativeAlgebra-test.jl index c478a6d0fe..5ee356686d 100644 --- a/test/generic/FreeAssociativeAlgebra-test.jl +++ b/test/generic/FreeAssociativeAlgebra-test.jl @@ -140,10 +140,6 @@ end @test !occursin("\n", sprint(show, R)) end -function test_elem(R::Generic.FreeAssociativeAlgebra{elem_type(ZZ)}) - return rand(R, 0:4, 0:5, -10:10) -end - @testset "Generic.FreeAssociativeAlgebra.change_base_ring" begin F5, = residue_ring(ZZ, 5) R, varsR = polynomial_ring(F5, ["x"]) @@ -221,6 +217,11 @@ end end @testset "Generic.FreeAssociativeAlgebra.NCRing_interface" begin - test_NCRing_interface(free_associative_algebra(ZZ, 3)[1]) + S, = free_associative_algebra(ZZ, 3) + test_NCRing_interface(S) + + R, = QQ[:x, :y] + S, = free_associative_algebra(R, :z => 1:3) + test_NCRing_interface(S) end diff --git a/test/generic/MPoly-test.jl b/test/generic/MPoly-test.jl index 961ba23121..16045227a4 100644 --- a/test/generic/MPoly-test.jl +++ b/test/generic/MPoly-test.jl @@ -1704,3 +1704,18 @@ end @test leading_coefficient(f, 1) == coefficients(f, 1)[end] @test content(f, 1) == y end + +@testset "Generic.MPoly.Ring_interface" begin + S, = polynomial_ring(QQ, 0) + test_Ring_interface_recursive(S) + + S, = polynomial_ring(QQ, 1) + test_Ring_interface_recursive(S) + + S, = polynomial_ring(ZZ, 2) + test_Ring_interface_recursive(S) + + R, = QQ[:x] + S, = polynomial_ring(R, :z => 1:3) + test_Ring_interface(S) # _recursive needs too many ressources +end