Some more fixes

experimental/DoubleAndHyperComplexes/src/Morphisms/simplified_complexes.jl

@@ -509,12 +509,22 @@ function _simplify_matrix!(A::SMat; find_pivot=nothing)
     end
 
     # Adjust Sinv_transp
-    inc_row = sum(a*Sinv_transp[i] for (i, a) in a_col_del; init=sparse_row(R))
+    #inc_row = sum(a*Sinv_transp[i] for (i, a) in a_col_del; init=sparse_row(R))
+    inc_row = sparse_row(R)
+    for (i, a) in a_col_del
+      #is_zero(Sinv_transp[i]) && continue # can not be true since S is invertible
+      Hecke.add_scaled_row!(Sinv_transp[i], inc_row, a)
+    end
     Hecke.add_scaled_row!(inc_row, Sinv_transp[p], uinv)
 
     # Adjust T
     #T[q] = T[q] + sum(uinv*a*T[i] for (i, a) in a_row_del; init=sparse_row(R))
-    Hecke.add_scaled_row!(sum(a*T[i] for (i, a) in a_row_del if !iszero(T[i]); init=sparse_row(R)), T[q], uinv)
+    inc_row = sparse_row(R)
+    for (i, a) in a_row_del
+      # is_zero(T[i]) && continue # can not be true since T is invertible
+      Hecke.add_scaled_row!(T[i], inc_row, a)
+    end
+    Hecke.add_scaled_row!(inc_row, T[q], uinv)
 
     # Adjust Tinv_transp
     v = deepcopy(Tinv_transp[q])

experimental/Schemes/src/elliptic_surface.jl

@@ -847,7 +847,7 @@ Output a list of tuples with each tuple as follows
 - gram matrix of the intersection of [F0,...,Fn], it is an extended ADE-lattice.
 """
 function standardize_fiber(S::EllipticSurface, f::Vector{<:WeilDivisor})
-  @req all(is_prime(i) for i in f) "not a vector of prime divisors"
+  @hassert :EllipticSurface 2 all(is_prime(i) for i in f)
   f = copy(f)
   O = components(zero_section(S))[1]
   local f0
@@ -972,7 +972,7 @@ function fiber_components(S::EllipticSurface, P; algorithm=:exceptional_divisors
   end
   F = pullback(S.inc_Y, F)
   F = weil_divisor(F, ZZ)
-  fiber_components = [weil_divisor(E, ZZ) for E in EP]
+  fiber_components = [weil_divisor(E, ZZ; check=false) for E in EP]
   push!(fiber_components, F)
   return fiber_components
 end

src/Rings/mpoly-ideals.jl

@@ -572,7 +572,23 @@ Ideal generated by
   3*a*b*c
 ```
 """
-@attr T function radical(I::T) where {T <: MPolyIdeal}
+@attr MPolyIdeal{T} function radical(I::MPolyIdeal{T}) where {T <: MPolyRingElem}
+  is_known_to_be_radical(I) && return I
+  # Calling `elimpart` (within `simplify`) turns out to significantly speed things up in many cases.
+  R = base_ring(I)
+  Q, pr = quo(R, I)
+  S, iso, iso_inv = simplify(Q)
+  is_zero(ngens(S)) && return I # This only happens when all variables can be eliminated. 
+                                # Then the ideal defines a reduced point.
+  J = modulus(S)
+  pre_res = _compute_radical(J)
+  pre_res = ideal(R, elem_type(R)[lift(iso_inv(S(g))) for g in gens(pre_res)])
+  res = pre_res + ideal(R, [g for g in gens(I) if !(g in pre_res)])
+  set_attribute!(res, :is_radical=>true)
+  return res
+end
+
+function _compute_radical(I::T) where {T <: MPolyIdeal}
   R = base_ring(I)
   if isa(base_ring(R), NumField) && !isa(base_ring(R), AbsSimpleNumField)
     A, mA = absolute_simple_field(base_ring(R))
@@ -597,6 +613,7 @@ function radical(
     I::MPolyIdeal{T};
     factor_generators::Bool=true
   ) where {U<:Union{AbsSimpleNumFieldElem, <:Hecke.RelSimpleNumFieldElem}, T<:MPolyRingElem{U}}
+  is_known_to_be_radical(I) && return I
   get_attribute!(I, :radical) do
     is_one(I) && return I
     R = base_ring(I)
@@ -629,7 +646,7 @@ end
 
 function map_coefficients(mp, I::MPolyIdeal; parent = nothing)
   if parent === nothing
-    parent = Oscar.parent(map_coefficients(mp, gen(I, 1)))
+    parent = Oscar.parent(map_coefficients(mp, zero(base_ring(I))))
   end
   return ideal(parent, [map_coefficients(mp, g, parent = parent) for g = gens(I)])
 end
@@ -642,12 +659,17 @@ Return whether `I` is a radical ideal.
 Computes the radical.
 """
 @attr Bool function is_radical(I::MPolyIdeal)
-  if has_attribute(I, :is_prime) && is_prime(I)
-    return true
-  end
+  has_attribute(I, :is_prime) && is_prime(I) && return true
+  return is_subset(radical(I), I)
+end
 
-  return I == radical(I)
+function is_known_to_be_radical(I::Ideal)
+  has_attribute(I, :is_radical) && get_attribute(I, :is_radical) === true && return true
+  has_attribute(I, :is_prime) && get_attribute(I, :is_prime) === true && return true
+  return false
 end
+
+
 #######################################################
 @doc raw"""
     primary_decomposition(I::MPolyIdeal; algorithm = :GTZ, cache=true)
@@ -1149,7 +1171,14 @@ function minimal_primes(
       end
       J = K
     end
-    result = unique!(filter!(!is_one, vcat([minimal_primes(j; algorithm, factor_generators=false) for j in J]...)))
+    # unique! seems to fail here. We have to do it manually.
+    pre_result = filter!(!is_one, vcat([minimal_primes(j; algorithm, factor_generators=false) for j in J]...))
+    result = typeof(I)[]
+    for p in pre_result
+      p in result && continue
+      push!(result, p)
+    end
+
     # The list might not consist of minimal primes only. We have to discard the embedded ones
     final_list = typeof(I)[]
     for p in result
@@ -1209,7 +1238,7 @@ julia> L = equidimensional_decomposition_weak(I)
  Ideal with 1 generator
 ```
 """
-@attr Any function equidimensional_decomposition_weak(I::MPolyIdeal)
+@attr Vector{typeof(I)} function equidimensional_decomposition_weak(I::MPolyIdeal)
   R = base_ring(I)
   iszero(I) && return [I]
   @req coefficient_ring(R) isa AbstractAlgebra.Field "The coefficient ring must be a field"
@@ -1235,7 +1264,7 @@ julia> L = equidimensional_decomposition_weak(I)
   return V
 end
 
-@attr Any function equidimensional_decomposition_weak(
+@attr Vector{typeof(I)} function equidimensional_decomposition_weak(
     I::MPolyIdeal{T}
   ) where {U<:Union{AbsSimpleNumFieldElem, <:Hecke.RelSimpleNumFieldElem}, T<:MPolyRingElem{U}}
   R = base_ring(I)
@@ -2249,3 +2278,11 @@ function flag_pluecker_ideal(ring::MPolyRing{<: FieldElem}, dimensions::Vector{I
                    isReduced=true,
                    isGB=true))
 end
+
+# Since most ideals implement `==`, they have to implement the hash function.
+# See issue #4143 for problems entailed. Interestingly, this does not yet fix 
+# the failure of unique! on lists of ideals.
+function hash(I::Ideal, c::UInt)
+  return hash(typeof(I), hash(base_ring(I), c))
+end
+