From b9298e553c31f2540a02f8d33012b38d3bbf936c Mon Sep 17 00:00:00 2001 From: HechtiDerLachs Date: Thu, 19 Sep 2024 11:18:06 +0200 Subject: [PATCH 1/4] WIP on using sparse matrices for orderings. --- src/Rings/orderings.jl | 46 +++++++++++++++++++++++++++--------------- 1 file changed, 30 insertions(+), 16 deletions(-) diff --git a/src/Rings/orderings.jl b/src/Rings/orderings.jl index 42fecc6e41ad..81ca9447aece 100644 --- a/src/Rings/orderings.jl +++ b/src/Rings/orderings.jl @@ -71,7 +71,7 @@ end # if fullrank is true, then the matrix should be invertible struct MatrixOrdering <: AbsGenOrdering vars::Vector{Int} - matrix::ZZMatrix + matrix::SMat{ZZRingElem} fullrank::Bool function MatrixOrdering(v, m::ZZMatrix, fullrank::Bool) @req length(v) == ncols(m) "number of variables should match the number of columns" @@ -87,29 +87,43 @@ struct MatrixOrdering <: AbsGenOrdering end end -function _canonical_matrix(w) - ww = matrix(ZZ, 0, ncols(w), []) +# TODO: Move the following two methods to Hecke +function content(v::SRow{ZZRingElem}) + return reduce(gcd, v.values) # TODO: Make this clean! +end + +function _canonical_matrix(w_in::ZZMatrix) + # The matrix coming in from polynomial rings is not sparse. + # We first convert it, but eventually we should consider using sparse matrices also there. + w = sparse_matrix(w_in) + ww = sparse_matrix(ZZ, 0, ncols(w)) for i in 1:nrows(w) - if is_zero_row(w, i) - continue - end - nw = w[i:i, :] + is_zero(w[i]) && continue + # if is_zero_row(w, i) + # continue + # end + # nw = w[i:i, :] + nw = w[i] c = content(nw) if !isone(c) nw = divexact(nw, c) end for j in 1:nrows(ww) - h = findfirst(x->ww[j, x] != 0, 1:ncols(w)) - if !iszero(nw[1, h]) - nw = abs(ww[j, h])*nw - sign(ww[j, h])*nw[1, h]*ww[j:j, :] - end + is_zero(ww[j]) && continue + (h, ent) = first(ww[j]) + nw = abs(ww[j, h])*nw - sign(ww[j, h])*ent*ww[j] + #h = findfirst(x->ww[j, x] != 0, 1:ncols(w)) + #if !iszero(nw[1, h]) + # nw = abs(ww[j, h])*nw - sign(ww[j, h])*nw[1, h]*ww[j:j, :] + #end end if !iszero(nw) c = content(nw) if !isone(c) nw = divexact(nw, c) end - ww = vcat(ww, nw) + push!(ww, nw) + #ww = vcat(ww, nw) end end @assert nrows(ww) <= ncols(ww) @@ -175,7 +189,7 @@ mutable struct MonomialOrdering{S} o::AbsGenOrdering is_total::Bool is_total_is_known::Bool - canonical_matrix::ZZMatrix + canonical_matrix::SMat{ZZRingElem} function MonomialOrdering(R::S, o::AbsGenOrdering) where {S} return new{S}(R, o, false, false) @@ -1367,7 +1381,7 @@ function Hecke.simplify(M::MonomialOrdering) end function canonical_matrix(nvars::Int, M::AbsOrdering) - return _canonical_matrix(_matrix(nvars, M)) + return _canonical_matrix(_matrix(nvars, M))::sparse_matrix_type(ZZRingElem) end @doc raw""" @@ -1809,7 +1823,7 @@ end function _canonical_matrix_intern(o::ModuleOrdering) nvrs = ngens(o.M) + ngens(base_ring(o.M)) eo = _embedded_ring_ordering(o) - return canonical_matrix(nvrs, eo) + return canonical_matrix(nvrs, eo)::sparse_matrix_type(ZZRingElem) end function canonical_matrix(o::ModuleOrdering) @@ -1820,7 +1834,7 @@ function canonical_matrix(o::ModuleOrdering) o.canonical_matrix = _canonical_matrix_intern(o) end end - return o.canonical_matrix + return o.canonical_matrix::sparse_matrix_type(ZZRingElem) end # is_obviously_equal checks whether the two orderings are defined in the exact From af7ed5b7fe07e267f4794fe1101868bb14e0c071 Mon Sep 17 00:00:00 2001 From: Johannes Schmitt Date: Thu, 19 Sep 2024 14:33:09 +0200 Subject: [PATCH 2/4] Make it work --- src/Rings/orderings.jl | 26 +++++++++++--------------- 1 file changed, 11 insertions(+), 15 deletions(-) diff --git a/src/Rings/orderings.jl b/src/Rings/orderings.jl index 81ca9447aece..e387e32e7cfd 100644 --- a/src/Rings/orderings.jl +++ b/src/Rings/orderings.jl @@ -88,7 +88,11 @@ struct MatrixOrdering <: AbsGenOrdering end # TODO: Move the following two methods to Hecke -function content(v::SRow{ZZRingElem}) +function Hecke.content(v::SRow{ZZRingElem}) + isempty(v) && return ZZ(0) + if length(v) == 1 + return abs(v.values[1]) + end return reduce(gcd, v.values) # TODO: Make this clean! end @@ -98,39 +102,31 @@ function _canonical_matrix(w_in::ZZMatrix) w = sparse_matrix(w_in) ww = sparse_matrix(ZZ, 0, ncols(w)) for i in 1:nrows(w) - is_zero(w[i]) && continue - # if is_zero_row(w, i) - # continue - # end - # nw = w[i:i, :] + is_empty(w[i]) && continue nw = w[i] c = content(nw) if !isone(c) nw = divexact(nw, c) end for j in 1:nrows(ww) - is_zero(ww[j]) && continue + is_empty(ww[j]) && continue (h, ent) = first(ww[j]) - nw = abs(ww[j, h])*nw - sign(ww[j, h])*ent*ww[j] - #h = findfirst(x->ww[j, x] != 0, 1:ncols(w)) - #if !iszero(nw[1, h]) - # nw = abs(ww[j, h])*nw - sign(ww[j, h])*nw[1, h]*ww[j:j, :] - #end + hnw = findfirst(isequal(h), nw.pos) + hnw == nothing && continue + nw = abs(ent)*nw - sign(ent)*nw.values[hnw]*ww[j] end - if !iszero(nw) + if !is_empty(nw) c = content(nw) if !isone(c) nw = divexact(nw, c) end push!(ww, nw) - #ww = vcat(ww, nw) end end @assert nrows(ww) <= ncols(ww) return ww end - # convert (y,x,z) => (2,1,3) and check uniqueness function _unique_var_indices(a::AbstractVector{<:MPolyRingElem}) !isempty(a) || error("need at least one variable") From dab90171b93b29e1959fe37433511296ea127910 Mon Sep 17 00:00:00 2001 From: Johannes Schmitt Date: Thu, 19 Sep 2024 15:06:02 +0200 Subject: [PATCH 3/4] Make it faster --- src/Rings/orderings.jl | 11 +++++++---- 1 file changed, 7 insertions(+), 4 deletions(-) diff --git a/src/Rings/orderings.jl b/src/Rings/orderings.jl index e387e32e7cfd..a0b38afd7872 100644 --- a/src/Rings/orderings.jl +++ b/src/Rings/orderings.jl @@ -101,21 +101,24 @@ function _canonical_matrix(w_in::ZZMatrix) # We first convert it, but eventually we should consider using sparse matrices also there. w = sparse_matrix(w_in) ww = sparse_matrix(ZZ, 0, ncols(w)) + ent = ZZ() for i in 1:nrows(w) - is_empty(w[i]) && continue + iszero(w[i]) && continue nw = w[i] c = content(nw) if !isone(c) nw = divexact(nw, c) end for j in 1:nrows(ww) - is_empty(ww[j]) && continue - (h, ent) = first(ww[j]) + iszero(ww[j]) && continue + h = ww[j].pos[1] + # TODO: Provide an "in place" getindex for ZZRingElem_Array + ccall((:fmpz_set, Nemo.libflint), Nothing, (Ref{ZZRingElem}, Ptr{Int}), ent, Hecke.get_ptr(ww[j].values, 1)) hnw = findfirst(isequal(h), nw.pos) hnw == nothing && continue nw = abs(ent)*nw - sign(ent)*nw.values[hnw]*ww[j] end - if !is_empty(nw) + if !iszero(nw) c = content(nw) if !isone(c) nw = divexact(nw, c) From 0002d7878ab9c38bd5311739ddfb8e4ad897d18e Mon Sep 17 00:00:00 2001 From: Johannes Schmitt Date: Thu, 19 Sep 2024 15:09:53 +0200 Subject: [PATCH 4/4] Adjust type --- src/Rings/orderings.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/Rings/orderings.jl b/src/Rings/orderings.jl index a0b38afd7872..63699798ba04 100644 --- a/src/Rings/orderings.jl +++ b/src/Rings/orderings.jl @@ -71,7 +71,7 @@ end # if fullrank is true, then the matrix should be invertible struct MatrixOrdering <: AbsGenOrdering vars::Vector{Int} - matrix::SMat{ZZRingElem} + matrix::ZZMatrix fullrank::Bool function MatrixOrdering(v, m::ZZMatrix, fullrank::Bool) @req length(v) == ncols(m) "number of variables should match the number of columns" @@ -1770,7 +1770,7 @@ mutable struct ModuleOrdering{S} M::S o::AbsModOrdering # must allow gen*mon or mon*gen product ordering - canonical_matrix::ZZMatrix + canonical_matrix::SMat{ZZRingElem, Hecke.ZZRingElem_Array_Mod.ZZRingElem_Array} function ModuleOrdering(M::S, o::AbsModOrdering) where {S} return new{S}(M, o)