Skip to content

Commit

Permalink
removes more singular_assure calls (#2817)
Browse files Browse the repository at this point in the history
  • Loading branch information
@ederc @fieker
ederc authored and fieker committed Sep 29, 2023
1 parent 9ccf6a9 commit 2915158
Show file tree
Hide file tree
Showing 3 changed files with 48 additions and 65 deletions.
10 changes: 10 additions & 0 deletions src/Rings/groebner.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -54,6 +54,16 @@ function groebner_assure(I::MPolyIdeal, ordering::MonomialOrdering, complete_red
end
end

function singular_groebner_generators(I::MPolyIdeal, ordering::MonomialOrdering = default_ordering(base_ring(I)), complete_reduction::Bool = false)
standard_basis(I, ordering=ordering, complete_reduction = complete_reduction)
return singular_generators(I.gb[ordering], ordering)
end

function singular_groebner_generators(I::MPolyIdeal, complete_reduction::Bool, need_global::Bool)
G = groebner_assure(I, complete_reduction, need_global)
return singular_generators(G, G.ord)
end

@doc raw"""
_compute_standard_basis(B::IdealGens; ordering::MonomialOrdering,
complete_reduction::Bool = false)
Expand Down
99 changes: 34 additions & 65 deletions src/Rings/mpoly-ideals.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -94,8 +94,7 @@ ideal(x^3, x^2*y, x*y^2, y^3)
```
"""
function Base.:^(I::MPolyIdeal, m::Int)
singular_assure(I)
return MPolyIdeal(base_ring(I), I.gens.S^m)
return MPolyIdeal(base_ring(I), (singular_generators(I))^m)
end

@doc raw"""
Expand Down Expand Up @@ -189,12 +188,8 @@ ideal(x^3*y - x*y - y^3, x^4 - x^2 - x*y^2)
```
"""
function Base.intersect(I::MPolyIdeal{T}, Js::MPolyIdeal{T}...) where T
singular_assure(I)
si = I.gens.S
for J in Js
singular_assure(J)
end
si = Singular.intersection(si, [J.gens.S for J in Js]...)
si = singular_generators(I)
si = Singular.intersection(si, [singular_generators(J) for J in Js]...)
return MPolyIdeal(base_ring(I), si)
end

Expand Down Expand Up @@ -240,9 +235,7 @@ ideal(x^2*y + y^3, x^3*z + x*y^2*z + y^4, x^2*z^2 + x*y^3 - x*y^2*z + y^2*z^2, x
"""
function quotient(I::MPolyIdeal{T}, J::MPolyIdeal{T}) where T
@assert base_ring(I) == base_ring(J)
singular_assure(I)
singular_assure(J)
return MPolyIdeal(base_ring(I), Singular.quotient(I.gens.S, J.gens.S))
return MPolyIdeal(base_ring(I), Singular.quotient(singular_generators(I), singular_generators(J)))
end

function quotient(I::MPolyIdeal{T}, f::T) where T
Expand Down Expand Up @@ -277,18 +270,14 @@ ideal(z, x*y)
```
"""
function saturation(I::MPolyIdeal{T}, J::MPolyIdeal{T}) where T
singular_assure(I)
singular_assure(J)
K, _ = Singular.saturation(I.gens.S, J.gens.S)
K, _ = Singular.saturation(singular_generators(I), singular_generators(J))
return MPolyIdeal(base_ring(I), K)
end

# the following is corresponding to saturation2 from Singular
# TODO: think about how to use use this properly/automatically
function _saturation2(I::MPolyIdeal{T}, J::MPolyIdeal{T}) where T
singular_assure(I)
singular_assure(J)
K, _ = Singular.saturation2(I.gens.S, J.gens.S)
K, _ = Singular.saturation2(singular_generators(I), singular_generators(J))
return MPolyIdeal(base_ring(I), K)
end

Expand All @@ -314,9 +303,7 @@ julia> K, m = saturation_with_index(I, J)
```
"""
function saturation_with_index(I::MPolyIdeal{T}, J::MPolyIdeal{T}) where T
singular_assure(I)
singular_assure(J)
K, k = Singular.saturation(I.gens.S, J.gens.S)
K, k = Singular.saturation(singular_generators(I), singular_generators(J))
return (MPolyIdeal(base_ring(I), K), k)
end

Expand Down Expand Up @@ -363,9 +350,8 @@ Multivariate polynomial ring in 4 variables t, x, y, z
```
"""
function eliminate(I::MPolyIdeal{T}, l::Vector{T}) where T <: MPolyRingElem
singular_assure(I)
S = base_ring(I.gens.S)
s = Singular.eliminate(I.gens.S, [S(x) for x = l]...)
S = singular_polynomial_ring(I)
s = Singular.eliminate(singular_generators(I), [S(x) for x = l]...)
return MPolyIdeal(base_ring(I), s)
end
function eliminate(I::MPolyIdeal, l::AbstractVector{Int})
Expand Down Expand Up @@ -433,7 +419,6 @@ ideal(102*b*d, 78*a*d, 51*b*c, 39*a*c, 6*a*b*d, 3*a*b*c)
```
"""
@attr T function radical(I::T) where {T <: MPolyIdeal}
singular_assure(I)
R = base_ring(I)
if isa(base_ring(R), NumField) && !isa(base_ring(R), AnticNumberField)
A, mA = absolute_simple_field(base_ring(R))
Expand All @@ -442,9 +427,9 @@ ideal(102*b*d, 78*a*d, 51*b*c, 39*a*c, 6*a*b*d, 3*a*b*c)
set_attribute!(Irad, :is_radical => true)
return Irad
elseif elem_type(base_ring(R)) <: FieldElement
J = Singular.LibPrimdec.radical(I.gens.Sx, I.gens.S)
elseif base_ring(I.gens.Sx) isa Singular.Integers
J = Singular.LibPrimdecint.radicalZ(I.gens.Sx, I.gens.S)
J = Singular.LibPrimdec.radical(singular_polynomial_ring(I), singular_generators(I))
elseif base_ring(singular_polynomial_ring(I)) isa Singular.Integers
J = Singular.LibPrimdecint.radicalZ(singular_polynomial_ring(I), singular_generators(I))
else
error("not implemented for base ring")
end
Expand Down Expand Up @@ -556,17 +541,16 @@ function _compute_primary_decomposition(I::MPolyIdeal; algorithm::Symbol=:GTZ)
end
return Tuple{typeof(I), typeof(I)}[(map_coefficients(mA, x[1], parent = R), map_coefficients(mA, x[2], parent = R)) for x = pd]
end
singular_assure(I)
if elem_type(base_ring(R)) <: FieldElement
if algorithm == :GTZ
L = Singular.LibPrimdec.primdecGTZ(I.gens.Sx, I.gens.S)
L = Singular.LibPrimdec.primdecGTZ(singular_polynomial_ring(I), singular_generators(I))
elseif algorithm == :SY
L = Singular.LibPrimdec.primdecSY(I.gens.Sx, I.gens.S)
L = Singular.LibPrimdec.primdecSY(singular_polynomial_ring(I), singular_generators(I))
else
error("algorithm invalid")
end
elseif base_ring(I.gens.Sx) isa Singular.Integers
L = Singular.LibPrimdecint.primdecZ(I.gens.Sx, I.gens.S)
elseif base_ring(singular_polynomial_ring(I)) isa Singular.Integers
L = Singular.LibPrimdecint.primdecZ(singular_polynomial_ring(I), singular_generators(I))
else
error("base ring not implemented")
end
Expand Down Expand Up @@ -661,8 +645,7 @@ Multivariate polynomial ring in 2 variables over number field graded by
if is_zero(I)
return [(ideal(R, zero(R)), ideal(R, zero(R)), ideal(R, zero(R)), 1)]
end
singular_assure(I)
(S, d) = Singular.LibPrimdec.absPrimdecGTZ(I.gens.Sx, I.gens.S)
(S, d) = Singular.LibPrimdec.absPrimdecGTZ(singular_polynomial_ring(I), singular_generators(I))
decomp = d[:primary_decomp]
absprimes = d[:absolute_primes]
@assert length(decomp) == length(absprimes)
Expand Down Expand Up @@ -773,17 +756,16 @@ function minimal_primes(I::MPolyIdeal; algorithm::Symbol = :GTZ)
mp = minimal_primes(map_coefficients(pseudo_inv(mA), I); algorithm = algorithm)
return typeof(I)[map_coefficients(mA, x) for x = mp]
end
singular_assure(I)
if elem_type(base_ring(R)) <: FieldElement
if algorithm == :GTZ
l = Singular.LibPrimdec.minAssGTZ(I.gens.Sx, I.gens.S)
l = Singular.LibPrimdec.minAssGTZ(singular_polynomial_ring(I), singular_generators(I))
elseif algorithm == :charSets
l = Singular.LibPrimdec.minAssChar(I.gens.Sx, I.gens.S)
l = Singular.LibPrimdec.minAssChar(singular_polynomial_ring(I), singular_generators(I))
else
error("algorithm invalid")
end
elseif base_ring(I.gens.Sx) isa Singular.Integers
l = Singular.LibPrimdecint.minAssZ(I.gens.Sx, I.gens.S)
elseif base_ring(singular_polynomial_ring(I)) isa Singular.Integers
l = Singular.LibPrimdecint.minAssZ(singular_polynomial_ring(I), singular_generators(I))
else
error("base ring not implemented")
end
Expand Down Expand Up @@ -833,8 +815,7 @@ julia> L = equidimensional_decomposition_weak(I)
eq = equidimensional_decomposition_weak(map_coefficients(pseudo_inv(mA), I))
return typeof(I)[map_coefficients(mA, x, parent = R) for x = eq]
end
singular_assure(I)
l = Singular.LibPrimdec.equidim(I.gens.Sx, I.gens.S)
l = Singular.LibPrimdec.equidim(singular_polynomial_ring(I), singular_generators(I))
V = [ideal(R, i) for i in l]
if length(V) == 1 && is_one(gen(V[1], 1))
return typeof(I)[]
Expand Down Expand Up @@ -878,8 +859,7 @@ julia> L = equidimensional_decomposition_radical(I)
eq = equidimensional_decomposition_radical(map_coefficients(pseudo_inv(mA), I))
return typeof(I)[map_coefficients(mA, x) for x = eq]
end
singular_assure(I)
l = Singular.LibPrimdec.prepareAss(I.gens.Sx, I.gens.S)
l = Singular.LibPrimdec.prepareAss(singular_polynomial_ring(I), singular_generators(I))
V = [ideal(R, i) for i in l]
if length(V) == 1 && is_one(gen(V[1], 1))
return typeof(I)[]
Expand Down Expand Up @@ -939,11 +919,10 @@ function equidimensional_hull(I::MPolyIdeal)
eq = equidimensional_hull(map_coefficients(pseudo_inv(mA), I))
return map_coefficients(mA, eq)
end
singular_assure(I)
if elem_type(base_ring(R)) <: FieldElement
i = Singular.LibPrimdec.equidimMax(I.gens.Sx, I.gens.S)
elseif base_ring(I.gens.Sx) isa Singular.Integers
i = Singular.LibPrimdecint.equidimZ(I.gens.Sx, I.gens.S)
i = Singular.LibPrimdec.equidimMax(singular_polynomial_ring(I), singular_generators(I))
elseif base_ring(singular_polynomial_ring(I)) isa Singular.Integers
i = Singular.LibPrimdecint.equidimZ(singular_polynomial_ring(I), singular_generators(I))
else
error("base ring not implemented")
end
Expand Down Expand Up @@ -984,8 +963,7 @@ function equidimensional_hull_radical(I::MPolyIdeal)
eq = equidimensional_hull_radical(map_coefficients(pseudo_inv(mA), I))
return map_coefficients(mA, eq)
end
singular_assure(I)
i = Singular.LibPrimdec.equiRadical(I.gens.Sx, I.gens.S)
i = Singular.LibPrimdec.equiRadical(singular_polynomial_ring(I), singular_generators(I))
return ideal(R, i)
end

Expand Down Expand Up @@ -1040,10 +1018,7 @@ true
```
"""
function is_subset(I::MPolyIdeal{T}, J::MPolyIdeal{T}) where T
singular_assure(I)
G = groebner_assure(J)
singular_assure(G)
return Singular.iszero(Singular.reduce(I.gens.S, G.S))
return Singular.iszero(Singular.reduce(singular_generators(I), singular_groebner_generators(J)))
end

### todo: wenn schon GB's bekannt ...
Expand Down Expand Up @@ -1077,10 +1052,8 @@ false
```
"""
function ideal_membership(f::T, I::MPolyIdeal{T}; ordering::MonomialOrdering = default_ordering(base_ring(I))) where T
GI = standard_basis(I, ordering=ordering, complete_reduction=false)
singular_assure(GI)
Sx = base_ring(GI.S)
return Singular.iszero(Singular.reduce(Sx(f), GI.S))
Sx = singular_polynomial_ring(I, ordering)
return Singular.iszero(Singular.reduce(Sx(f), singular_groebner_generators(I, ordering)))
end
Base.:in(f::MPolyRingElem, I::MPolyIdeal) = ideal_membership(f,I)
#######################################################
Expand Down Expand Up @@ -1112,9 +1085,8 @@ false
```
"""
function radical_membership(f::T, I::MPolyIdeal{T}) where T
singular_assure(I)
Sx = base_ring(I.gens.S)
return Singular.LibPolylib.rad_con(Sx(f), I.gens.S) == 1
Sx = singular_polynomial_ring(I)
return Singular.LibPolylib.rad_con(Sx(f), singular_generators(I)) == 1
end
inradical(f::MPolyRingElem, I::MPolyIdeal) = radical_membership(f,I)
################################################################################
Expand Down Expand Up @@ -1286,9 +1258,7 @@ julia> dim(I)
if I.dim > -1
return I.dim
end
G = groebner_assure(I, false, true)
singular_assure(G)
I.dim = Singular.dimension(G.S)
I.dim = Singular.dimension(singular_groebner_generators(I, false, true))
return I.dim
end

Expand Down Expand Up @@ -1491,8 +1461,7 @@ function small_generating_set(I::MPolyIdeal)
@req coefficient_ring(R) isa Field "The coefficient ring must be a field"

# in the ungraded case, mstd's heuristic returns smaller gens when recomputing gb
singular_assure(I)
sing_gb, sing_min = Singular.mstd(I.gens.gens.S)
sing_gb, sing_min = Singular.mstd(singular_generators(I))
ring = I.gens.Ox
computed_gb = IdealGens(ring, sing_gb, true)
if !haskey(I.gb,computed_gb.ord)
Expand Down
4 changes: 4 additions & 0 deletions src/Rings/mpoly.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -709,6 +709,10 @@ function singular_generators(I::MPolyIdeal, monorder::MonomialOrdering=default_o
return singular_generators(generating_system(I), monorder)
end

function singular_polynomial_ring(I::MPolyIdeal, monorder::MonomialOrdering=default_ordering(base_ring(I)))
return (singular_generators(I, monorder)).base_ring
end

function singular_assure(I::BiPolyArray)
if !isdefined(I, :S)
I.Sx = singular_poly_ring(I.Ox)
Expand Down

0 comments on commit 2915158

Please sign in to comment.