From b5f7a7bfa7cac4e0c540e3d5156774a0b609a4a0 Mon Sep 17 00:00:00 2001
From: Wolfram Decker <decker@mathematik.uni-kl.de>
Date: Sun, 27 Oct 2024 12:11:22 +0100
Subject: [PATCH 1/4] Saturation for modules

---
 .../subquotients.md                           |  11 +-
 .../src/CommutativeAlgebra/affine_algebras.md |  29 +-
 src/Modules/UngradedModules/SubquoModule.jl   | 283 +++++++++++++++++-
 3 files changed, 305 insertions(+), 18 deletions(-)

diff --git a/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md b/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md
index 695bf907fe8..5f2c5c7de34 100644
--- a/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md
+++ b/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md
@@ -363,9 +363,18 @@ quotient(M::SubquoModule{T}, N::SubquoModule{T}) where T
 ```
 
 ```@docs
-quotient(M::SubquoModule{T}, J::MPolyIdeal{T}) where T
+quotient(M:: SubquoModule, J::Ideal)
 ```
 
+```@docs
+saturation(M:: SubquoModule, J::Ideal)
+```
+
+```@docs
+saturation_with_index(M:: SubquoModule, J::Ideal)
+```
+
+
 ## Submodules and Quotients
 
 ```@docs
diff --git a/docs/src/CommutativeAlgebra/affine_algebras.md b/docs/src/CommutativeAlgebra/affine_algebras.md
index 9e65ec7002b..3d2a110c1b8 100644
--- a/docs/src/CommutativeAlgebra/affine_algebras.md
+++ b/docs/src/CommutativeAlgebra/affine_algebras.md
@@ -294,11 +294,38 @@ minimal_generating_set(I::MPolyQuoIdeal{<:MPolyDecRingElem})
 
 #### Intersection of Ideals
 
-```@docs
+```@julia
 intersect(a::MPolyQuoIdeal{T}, bs::MPolyQuoIdeal{T}...) where T
 intersect(V::Vector{MPolyQuoIdeal{T}}) where T
 ```
 
+Return the intersection of two or more ideals.
+
+###### Examples
+
+```jldoctest
+julia> R, (x, y) = polynomial_ring(QQ, [:x, :y]);
+
+julia> A, _ = quo(R, ideal(R, [x^2-y^3, x-y]));
+
+julia> a = ideal(A, [y^2])
+Ideal generated by
+  y^2
+
+julia> b = ideal(A, [x])
+Ideal generated by
+  x
+
+julia> intersect(a,b)
+Ideal generated by
+  x*y
+
+julia> intersect([a,b])
+Ideal generated by
+x*y
+
+```
+
 #### Ideal Quotients
 
 ```@docs
diff --git a/src/Modules/UngradedModules/SubquoModule.jl b/src/Modules/UngradedModules/SubquoModule.jl
index e44268dc784..1721a1091e1 100644
--- a/src/Modules/UngradedModules/SubquoModule.jl
+++ b/src/Modules/UngradedModules/SubquoModule.jl
@@ -1771,7 +1771,7 @@ Ideal generated by
 function quotient(M::SubquoModule{T}, N::SubquoModule{T}) where T
   @assert base_ring(M) == base_ring(N)
   @assert ambient_module(M) == ambient_module(N)
-  MplusN, iM, _ = sum(M, N)
+  MplusN = M+N
   Q, _ = quo(MplusN, [iM(x) for x in gens(M)])
   return annihilator(Q)
 end
@@ -1779,8 +1779,7 @@ end
 (::Colon)(M::SubquoModule{T}, N::SubquoModule{T}) where T = quotient(M, N)
 
 @doc raw"""
-     quotient(M::SubquoModule{T}, J::MPolyIdeal{T}) where T
-
+     quotient(M::SubquoModule, J::Ideal)
 
 Return the quotient of `M` by `J`.
 
@@ -1819,12 +1818,10 @@ Ideal generated by
   z^2
 
 julia> L = quotient(M, J)
-Subquotient of Submodule with 5 generators
+Subquotient of Submodule with 3 generators
 1 -> x*e[1]
-2 -> y^3*e[1]
-3 -> z^4*e[1]
-4 -> y*z^3*e[1]
-5 -> y^2*z^2*e[1]
+2 -> y*z^3*e[1]
+3 -> y^2*z^2*e[1]
 by Submodule with 3 generators
 1 -> x^2*e[1]
 2 -> y^3*e[1]
@@ -1833,32 +1830,286 @@ by Submodule with 3 generators
 julia> ambient_free_module(L) == ambient_free_module(M)
 true
 
+```
+
+```jldoctest
+julia> S, (x, y, z) = polynomial_ring(QQ, [:x, :y, :z]);
+
+julia> R, _ = quo(S, ideal(S, [x+y+z]));
+
+julia> F = free_module(R, 1);
+
+julia> AM = R[x;];
+
+julia> BM = R[x^2; y^3; z^4];
+
+julia> M = SubquoModule(F, AM, BM)
+Subquotient of Submodule with 1 generator
+1 -> (-y - z)*e[1]
+by Submodule with 3 generators
+1 -> (y^2 + 2*y*z + z^2)*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+julia> J = ideal(R, [x, y, z])^2
+Ideal generated by
+  x^2
+  x*y
+  x*z
+  y^2
+  y*z
+  z^2
+
+julia> quotient(M, J)
+Subquotient of Submodule with 2 generators
+1 -> z*e[1]
+2 -> y*e[1]
+by Submodule with 3 generators
+1 -> (y^2 + 2*y*z + z^2)*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
 ```
 """
-function quotient(M::SubquoModule{T}, J::MPolyIdeal{T}) where T
+function quotient(M::SubquoModule, J::Ideal)
+  @assert base_ring(M) == base_ring(J)
+  M_quo = isdefined(M, :quo) ? M.quo : SubModuleOfFreeModule(ambient_free_module(M), Vector{elem_type(ambient_free_module(M))}())
+  U = M.sub+M.quo
+  UF = _quotient(U, J)
+  res = SubquoModule(UF)
+  res.quo = M.quo
+  return simplify_light(res)[1]
+end
+
+(::Colon)(M::SubquoModule, J::Ideal) = quotient(M, N)
+
+function _quotient(U::SubModuleOfFreeModule, J::Ideal) ### TODO Replace by generic method
   error("not implemented for the given types of modules.")
 end
 
-(::Colon)(M::SubquoModule{T}, J::MPolyIdeal{T}) where T = quotient(M, N)
+function _quotient(U::SubModuleOfFreeModule, J::Ideal{T}) where T <: Union{MPolyRingElem, MPolyQuoRingElem}
+  F = ambient_free_module(U)
+  SgU = singular_generators(U.gens)
+  SgJ = singular_generators(J.gens)
+  SQ = Singular.quotient(SgU, SgJ)
+  MG = ModuleGens(F, SQ)
+  return SubModuleOfFreeModule(F, MG)
+end
 
-function quotient(M::SubquoModule{T}, J::MPolyIdeal{T}) where T <: MPolyRingElem
+########################################
+### saturation for modules
+########################################
+
+@doc raw"""
+     saturation(M::SubquoModule,
+               J::Ideal = ideal(base_ring(M), gens(base_ring(M)));
+               iteration::Bool = false)
+
+Return the saturation $M:J^{\infty}$ of `M` with respect to `J`.
+
+If the ideal `J` is not given, the ideal generated by the generators (variables) of `base_ring(M)` is used.
+
+Setting `iteration` to `true` only has an effect over rings of type `MPolyRing` or `MPolyQuoRing`. Over such rings,
+if `iteration` is set to `true`, the saturation is done by carrying out successive ideal quotient computations as
+suggested by the definition of saturation. Otherwise, a more sophisticated Gröbner basis approach is used which is typically
+faster. Applying the two approaches may lead to different generating sets of the saturation.
+
+!!! note
+    By definition,
+    $M:J^{\infty} = \{ a \in A \mid J^ka \subset M \text{ for some } k\geq 1 \} = \bigcup\limits_{k=1}^{\infty} (M:J^k).$
+    Here, $A$ is the ambient module of $M$.
+
+# Examples
+```jldoctest
+julia> R, (x, y, z) = polynomial_ring(QQ, [:x, :y, :z]);
+
+julia> F = free_module(R, 1);
+
+julia> AM = R[x;];
+
+julia> BM = R[x^2; y^3; z^4];
+
+julia> M = SubquoModule(F, AM, BM)
+Subquotient of Submodule with 1 generator
+1 -> x*e[1]
+by Submodule with 3 generators
+1 -> x^2*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+julia> J = ideal(R, [x, y, z])^2
+Ideal generated by
+  x^2
+  x*y
+  x*z
+  y^2
+  y*z
+  z^2
+
+julia> saturation(M, J)
+Subquotient of Submodule with 1 generator
+1 -> e[1]
+by Submodule with 3 generators
+1 -> x^2*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+```
+"""
+function saturation(M::SubquoModule, J::Ideal = ideal(base_ring(M), gens(base_ring(M))); iteration::Bool = false)
   @assert base_ring(M) == base_ring(J)
-  F = ambient_free_module(M)
+  M_quo = isdefined(M, :quo) ? M.quo : SubModuleOfFreeModule(ambient_free_module(M), Vector{elem_type(ambient_free_module(M))}())
   U = M.sub+M.quo
+  UF = _saturation(U, J; iteration = iteration)
+  res = SubquoModule(UF)
+  res.quo = M.quo
+  return simplify_light(res)[1]
+end
+
+function _saturation(U::SubModuleOfFreeModule, J::Ideal) ### TODO Replace by generic method
+  error("not implemented for the given types of modules.")
+end
+
+function _saturation(U::SubModuleOfFreeModule, J::Ideal{T}; iteration::Bool = false) where T <: Union{MPolyRingElem, MPolyQuoRingElem}
+  F = ambient_free_module(U)
   SgU = singular_generators(U.gens)
-  singular_assure(J)
-  SQ = Singular.quotient(SgU, J.gens.S)
+  SgJ = singular_generators(J.gens)
+  SQ, _ = Singular.saturation(SgU, SgJ)
   MG = ModuleGens(F, SQ)
-  UF = SubModuleOfFreeModule(F, MG)
+  return SubModuleOfFreeModule(F, MG)
+end
+
+@doc raw"""
+     saturation_with_index(M::SubquoModule,
+               J::MPolyIdeal = ideal(base_ring(M), gens(base_ring(M)));
+               iteration::Bool = false)
+
+
+Return the saturation $M:J^{\infty}$ of $M$ with respect to $J$ and the smallest integer $k$ such that $I:J^k = I:J^{\infty}$ (saturation index).
+
+If the ideal `J` is not given, the ideal generated by the generators (variables) of `base_ring(M)` is used.
+
+Setting `iteration` to `true` only has an effect over rings of type `MPolyRing` or `MPolyQuoRing`. Over such rings,
+if `iteration` is set to `true`, the saturation is done by carrying out successive module quotient computations as
+suggested by the definition of saturation. Otherwise, a more sophisticated Gröbner basis approach is used which is typically
+faster. Applying the two approaches may lead to different generating sets of the saturation.
+
+!!! note
+    By definition,
+    $M:J^{\infty} = \{ a \in A \mid J^ka \subset M \text{ for some } k\geq 1 \} = \bigcup\limits_{k=1}^{\infty} (M:J^k).$
+    Here, $A$ is the ambient module of $M$.
+
+# Examples
+```jldoctest
+julia> R, (x, y, z) = polynomial_ring(QQ, [:x, :y, :z]);
+
+julia> F = free_module(R, 1);
+
+julia> AM = R[x;];
+
+julia> BM = R[x^2; y^3; z^4];
+
+julia> M = SubquoModule(F, AM, BM)
+Subquotient of Submodule with 1 generator
+1 -> x*e[1]
+by Submodule with 3 generators
+1 -> x^2*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+julia> J = ideal(R, [x, y, z])^2
+Ideal generated by
+  x^2
+  x*y
+  x*z
+  y^2
+  y*z
+  z^2
+
+julia> L = saturation_with_index(M, J);
+
+julia> L[1]
+Subquotient of Submodule with 1 generator
+1 -> e[1]
+by Submodule with 3 generators
+1 -> x^2*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+julia> L[2]
+3
+
+```
+
+```jldoctest
+julia> S, (x, y, z) = polynomial_ring(QQ, [:x, :y, :z]);
+
+julia> R, _ = quo(S, ideal(S, [x+y+z]));
+
+julia> F = free_module(R, 1);
+
+julia> AM = R[x;];
+
+julia> BM = R[x^2; y^3; z^4];
+
+julia> M = SubquoModule(F, AM, BM)
+Subquotient of Submodule with 1 generator
+1 -> (-y - z)*e[1]
+by Submodule with 3 generators
+1 -> (y^2 + 2*y*z + z^2)*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+julia> J = ideal(R, [x, y, z])^2
+Ideal generated by
+  x^2
+  x*y
+  x*z
+  y^2
+  y*z
+  z^2
+
+julia> L = saturation_with_index(M, J);
+
+julia> L[1]
+Subquotient of Submodule with 1 generator
+1 -> e[1]
+by Submodule with 3 generators
+1 -> (y^2 + 2*y*z + z^2)*e[1]
+2 -> y^3*e[1]
+3 -> z^4*e[1]
+
+julia> L[2]
+2
+
+```
+"""
+function saturation_with_index(M::SubquoModule, J::Ideal = ideal(base_ring(M), gens(base_ring(M))); iteration::Bool = false)
+@assert base_ring(M) == base_ring(J)
+  M_quo = isdefined(M, :quo) ? M.quo : SubModuleOfFreeModule(ambient_free_module(M), Vector{elem_type(ambient_free_module(M))}())
+  U = M.sub+M.quo
+  UF, k = _saturation_with_index(U, J; iteration = iteration)
   res = SubquoModule(UF)
   res.quo = M.quo
-  return res
+  return simplify_light(res)[1], k
 end
 
+function _saturation_with_index(U::SubModuleOfFreeModule, J::Ideal) ### TODO Replace by generic method
+  error("not implemented for the given types of modules.")
+end
 
+function _saturation_with_index(U::SubModuleOfFreeModule, J::Ideal{T}; iteration::Bool = false) where T <: Union{MPolyRingElem, MPolyQuoRingElem}
+  F = ambient_free_module(U)
+  SgU = singular_generators(U.gens)
+  SgJ = singular_generators(J.gens)
+  SQ, k = Singular.saturation(SgU, SgJ)
+  MG = ModuleGens(F, SQ)
+  return SubModuleOfFreeModule(F, MG), k
+end
 
 ########################################
 
+
 @doc raw"""
     represents_element(a::FreeModElem, SQ::SubquoModule)
 

From e08f09c2b099bc19293e6bd0ba061c6a0030457c Mon Sep 17 00:00:00 2001
From: Wolfram Decker <decker@mathematik.uni-kl.de>
Date: Sun, 27 Oct 2024 12:17:03 +0100
Subject: [PATCH 2/4] remove whitespace

---
 .../ModulesOverMultivariateRings/subquotients.md                | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md b/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md
index 5f2c5c7de34..5dfd40282bf 100644
--- a/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md
+++ b/docs/src/CommutativeAlgebra/ModulesOverMultivariateRings/subquotients.md
@@ -363,7 +363,7 @@ quotient(M::SubquoModule{T}, N::SubquoModule{T}) where T
 ```
 
 ```@docs
-quotient(M:: SubquoModule, J::Ideal)
+quotient(M::SubquoModule, J::Ideal)
 ```
 
 ```@docs

From bdf8898fb1f8ea792366794a14f1e89be7fdbbea Mon Sep 17 00:00:00 2001
From: Wolfram Decker <decker@mathematik.uni-kl.de>
Date: Sun, 27 Oct 2024 12:33:12 +0100
Subject: [PATCH 3/4] Correction

---
 src/Modules/UngradedModules/SubquoModule.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/src/Modules/UngradedModules/SubquoModule.jl b/src/Modules/UngradedModules/SubquoModule.jl
index 1721a1091e1..6ac9ebeccc0 100644
--- a/src/Modules/UngradedModules/SubquoModule.jl
+++ b/src/Modules/UngradedModules/SubquoModule.jl
@@ -1771,7 +1771,7 @@ Ideal generated by
 function quotient(M::SubquoModule{T}, N::SubquoModule{T}) where T
   @assert base_ring(M) == base_ring(N)
   @assert ambient_module(M) == ambient_module(N)
-  MplusN = M+N
+  MplusN, iM, _ = sum(M, N)
   Q, _ = quo(MplusN, [iM(x) for x in gens(M)])
   return annihilator(Q)
 end

From 6faf888e33d756db703e3431d7f77f4f6f8a68a2 Mon Sep 17 00:00:00 2001
From: Wolfram Decker <decker@mathematik.uni-kl.de>
Date: Mon, 28 Oct 2024 11:40:12 +0100
Subject: [PATCH 4/4] reacting to the discussion

---
 .../src/CommutativeAlgebra/affine_algebras.md | 30 +------------------
 1 file changed, 1 insertion(+), 29 deletions(-)

diff --git a/docs/src/CommutativeAlgebra/affine_algebras.md b/docs/src/CommutativeAlgebra/affine_algebras.md
index 3d2a110c1b8..4c48d54b1aa 100644
--- a/docs/src/CommutativeAlgebra/affine_algebras.md
+++ b/docs/src/CommutativeAlgebra/affine_algebras.md
@@ -294,36 +294,8 @@ minimal_generating_set(I::MPolyQuoIdeal{<:MPolyDecRingElem})
 
 #### Intersection of Ideals
 
-```@julia
+```@docs
 intersect(a::MPolyQuoIdeal{T}, bs::MPolyQuoIdeal{T}...) where T
-intersect(V::Vector{MPolyQuoIdeal{T}}) where T
-```
-
-Return the intersection of two or more ideals.
-
-###### Examples
-
-```jldoctest
-julia> R, (x, y) = polynomial_ring(QQ, [:x, :y]);
-
-julia> A, _ = quo(R, ideal(R, [x^2-y^3, x-y]));
-
-julia> a = ideal(A, [y^2])
-Ideal generated by
-  y^2
-
-julia> b = ideal(A, [x])
-Ideal generated by
-  x
-
-julia> intersect(a,b)
-Ideal generated by
-  x*y
-
-julia> intersect([a,b])
-Ideal generated by
-x*y
-
 ```
 
 #### Ideal Quotients