From caaed2c1d9786d8fc9badba07e731ccd3c6eb077 Mon Sep 17 00:00:00 2001
From: Matthias Zach <85350711+HechtiDerLachs@users.noreply.github.com>
Date: Wed, 2 Oct 2024 23:01:24 +0200
Subject: [PATCH] Some more fixes (#4142)

* Speed up simplification of complexes.

* Make radical computation actually work over extensions of QQ.

* Add is_known_to_be_radical.

* Some improvements on elliptic surfaces.

* Implement base change for divisors.

* Extract types for resolutions of singularities to its own file.

* Make elimination an option in radical computation and update docstring.

* Speed up arithmetic in localized rings.

* Speed up restriction of maps by introducing _images_of_generators.

* Allow slightly more flexibility in the strict transform of closed embeddings.

* Disable isomorphism_on_open_subset.
---
 .../src/Morphisms/simplified_complexes.jl     |  14 +-
 experimental/Schemes/src/BlowupMorphism.jl    | 187 +++--------
 .../Schemes/src/BlowupMorphismTypes.jl        | 306 +++++++++++++++++
 experimental/Schemes/src/Schemes.jl           |   1 +
 experimental/Schemes/src/Types.jl             | 314 ------------------
 experimental/Schemes/src/elliptic_surface.jl  | 180 ++++++++--
 .../AffineSchemes/Morphisms/Attributes.jl     |  14 -
 .../AffineSchemes/Morphisms/Methods.jl        |  34 +-
 .../AffineSchemes/Objects/Attributes.jl       |   2 +-
 .../Schemes/AffineSchemes/Objects/Methods.jl  |  10 +-
 .../CoveredSchemes/Morphisms/Attributes.jl    |  16 +-
 .../Schemes/Covering/Morphisms/Methods.jl     |  10 +-
 .../Schemes/Covering/Objects/Methods.jl       |  12 +
 .../Schemes/Divisors/base_change.jl           |  19 ++
 .../Schemes/FunctionField/FunctionFields.jl   |  16 -
 .../Schemes/Gluing/Methods.jl                 |   2 +
 .../PrincipalOpenSubset/Objects/Methods.jl    |   1 +
 src/AlgebraicGeometry/Schemes/main.jl         |   1 +
 src/Rings/MPolyMap/MPolyAnyMap.jl             |  72 ++--
 src/Rings/MPolyMap/Types.jl                   |  47 +++
 src/Rings/MPolyQuo.jl                         |  10 +-
 src/Rings/Rings.jl                            |   1 +
 src/Rings/mpoly-ideals.jl                     | 104 +++++-
 src/Rings/mpoly-localizations.jl              |   4 +-
 .../Schemes/BlowupMorphism.jl                 |  11 -
 .../Schemes/CoveredProjectiveSchemes.jl       |   2 +-
 test/Rings/MPolyAnyMap/MPolyRing.jl           |   6 +-
 27 files changed, 773 insertions(+), 623 deletions(-)
 create mode 100644 experimental/Schemes/src/BlowupMorphismTypes.jl
 create mode 100644 src/AlgebraicGeometry/Schemes/Divisors/base_change.jl
 create mode 100644 src/Rings/MPolyMap/Types.jl

diff --git a/experimental/DoubleAndHyperComplexes/src/Morphisms/simplified_complexes.jl b/experimental/DoubleAndHyperComplexes/src/Morphisms/simplified_complexes.jl
index 7cd414decd53..f37708b45121 100644
--- a/experimental/DoubleAndHyperComplexes/src/Morphisms/simplified_complexes.jl
+++ b/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])
diff --git a/experimental/Schemes/src/BlowupMorphism.jl b/experimental/Schemes/src/BlowupMorphism.jl
index 92da05d1d4fa..166d8491aea2 100644
--- a/experimental/Schemes/src/BlowupMorphism.jl
+++ b/experimental/Schemes/src/BlowupMorphism.jl
@@ -3,7 +3,18 @@ export center
 export exceptional_divisor
 export projection
 
+########################################################################
+# `AbsDesingMor` and `AbsBlowupMorphism`
+#
+# An abstract type for blowup morphisms.
+#
+# This should also comprise sequences of simple blowups leading
+# to a partial or full resolution of singularities. The interface
+# is specified below.
+########################################################################
 
+### See experimental/Schemes/src/BlowupMorphismTypes.jl for the definition of 
+# `AbsDesingMor` and `AbsBlowupMorphism`.
 
 # The interface inherits all functionality from AbsCoveredSchemeMorphism.
 # This is extended by the following:
@@ -70,6 +81,16 @@ function total_transform(f::AbsBlowupMorphism, a::Any)
   error("not implemented")
 end
 
+########################################################################
+# `AbsSimpleBlowupMorphism`
+# An abstract type for classical blowups of ideal sheaves.
+#
+# This can either be a BlowupMorphism as below, but also a special
+# toric morphism induced by fan subdivisions.
+########################################################################
+
+### See BlowupMorphismTypes.jl for the definition of 
+# `AbsSimpleBlowupMorphism`
 
 @doc raw"""
     exceptional_divisor(f::AbsSimpleBlowupMorphism)
@@ -142,7 +163,16 @@ function controlled_transform(f::AbsSimpleBlowupMorphism, a::Any, k::Int)
   error("not implemented")
 end
 
+########################################################################
+# `BlowupMorphism`
+#
+# A datastructure to maintain all information necessary to effectively
+# handle blowups. This is work in progress and will one day serve as
+# a building block for sequences of blowups
+########################################################################
 
+### See experimental/Schemes/src/BlowupMorphismTypes.jl for the concrete type 
+# `BlowupMorphism`.
 
 ### Forward the essential functionality
 underlying_morphism(phi::BlowupMorphism) = projection(phi)
@@ -210,65 +240,8 @@ function strict_transform(p::AbsSimpleBlowupMorphism, inc::CoveredClosedEmbeddin
   inc_cod_cov = covering_morphism(inc)
 
   Z_trans = domain(inc_Z_trans)
-  pr_res = restrict(projection(p), inc_Z_trans, inc)
-
-  if has_attribute(p, :isomorphism_on_open_subset) # will only happen when p is a BlowupMorphism
-    OOX = OO(X)
-    OOY = OO(Y)
-    p_iso = isomorphism_on_open_subset(p)
-    U = domain(p_iso)::PrincipalOpenSubset
-    V = codomain(p_iso)::PrincipalOpenSubset
-    V_amb = ambient_scheme(V)
-    U_amb = ambient_scheme(U)
-
-    # We have the following diagram:
-    #             inc_dom
-    #     Z_trans    ↪     Y   ⊃ U_amb ⊃ U
-    #
-    #  pr_res↓             ↓ p           ↓ p_iso
-    #        Z       ↪     X   ⊃ V_amb ⊃ V
-    #             inc_cod
-    #
-    # Given all the refinements that potentially had to be done, we have
-    # to do the following.
-    #  1. Find a `patch` `U_sub` of the `domain` of `inc_dom_cov` for which
-    #    inc_dom : U_sub -> W has a codomain `W` which has `U_amb` as an
-    #    ancestor. Since `U_amb` is one of the `affine_charts` of `Y`, this will work.
-
-    k = findfirst(x->has_ancestor(y->y===U_amb, codomain(inc_dom_cov[x])), patches(domain(inc_dom_cov)))
-    U_sub = patches(domain(inc_dom_cov))[k]
-
-    #  2. pr_res : U_amb -> V_amb has some codomain such that there exists
-    #    an ancestor `VV` in the `domain` of `inc_dom_cov` such that
-    #    inc_dom : VV -> W' has a codomain with `V_amb` as an ancestor.
-    #  3. outside the `complement_equation`s of `U` and `V` the projection
-    #    was an isomorphism. Then the restriction of `pr_res` to those
-    #    complements in `U_sub` and `V_sub` also is.
-    U_sub_res = PrincipalOpenSubset(U_sub,
-                   pullback(inc_dom_cov[U_sub])(
-                       OOX(U_amb, codomain(inc_dom_cov[U_sub]))(
-                           complement_equation(U)
-                         )
-                     )
-                 )
-
-    pr_res_cov = covering_morphism(pr_res)
-    pr_sub = pr_res_cov[U_sub]
-
-    V_sub = codomain(pr_sub)
-    @assert has_ancestor(x->any(y->y===x, patches(domain(inc_cod_cov))), V_sub)
-    V_sub_inc, comp_eqns = _find_chart(V_sub, domain(inc_cod_cov))
-    OOZ = OO(Z)
-    dummy_dom = codomain(V_sub_inc)
-    dummy_map = inc_cod_cov[dummy_dom]
-    dummy_cod = codomain(dummy_map)
-    V_sub_res = PrincipalOpenSubset(V_sub,
-                                    OOZ(dummy_dom, V_sub)(pullback(dummy_map)(OOY(V_amb, dummy_cod)(complement_equation(V)))))
-    result = restrict(pr_sub, U_sub_res, V_sub_res)
-    # TODO: Obtain the inverse another way?
-    #@assert is_isomorphism(result)
-    set_attribute!(pr_res, :isomorphism_on_open_subset, result)
-  end
+  pr_res = restrict(projection(p), inc_Z_trans, inc; check=false)
+
   return Z_trans, inc_Z_trans, pr_res
 end
 
@@ -494,29 +467,35 @@ function restrict(f::AbsCoveredSchemeMorphism,
   inc_dom_cov = covering_morphism(inc_dom)
   inc_cod_cov = covering_morphism(inc_cod)
 
-  # We need to do the following.
-  # - Pass to a common refinement ref_cod in X that both
-  #   f and inc_cod can restrict to.
-  # - Pass to a common refinement in Y
-  ref_cod, a, b = _register!(common_refinement(codomain(f_cov), codomain(inc_cod_cov)), codomain(f))
-  inc_cod_ref = restrict(inc_cod, ref_cod)
-  f_res = restrict(f, ref_cod)
-  ref_dom, aa, bb = _register!(common_refinement(domain(f_res), codomain(inc_dom_cov)), domain(f))
-  inc_dom_ref = restrict(inc_dom, ref_dom)
-  inc_dom_ref = compose(inc_dom_ref, aa)
+  # Build up the common refinement 
+  success, dom_ref = is_refinement(codomain(inc_dom_cov), domain(f_cov))
+  @assert success "restriction not implemented for this constellation of refinements"
+  inc_dom_f = codomain(inc_dom_cov) === domain(f_cov) ? compose(inc_dom_cov, f_cov) : compose(compose(inc_dom_cov, dom_ref), f_cov)
+
+  success, cod_ref_map = is_refinement(codomain(inc_dom_f), codomain(inc_cod_cov))
+  inc_dom_f_ref = inc_dom_f # initialize the variable
+  if !success
+    success, cod_ref_map = is_refinement(codomain(inc_cod_cov), codomain(inc_dom_f))
+    @assert success "restriction not implemented for this constellation of refinements"
+    dom_ref2, to_inc_dom_f, to_inc_cod_cov = fiber_product(inc_dom_f, cod_ref_map)
+    inc_dom_f_ref = to_inc_cod_cov
+  else
+    inc_dom_f_ref = compose(inc_dom_f, cod_ref_map)
+  end
+
   # Collecting the maps for the restricted projection here
   map_dict = IdDict{AbsAffineScheme, AbsAffineSchemeMor}()
-  for U in patches(domain(inc_dom_ref))
-    q_res = compose(inc_dom_ref[U], f_res[codomain(inc_dom_ref[U])])
+  for U in patches(domain(inc_dom_f_ref))
+    q_res = inc_dom_f_ref[U]
     V = codomain(q_res)
-    g = maps_with_given_codomain(inc_cod_ref, V)
+    g = maps_with_given_codomain(inc_cod_cov, V)
     if !isone(length(g))
       error()
     end
     pre_V = domain(first(g))
-    map_dict[U] = restrict(q_res, domain(q_res), pre_V, check=false)
+    map_dict[U] = restrict(q_res, domain(q_res), pre_V; check)
   end
-  psi = CoveringMorphism(domain(inc_dom_ref), domain(inc_cod_ref), map_dict, check=false)
+  psi = CoveringMorphism(domain(inc_dom_f_ref), domain(inc_cod_cov), map_dict; check)
   return CoveredSchemeMorphism(domain(inc_dom), domain(inc_cod), psi)
 end
 
@@ -615,61 +594,6 @@ end
   return p_res
 end
 
-### Some functionality used by function fields
-
-# If set this attribute shall return some isomorphism on some open subsets of the domain
-# and codomain of phi. If both are irreducible, this automatically implies that both are
-# dense, but we do not check for this.
-@attr AbsAffineSchemeMor function isomorphism_on_open_subset(f::BlowupMorphism)
-  X = domain(f)
-  Y = codomain(f)
-  iso_dict = get_attribute(f, :isos_on_complement_of_center)
-  pr_res = first(values(iso_dict))
-  U = domain(pr_res)
-  V = codomain(pr_res)
-  iso_U = _flatten_open_subscheme(U, default_covering(X))
-  U_flat = codomain(iso_U)
-  iso_V = _flatten_open_subscheme(V, default_covering(Y))
-  V_flat = codomain(iso_V)
-  phi = morphism(U_flat, V_flat, pullback(inverse(iso_U)).(pullback(pr_res).(pullback(iso_V).(gens(OO(V_flat))))), check=false)
-  phi_inv = morphism(V_flat, U_flat, pullback(inverse(iso_V)).(pullback(inverse(pr_res)).(pullback(iso_U).(gens(OO(U_flat))))), check=false)
-  set_attribute!(phi, :inverse, phi_inv)
-  set_attribute!(phi_inv, :inverse, phi)
-  return phi
-end
-
-function pullback(f::BlowupMorphism, g::VarietyFunctionFieldElem)
-  X = domain(projection(f))
-  Y = codomain(projection(f))
-  FX = function_field(X)
-  FY = function_field(Y)
-  phi = isomorphism_on_open_subset(f)
-  U = ambient_scheme(domain(phi))
-  V = ambient_scheme(codomain(phi))
-  parent(g) === FY || error("element does not belong to the correct field")
-  h = g[V]
-  pbg = fraction(pullback(phi)(OO(V)(numerator(h))))//fraction(pullback(phi)(OO(V)(denominator(h))))
-  return FX.(pbg)
-end
-
-function pushforward(f::BlowupMorphism, g::VarietyFunctionFieldElem)
-  X = domain(projection(f))
-  Y = codomain(projection(f))
-  FX = function_field(X)
-  FY = function_field(Y)
-  phi = inverse(isomorphism_on_open_subset(f))
-  U = ambient_scheme(domain(phi))
-  V = ambient_scheme(codomain(phi))
-  parent(g) === FX || error("element does not belong to the correct field")
-  h = g[V]
-  pfg = fraction(pullback(phi)(OO(V)(numerator(h))))//fraction(pullback(phi)(OO(V)(denominator(h))))
-  return FY.(pfg)
-end
-
-@attr AbsAffineSchemeMor function isomorphism_on_open_subset(phi::AbsCoveredSchemeMorphism)
-  error("attribute not found; this needs to be set manually in general")
-end
-
 function compose(f::AbsSimpleBlowupMorphism, g::AbsSimpleBlowupMorphism)
   return composite_map(f, g)
 end
@@ -682,9 +606,7 @@ function compose(f::AbsCoveredSchemeMorphism, g::AbsSimpleBlowupMorphism)
   return composite_map(f, g)
 end
 
-
-
-
+### See experimental/Schemes/src/BlowupMorphismTypes.jl for declaration of `BlowUpSequence`.
 
 ########################################################################
 # Convenience method                                                   #
@@ -696,6 +618,7 @@ function blow_up(m::AbsCoveredScheme, I::AbsIdealSheaf; coordinate_name = "e")
 end
 
 
+### See experimental/Schemes/src/BlowupMorphismTypes.jl for the declaration of `StrictTransformIdealSheaf`.
 
 morphism(I::StrictTransformIdealSheaf) = I.morphism
 original_ideal_sheaf(I::StrictTransformIdealSheaf) = I.orig
diff --git a/experimental/Schemes/src/BlowupMorphismTypes.jl b/experimental/Schemes/src/BlowupMorphismTypes.jl
new file mode 100644
index 000000000000..63ace393a85f
--- /dev/null
+++ b/experimental/Schemes/src/BlowupMorphismTypes.jl
@@ -0,0 +1,306 @@
+@doc raw"""
+    AbsDesingMor{
+                                  DomainType<:AbsCoveredScheme,
+                                  CodomainType<:AbsCoveredScheme,
+                                  BlowupMorphismType
+       } <: AbsCoveredSchemeMorphism{
+                                 DomainType,
+                                 CodomainType,
+                                 Nothing,
+                                 BlowupMorphismType
+                                }
+Abstract type for desingularizations ``f : X -> Y `` of schemes where
+
+  * ``Y`` is the scheme of which the singularities are to be resolved
+  * ``f`` is a birational proper map 
+          may for instance be BlowUpSequence or Lipman-style combination of blow-ups and normalization
+  * ``Y`` is a regular scheme
+"""
+abstract type AbsDesingMor{
+                           DomainType<:AbsCoveredScheme,
+                           CodomainType<:AbsCoveredScheme,
+                           BlowupMorphismType
+                          } <: AbsCoveredSchemeMorphism{
+                                                        DomainType,
+                                                        CodomainType,
+                                                        Nothing,
+                                                        BlowupMorphismType
+                                                       }
+end
+
+########################################################################
+# An abstract type for blowup morphisms.
+#
+# This should also comprise sequences of simple blowups leading
+# to a partial or full resolution of singularities. The interface
+# is specified below.
+########################################################################
+abstract type AbsBlowupMorphism{DomainType<:AbsCoveredScheme,
+                                  CodomainType<:AbsCoveredScheme,
+                                  BlowupMorphismType
+                                 } <: AbsDesingMor{DomainType, CodomainType, BlowupMorphismType}
+end
+
+abstract type AbsSimpleBlowupMorphism{DomainType<:AbsCoveredScheme,
+                                     CodomainType<:AbsCoveredScheme,
+                                     BlowupMorphismType
+    } <: AbsBlowupMorphism{
+                             DomainType,
+                             CodomainType,
+                             BlowupMorphismType
+                            }
+end
+
+
+########################################################################
+# BlowupMorphism
+#
+# A datastructure to maintain all information necessary to effectively
+# handle blowups. This is work in progress and will one day serve as
+# a building block for sequences of blowups
+########################################################################
+
+@doc raw"""
+    BlowupMorphism
+
+A datastructure to encode blowups of covered schemes in some sheaves of ideals.
+
+It is described as a morphism from the new scheme to the blown-up scheme, with
+information about its center (i.e. the ideal sheaves blown-up in the bottom
+scheme) and its exceptional locus (i.e. the preimage of the center under the
+blowup).
+
+# Examples
+```jldoctest
+julia> R, (x,y,z) = QQ[:x, :y, :z];
+
+julia> A3 = spec(R)
+Spectrum
+  of multivariate polynomial ring in 3 variables x, y, z
+    over rational field
+
+julia> I = ideal(R, [x,y,z])
+Ideal generated by
+  x
+  y
+  z
+
+julia> bl = blow_up(A3, I)
+Blowup
+  of scheme over QQ covered with 1 patch
+    1b: [x, y, z]   affine 3-space
+  in sheaf of ideals with restriction
+    1b: Ideal (x, y, z)
+with domain
+  scheme over QQ covered with 3 patches
+    1a: [(s1//s0), (s2//s0), x]   scheme(0, 0, 0)
+    2a: [(s0//s1), (s2//s1), y]   scheme(0, 0, 0)
+    3a: [(s0//s2), (s1//s2), z]   scheme(0, 0, 0)
+and exceptional divisor
+  effective cartier divisor defined by
+    sheaf of ideals with restrictions
+      1a: Ideal (x)
+      2a: Ideal (y)
+      3a: Ideal (z)
+
+julia> E = exceptional_divisor(bl)
+Effective cartier divisor
+  on scheme over QQ covered with 3 patches
+    1: [(s1//s0), (s2//s0), x]   scheme(0, 0, 0)
+    2: [(s0//s1), (s2//s1), y]   scheme(0, 0, 0)
+    3: [(s0//s2), (s1//s2), z]   scheme(0, 0, 0)
+defined by
+  sheaf of ideals with restrictions
+    1: Ideal (x)
+    2: Ideal (y)
+    3: Ideal (z)
+
+julia> Z = center(bl)
+Sheaf of ideals
+  on scheme over QQ covered with 1 patch
+    1: [x, y, z]   affine 3-space
+with restriction
+  1: Ideal (x, y, z)
+```
+"""
+@attributes mutable struct BlowupMorphism{
+     DomainType<:AbsCoveredScheme, # Not a concrete type in general because this is lazy
+     CodomainType<:AbsCoveredScheme,
+   } <: AbsSimpleBlowupMorphism{
+                                  DomainType,
+                                  CodomainType,
+                                  BlowupMorphism{DomainType, CodomainType}
+                                 }
+  projective_bundle::CoveredProjectiveScheme
+  codomain::CodomainType   # in general a CoveredScheme
+  center::AbsIdealSheaf      # on codomain
+  projection::AbsCoveredSchemeMorphism
+  domain::AbsCoveredScheme # in general a CoveredScheme
+  exceptional_divisor::EffectiveCartierDivisor
+
+  function BlowupMorphism(
+      IP::CoveredProjectiveScheme,
+      I::AbsIdealSheaf
+    )
+    X = base_scheme(IP)
+    X === scheme(I) || error("ideal sheaf not compatible with blown up variety")
+    return new{AbsCoveredScheme, typeof(X)}(IP, X, I)
+  end
+end
+
+########################################################################
+# Resolutions of singularities                                         #
+########################################################################
+
+@doc raw"""
+    BlowUpSequence{
+    DomainType<:AbsCoveredScheme,
+    CodomainType<:AbsCoveredScheme
+   } <: AbsDesingMor{
+                                 DomainType,
+                                 CodomainType,
+                                }
+
+
+"""
+@attributes mutable struct BlowUpSequence{
+                                          DomainType<:AbsCoveredScheme,
+                                          CodomainType<:AbsCoveredScheme
+                                         }<:AbsBlowupMorphism{
+                                                                DomainType, CodomainType, 
+                                                                BlowUpSequence{DomainType, CodomainType}
+                                                               }
+  maps::Vector{<:BlowupMorphism}                 # count right to left:
+                                                 # original scheme is codomain of map 1
+  
+  embeddings::Vector{<:AbsCoveredSchemeMorphism} # if set,
+                                                 # assert codomain(maps[i])===codomain(embeddings[i]) 
+  # boolean flags
+  is_embedded::Bool                              # do not set embeddings, ex_mult, controlled_transform etc
+                                                 #     if is_embedded == false
+  resolves_sing::Bool                            # domain(maps[end]) smooth?
+  is_trivial::Bool                               # codomain already smooth?
+  is_strong::Bool                                # snc divisors ensured?
+  transform_type::Symbol                         # can be :strict, :weak or :control
+                                                 #     only relevant for is_embedded == true
+
+  # fields for caching, may be filled during computation
+  ex_div::Vector{<:EffectiveCartierDivisor}      # list of exc. divisors arising from individual steps
+                                                 # lives in domain(maps[end])
+  control::Int                                   # value of control for controlled transform
+  ex_mult::Vector{Int}                           # multiplicities of exceptional divisors removed from
+                                                 # total transform, not set for is_embedded == false
+                                                 # or transform_type == strict
+  controlled_transform::AbsIdealSheaf            # holds weak or controlled transform according to transform_type
+  dont_meet::Vector{Tuple{Int,Int}}              # mostly for dim=2: intersections which cannot exist according
+                                                 # to intermediate computations
+  caution_multi_charts::Vector{Tuple{Int,Int}}   # only for dim=2: intersection of divisors not
+                                                 # entirely visible in a single chart
+
+
+  # fields for caching to be filled a posteriori (on demand, only if partial_res==false)
+  underlying_morphism::CompositeCoveredSchemeMorphism{DomainType, CodomainType}
+  exceptional_divisor::CartierDivisor            # exceptional divisor of composed_map
+  exceptional_locus::WeilDivisor                 # exceptional locus of composed map
+  exceptional_divisor_on_X::CartierDivisor          # exceptional divisor of composed_map
+                                                 # restricted to domain(embeddings[end])
+
+  function BlowUpSequence(maps::Vector{<:BlowupMorphism})
+    n = length(maps)
+    for i in 1:n-1
+      @assert domain(maps[i]) === codomain(maps[i+1]) "not a sequence of morphisms"
+    end
+    return new{typeof(domain(maps[end])),typeof(codomain(first(maps)))}(maps)
+  end
+end
+
+########################################################################
+# strict transforms of ideal sheaves                                   #
+########################################################################
+
+@attributes mutable struct StrictTransformIdealSheaf{SpaceType, OpenType, OutputType,
+                                                     RestrictionType
+                                                    } <: AbsIdealSheaf{
+                                                                       SpaceType, OpenType,
+                                                                       OutputType, RestrictionType
+                                                                      }
+  morphism::AbsSimpleBlowupMorphism
+  orig::AbsIdealSheaf
+  underlying_presheaf::AbsPreSheaf
+
+  function StrictTransformIdealSheaf(
+      f::AbsSimpleBlowupMorphism,
+      J::AbsIdealSheaf
+    )
+    @assert scheme(J) === codomain(f)
+    X = domain(f)
+    Ipre = PreSheafOnScheme(X,
+                      OpenType=AbsAffineScheme, OutputType=Ideal,
+                      RestrictionType=Map,
+                      is_open_func=_is_open_func_for_schemes_without_affine_scheme_open_subscheme(X)
+                     )
+    I = new{typeof(X), AbsAffineScheme, Ideal, Map}(f, J, Ipre)
+    return I
+  end
+end
+
+#####################################################################################################
+# Desingularization morphism: birational map between covered schemes with smooth domain
+#####################################################################################################
+@doc raw"""
+    MixedBlowUpSequence{
+              DomainType<:AbsCoveredScheme,
+              CodomainType<:AbsCoveredScheme
+            }<:AbsDesingMor{  DomainType,
+                              CodomainType, 
+                              MixedBlowUpSequence{DomainType, CodomainType}
+                                              }
+A datastructure to encode sequences of blow-ups and normalizations of covered schemes
+as needed for desingularization of non-embedded schemes by the approaches of Zariski and of
+Lipman. 
+"""
+
+@attributes mutable struct MixedBlowUpSequence{
+                                                DomainType<:AbsCoveredScheme,
+                                                CodomainType<:AbsCoveredScheme
+                                              }<:AbsDesingMor{  DomainType,
+                                                                CodomainType, 
+                                                                MixedBlowUpSequence{DomainType, CodomainType}
+                                              }
+  maps::Vector{Union{<:BlowupMorphism,<:NormalizationMorphism}}       # count right to left:
+                                                 # original scheme is codomain of map 1
+  # boolean flags
+  resolves_sing::Bool                            # domain not smooth yet?
+  is_trivial::Bool                               # codomain already smooth?
+  is_strong::Bool                                # snc divisors ensured?
+
+  # fields for caching, to be filled during desingularization
+  # always carried along to domain(maps[end])) using strict_transform
+  ex_div::Vector{AbsIdealSheaf}      # list of exc. divisors arising from individual steps
+  dont_meet::Vector{Tuple{Int,Int}}             # mostly for dim=2: intersections which cannot exist
+                                                # according to intermediate computations
+  caution_multi_charts::Vector{Tuple{Int,Int}}  # only for dim=2: intersection of divisors not
+                                                # entirely visible in a single chart
+
+  # keep track of the normalization steps
+  normalization_steps::Vector{Int}
+
+  # fields for caching to be filled a posteriori (on demand, only if partial_res==false)
+  underlying_morphism::CompositeCoveredSchemeMorphism{DomainType, CodomainType}
+  exceptional_divisor::AbsWeilDivisor
+  exceptional_locus::AbsAlgebraicCycle
+
+  function MixedBlowUpSequence(maps::Vector{<:AbsCoveredSchemeMorphism})
+    n = length(maps)
+    for i in 1:n
+      @assert all(x->((x isa BlowupMorphism) || (x isa NormalizationMorphism)), maps) "only blow-ups and normalizations allowed"
+    end
+    for i in 1:n-1
+      @assert domain(maps[i]) === codomain(maps[i+1]) "not a sequence of morphisms"
+    end
+    resi = new{typeof(domain(maps[end])),typeof(codomain(first(maps)))}(maps)
+    resi.normalization_steps = [i for i in 1:n if maps[i] isa NormalizationMorphism]
+    return resi
+  end
+
+end
diff --git a/experimental/Schemes/src/Schemes.jl b/experimental/Schemes/src/Schemes.jl
index b5e6a87a6f48..215ff86cd640 100644
--- a/experimental/Schemes/src/Schemes.jl
+++ b/experimental/Schemes/src/Schemes.jl
@@ -6,6 +6,7 @@ include("Tjurina.jl")
 include("CoveredProjectiveSchemes.jl")
 
 include("Auxiliary.jl")
+include("BlowupMorphismTypes.jl")
 include("BlowupMorphism.jl")
 include("elliptic_surface.jl")
 include("critical_locus.jl")
diff --git a/experimental/Schemes/src/Types.jl b/experimental/Schemes/src/Types.jl
index 447d5ab178c5..59396ce9c518 100644
--- a/experimental/Schemes/src/Types.jl
+++ b/experimental/Schemes/src/Types.jl
@@ -183,317 +183,3 @@ end
 end
 
 
-@doc raw"""
-    AbsDesingMor{
-                                  DomainType<:AbsCoveredScheme,
-                                  CodomainType<:AbsCoveredScheme,
-                                  BlowupMorphismType
-       } <: AbsCoveredSchemeMorphism{
-                                 DomainType,
-                                 CodomainType,
-                                 Nothing,
-                                 BlowupMorphismType
-                                }
-Abstract type for desingularizations ``f : X -> Y `` of schemes where
-
-  * ``Y`` is the scheme of which the singularities are to be resolved
-  * ``f`` is a birational proper map 
-          may for instance be BlowUpSequence or Lipman-style combination of blow-ups and normalization
-  * ``Y`` is a regular scheme
-"""
-abstract type AbsDesingMor{
-                           DomainType<:AbsCoveredScheme,
-                           CodomainType<:AbsCoveredScheme,
-                           BlowupMorphismType
-                          } <: AbsCoveredSchemeMorphism{
-                                                        DomainType,
-                                                        CodomainType,
-                                                        Nothing,
-                                                        BlowupMorphismType
-                                                       }
-end
-
-########################################################################
-# An abstract type for blowup morphisms.
-#
-# This should also comprise sequences of simple blowups leading
-# to a partial or full resolution of singularities. The interface
-# is specified below.
-########################################################################
-abstract type AbsBlowupMorphism{DomainType<:AbsCoveredScheme,
-                                  CodomainType<:AbsCoveredScheme,
-                                  BlowupMorphismType
-                                 } <: AbsDesingMor{DomainType, CodomainType, BlowupMorphismType}
-end
-
-########################################################################
-# An abstract type for classical blowups of ideal sheaves.
-#
-# This can either be a BlowupMorphism as below, but also a special
-# toric morphism induced by fan subdivisions.
-########################################################################
-abstract type AbsSimpleBlowupMorphism{DomainType<:AbsCoveredScheme,
-                                     CodomainType<:AbsCoveredScheme,
-                                     BlowupMorphismType
-    } <: AbsBlowupMorphism{
-                             DomainType,
-                             CodomainType,
-                             BlowupMorphismType
-                            }
-end
-
-
-########################################################################
-# BlowupMorphism
-#
-# A datastructure to maintain all information necessary to effectively
-# handle blowups. This is work in progress and will one day serve as
-# a building block for sequences of blowups
-########################################################################
-
-@doc raw"""
-    BlowupMorphism
-
-A datastructure to encode blowups of covered schemes in some sheaves of ideals.
-
-It is described as a morphism from the new scheme to the blown-up scheme, with
-information about its center (i.e. the ideal sheaves blown-up in the bottom
-scheme) and its exceptional locus (i.e. the preimage of the center under the
-blowup).
-
-# Examples
-```jldoctest
-julia> R, (x,y,z) = QQ[:x, :y, :z];
-
-julia> A3 = spec(R)
-Spectrum
-  of multivariate polynomial ring in 3 variables x, y, z
-    over rational field
-
-julia> I = ideal(R, [x,y,z])
-Ideal generated by
-  x
-  y
-  z
-
-julia> bl = blow_up(A3, I)
-Blowup
-  of scheme over QQ covered with 1 patch
-    1b: [x, y, z]   affine 3-space
-  in sheaf of ideals with restriction
-    1b: Ideal (x, y, z)
-with domain
-  scheme over QQ covered with 3 patches
-    1a: [(s1//s0), (s2//s0), x]   scheme(0, 0, 0)
-    2a: [(s0//s1), (s2//s1), y]   scheme(0, 0, 0)
-    3a: [(s0//s2), (s1//s2), z]   scheme(0, 0, 0)
-and exceptional divisor
-  effective cartier divisor defined by
-    sheaf of ideals with restrictions
-      1a: Ideal (x)
-      2a: Ideal (y)
-      3a: Ideal (z)
-
-julia> E = exceptional_divisor(bl)
-Effective cartier divisor
-  on scheme over QQ covered with 3 patches
-    1: [(s1//s0), (s2//s0), x]   scheme(0, 0, 0)
-    2: [(s0//s1), (s2//s1), y]   scheme(0, 0, 0)
-    3: [(s0//s2), (s1//s2), z]   scheme(0, 0, 0)
-defined by
-  sheaf of ideals with restrictions
-    1: Ideal (x)
-    2: Ideal (y)
-    3: Ideal (z)
-
-julia> Z = center(bl)
-Sheaf of ideals
-  on scheme over QQ covered with 1 patch
-    1: [x, y, z]   affine 3-space
-with restriction
-  1: Ideal (x, y, z)
-```
-"""
-@attributes mutable struct BlowupMorphism{
-     DomainType<:AbsCoveredScheme, # Not a concrete type in general because this is lazy
-     CodomainType<:AbsCoveredScheme,
-   } <: AbsSimpleBlowupMorphism{
-                                  DomainType,
-                                  CodomainType,
-                                  BlowupMorphism{DomainType, CodomainType}
-                                 }
-  projective_bundle::CoveredProjectiveScheme
-  codomain::CodomainType   # in general a CoveredScheme
-  center::AbsIdealSheaf      # on codomain
-  projection::AbsCoveredSchemeMorphism
-  domain::AbsCoveredScheme # in general a CoveredScheme
-  exceptional_divisor::EffectiveCartierDivisor
-
-  function BlowupMorphism(
-      IP::CoveredProjectiveScheme,
-      I::AbsIdealSheaf
-    )
-    X = base_scheme(IP)
-    X === scheme(I) || error("ideal sheaf not compatible with blown up variety")
-    return new{AbsCoveredScheme, typeof(X)}(IP, X, I)
-  end
-end
-
-########################################################################
-# Resolutions of singularities                                         #
-########################################################################
-
-@doc raw"""
-    BlowUpSequence{
-    DomainType<:AbsCoveredScheme,
-    CodomainType<:AbsCoveredScheme
-   } <: AbsDesingMor{
-                                 DomainType,
-                                 CodomainType,
-                                }
-
-
-"""
-@attributes mutable struct BlowUpSequence{
-                                          DomainType<:AbsCoveredScheme,
-                                          CodomainType<:AbsCoveredScheme
-                                         }<:AbsBlowupMorphism{
-                                                                DomainType, CodomainType, 
-                                                                BlowUpSequence{DomainType, CodomainType}
-                                                               }
-  maps::Vector{<:BlowupMorphism}                 # count right to left:
-                                                 # original scheme is codomain of map 1
-  
-  embeddings::Vector{<:AbsCoveredSchemeMorphism} # if set,
-                                                 # assert codomain(maps[i])===codomain(embeddings[i]) 
-  # boolean flags
-  is_embedded::Bool                              # do not set embeddings, ex_mult, controlled_transform etc
-                                                 #     if is_embedded == false
-  resolves_sing::Bool                            # domain(maps[end]) smooth?
-  is_trivial::Bool                               # codomain already smooth?
-  is_strong::Bool                                # snc divisors ensured?
-  transform_type::Symbol                         # can be :strict, :weak or :control
-                                                 #     only relevant for is_embedded == true
-
-  # fields for caching, may be filled during computation
-  ex_div::Vector{<:EffectiveCartierDivisor}      # list of exc. divisors arising from individual steps
-                                                 # lives in domain(maps[end])
-  control::Int                                   # value of control for controlled transform
-  ex_mult::Vector{Int}                           # multiplicities of exceptional divisors removed from
-                                                 # total transform, not set for is_embedded == false
-                                                 # or transform_type == strict
-  controlled_transform::AbsIdealSheaf            # holds weak or controlled transform according to transform_type
-  dont_meet::Vector{Tuple{Int,Int}}              # mostly for dim=2: intersections which cannot exist according
-                                                 # to intermediate computations
-  caution_multi_charts::Vector{Tuple{Int,Int}}   # only for dim=2: intersection of divisors not
-                                                 # entirely visible in a single chart
-
-
-  # fields for caching to be filled a posteriori (on demand, only if partial_res==false)
-  underlying_morphism::CompositeCoveredSchemeMorphism{DomainType, CodomainType}
-  exceptional_divisor::CartierDivisor            # exceptional divisor of composed_map
-  exceptional_locus::WeilDivisor                 # exceptional locus of composed map
-  exceptional_divisor_on_X::CartierDivisor          # exceptional divisor of composed_map
-                                                 # restricted to domain(embeddings[end])
-
-  function BlowUpSequence(maps::Vector{<:BlowupMorphism})
-    n = length(maps)
-    for i in 1:n-1
-      @assert domain(maps[i]) === codomain(maps[i+1]) "not a sequence of morphisms"
-    end
-    return new{typeof(domain(maps[end])),typeof(codomain(first(maps)))}(maps)
-  end
-end
-
-########################################################################
-# strict transforms of ideal sheaves                                   #
-########################################################################
-
-@attributes mutable struct StrictTransformIdealSheaf{SpaceType, OpenType, OutputType,
-                                                     RestrictionType
-                                                    } <: AbsIdealSheaf{
-                                                                       SpaceType, OpenType,
-                                                                       OutputType, RestrictionType
-                                                                      }
-  morphism::AbsSimpleBlowupMorphism
-  orig::AbsIdealSheaf
-  underlying_presheaf::AbsPreSheaf
-
-  function StrictTransformIdealSheaf(
-      f::AbsSimpleBlowupMorphism,
-      J::AbsIdealSheaf
-    )
-    @assert scheme(J) === codomain(f)
-    X = domain(f)
-    Ipre = PreSheafOnScheme(X,
-                      OpenType=AbsAffineScheme, OutputType=Ideal,
-                      RestrictionType=Map,
-                      is_open_func=_is_open_func_for_schemes_without_affine_scheme_open_subscheme(X)
-                     )
-    I = new{typeof(X), AbsAffineScheme, Ideal, Map}(f, J, Ipre)
-    return I
-  end
-end
-
-#####################################################################################################
-# Desingularization morphism: birational map between covered schemes with smooth domain
-#####################################################################################################
-@doc raw"""
-    MixedBlowUpSequence{
-              DomainType<:AbsCoveredScheme,
-              CodomainType<:AbsCoveredScheme
-            }<:AbsDesingMor{  DomainType,
-                              CodomainType, 
-                              MixedBlowUpSequence{DomainType, CodomainType}
-                                              }
-A datastructure to encode sequences of blow-ups and normalizations of covered schemes
-as needed for desingularization of non-embedded schemes by the approaches of Zariski and of
-Lipman. 
-"""
-
-@attributes mutable struct MixedBlowUpSequence{
-                                                DomainType<:AbsCoveredScheme,
-                                                CodomainType<:AbsCoveredScheme
-                                              }<:AbsDesingMor{  DomainType,
-                                                                CodomainType, 
-                                                                MixedBlowUpSequence{DomainType, CodomainType}
-                                              }
-  maps::Vector{Union{<:BlowupMorphism,<:NormalizationMorphism}}       # count right to left:
-                                                 # original scheme is codomain of map 1
-  # boolean flags
-  resolves_sing::Bool                            # domain not smooth yet?
-  is_trivial::Bool                               # codomain already smooth?
-  is_strong::Bool                                # snc divisors ensured?
-
-  # fields for caching, to be filled during desingularization
-  # always carried along to domain(maps[end])) using strict_transform
-  ex_div::Vector{AbsIdealSheaf}      # list of exc. divisors arising from individual steps
-  dont_meet::Vector{Tuple{Int,Int}}             # mostly for dim=2: intersections which cannot exist
-                                                # according to intermediate computations
-  caution_multi_charts::Vector{Tuple{Int,Int}}  # only for dim=2: intersection of divisors not
-                                                # entirely visible in a single chart
-
-  # keep track of the normalization steps
-  normalization_steps::Vector{Int}
-
-  # fields for caching to be filled a posteriori (on demand, only if partial_res==false)
-  underlying_morphism::CompositeCoveredSchemeMorphism{DomainType, CodomainType}
-  exceptional_divisor::AbsWeilDivisor
-  exceptional_locus::AbsAlgebraicCycle
-
-  function MixedBlowUpSequence(maps::Vector{<:AbsCoveredSchemeMorphism})
-    n = length(maps)
-    for i in 1:n
-      @assert all(x->((x isa BlowupMorphism) || (x isa NormalizationMorphism)), maps) "only blow-ups and normalizations allowed"
-    end
-    for i in 1:n-1
-      @assert domain(maps[i]) === codomain(maps[i+1]) "not a sequence of morphisms"
-    end
-    resi = new{typeof(domain(maps[end])),typeof(codomain(first(maps)))}(maps)
-    resi.normalization_steps = [i for i in 1:n if maps[i] isa NormalizationMorphism]
-    return resi
-  end
-
-end
-
-    
diff --git a/experimental/Schemes/src/elliptic_surface.jl b/experimental/Schemes/src/elliptic_surface.jl
index d76d945f7d6f..e9144342748c 100644
--- a/experimental/Schemes/src/elliptic_surface.jl
+++ b/experimental/Schemes/src/elliptic_surface.jl
@@ -25,6 +25,7 @@ For now functionality is restricted to $C = \mathbb{P}^1$.
   inc_Weierstrass::CoveredClosedEmbedding # inclusion of the weierstrass chart in its ambient projective bundle
   inc_Y::CoveredClosedEmbedding # inclusion of Y in its ambient blown up projective bundle
   euler_characteristic::Int
+  resolution_strategy::Symbol
   # the following are temporary until we have a dedicated type for
   # iterated blow ups
   blowup::AbsCoveredSchemeMorphism
@@ -35,7 +36,12 @@ For now functionality is restricted to $C = \mathbb{P}^1$.
   fibration::AbsCoveredSchemeMorphism # the projection to IP^1
   fibration_weierstrass_model::AbsCoveredSchemeMorphism # the projection from the Weierstrass model
 
-  function EllipticSurface(generic_fiber::EllipticCurve{F}, euler_characteristic::Int, mwl_basis::Vector{<:EllipticCurvePoint}) where F
+  function EllipticSurface(
+      generic_fiber::EllipticCurve{F}, 
+      euler_characteristic::Int, 
+      mwl_basis::Vector{<:EllipticCurvePoint};
+      resolution_strategy::Symbol=:iterative
+    ) where F
     B = typeof(coefficient_ring(base_ring(base_field(generic_fiber))))
     S = new{B,F}()
     S.E = generic_fiber
@@ -45,6 +51,7 @@ For now functionality is restricted to $C = \mathbb{P}^1$.
     set_attribute!(S, :is_reduced=>true)
     set_attribute!(S, :is_integral=>true)
     set_attribute!(S, :is_equidimensional=>true)
+    S.resolution_strategy = resolution_strategy
     return S
   end
 
@@ -111,10 +118,11 @@ with generic fiber
 function elliptic_surface(generic_fiber::EllipticCurve{BaseField},
                           euler_characteristic::Int,
                           mwl_gens::Vector{<:EllipticCurvePoint}=EllipticCurvePoint[];
+                          resolution_strategy::Symbol=:iterative,
                           is_basis::Bool=true) where {
                           BaseField <: FracFieldElem{<:PolyRingElem{<:FieldElem}}}
   @req all(parent(i)==generic_fiber for i in mwl_gens) "not a vector of points on $(generic_fiber)"
-  S = EllipticSurface(generic_fiber, euler_characteristic, mwl_gens)
+  S = EllipticSurface(generic_fiber, euler_characteristic, mwl_gens; resolution_strategy)
   if is_basis
     return S
   end
@@ -585,7 +593,141 @@ This triggers the computation of the `underlying_scheme` of ``X``
 as a blowup from its Weierstrass model. It may take a few minutes.
 """
 function weierstrass_contraction(X::EllipticSurface)
-  Y = X
+  algorithm = X.resolution_strategy
+  if algorithm == :iterative
+    return weierstrass_contraction_iterative(X)
+  elseif algorithm == :simultaneous
+    return weierstrass_contraction_simultaneous(X)
+  else
+    error("algorithm not recognized")
+  end
+end
+
+function weierstrass_contraction_simultaneous(Y::EllipticSurface)
+  if isdefined(Y, :blowup)
+    return Y.blowup
+  end
+  S, inc_S = weierstrass_model(Y)
+  @assert has_attribute(S, :is_equidimensional) && get_attribute(S, :is_equidimensional) === true
+  
+  X0 = codomain(inc_S)
+  Y0 = S
+  set_attribute!(Y0, :is_reduced=>true)
+  set_attribute!(Y0, :is_irreducible=>true)
+  set_attribute!(Y0, :is_equidimensional=>true)
+  inc_Y0 = inc_S
+  I_sing_Y0 = AbsIdealSheaf[ideal_sheaf_of_singular_locus(Y0)]
+  #I_sing_Y0 = AbsIdealSheaf[simplify(ideal_sheaf_of_singular_locus(Y0))]
+  I_sing_X0 = simplify.(pushforward(inc_Y0).(I_sing_Y0))
+
+  # Prepare a covering which has a permanent weierstrass chart
+  U0 = X0[1][1]
+  disc = numerator(discriminant(generic_fiber(Y)))
+  U = PrincipalOpenSubset(U0, evaluate(disc, gens(OO(U0))[3]))
+  _find_chart(U, default_covering(X0))
+  Cref = Covering(vcat([U], affine_charts(X0)))
+  inherit_gluings!(Cref, X0[1])
+  push!(X0.coverings, Cref)
+  @assert has_decomposition_info(default_covering(X0))
+  inherit_decomposition_info!(X0, Cref)
+  @assert has_decomposition_info(Cref)
+
+  ambient_exceptionals = EffectiveCartierDivisor[]
+  varnames = [:a,:b,:c,:d,:e,:f,:g,:h,:i,:j,:k,:l,:m,:n,:o,:p,:q,:r,:u,:v,:w]
+  projectionsX = BlowupMorphism[]
+  projectionsY = AbsCoveredSchemeMorphism[]
+  count = 0
+
+  @vprint :EllipticSurface 2 "Blowing up Weierstrass model simultaneously in all singular points\n"
+  while true
+    count = count+1
+    @vprint :EllipticSurface 1 "blowup number: $(count)\n"
+    @vprint :EllipticSurface 1 "number of ideal sheaves to be blown up: $(length(I_sing_X0))\n"
+    if length(I_sing_X0)==0
+      # stop if smooth
+      break
+    end
+    # make sure we have a Weierstrass chart kept.
+    if count == 1
+      cov = Cref
+    else
+      cov = simplified_covering(X0)
+    end
+    # take the first ideal sheaf and blow it up
+    J = SimplifiedIdealSheaf(I_sing_X0[1])
+    pr_X1 = blow_up(J, covering=cov, var_name=varnames[1+mod(count, length(varnames))])
+
+    X1 = domain(pr_X1)
+    @vprint :EllipticSurface 1 "$(X1)\n"
+    E1 = exceptional_divisor(pr_X1)
+
+    @vprint :EllipticSurface 2 "computing strict transforms\n"
+    # compute the exceptional divisors
+    ambient_exceptionals = EffectiveCartierDivisor[strict_transform(pr_X1, e) for e in ambient_exceptionals]
+    # move the divisors coming originally from S up to the next chart
+    push!(ambient_exceptionals, E1)
+
+    Y1, inc_Y1, pr_Y1 = strict_transform(pr_X1, inc_Y0)
+    # Speed up the computation of singular loci
+    set_attribute!(Y1, :is_irreducible=> true)
+    set_attribute!(Y1, :is_reduced=>true)
+    set_attribute!(Y1, :is_integral=>true)
+    set_attribute!(Y1, :is_equidimensional=>true)
+
+    # transform the singular loci
+    I_sing_X0 = AbsIdealSheaf[pullback(pr_X1, J) for J in I_sing_X0[2:end]]
+
+    # Add eventual new components
+    @vprint :EllipticSurface 2 "computing singular locus\n"
+    I_sing_new = ideal_sheaf_of_singular_locus(Y1; focus=pullback(inc_Y1, ideal_sheaf(E1)))
+    #I_sing_new = pushforward(inc_Y1, I_sing_new) + ideal_sheaf(E1) # new components only along the exc. set
+    I_sing_new = simplify(pushforward(inc_Y1, I_sing_new))
+
+    @vprint :EllipticSurface 2 "decomposing singular locus\n"
+    !is_one(I_sing_new) && push!(I_sing_X0, I_sing_new)
+
+    push!(projectionsX, pr_X1)
+    push!(projectionsY, pr_Y1)
+    simplify!(Y1)
+
+    # set up for the next iteration
+    Y0 = Y1
+    inc_Y0 = inc_Y1
+    X0 = X1
+    # Speed up the computation of singular loci
+    set_attribute!(Y0, :is_irreducible=> true)
+    set_attribute!(Y0, :is_reduced=>true)
+    set_attribute!(Y0, :is_integral=>true)
+    set_attribute!(Y0, :is_equidimensional=>true)
+    set_attribute!(X0, :is_irreducible=> true)
+    set_attribute!(X0, :is_reduced=>true)
+    set_attribute!(X0, :is_integral=>true)
+  end
+  Y.Y = Y0
+  Y.blowups = projectionsY
+
+  # We need to rewrap the last maps so that the domain is really Y
+  last_pr = pop!(projectionsY)
+  last_pr_wrap = CoveredSchemeMorphism(Y, codomain(last_pr), covering_morphism(last_pr))
+
+  push!(projectionsY, last_pr_wrap)
+  Y.ambient_blowups = projectionsX
+
+  Y.ambient_exceptionals = ambient_exceptionals
+  piY = CompositeCoveredSchemeMorphism(reverse(projectionsY))
+  Y.blowup = piY
+
+  inc_Y0_wrap = CoveredClosedEmbedding(Y, codomain(inc_Y0), covering_morphism(inc_Y0), check=false)
+  Y.inc_Y = inc_Y0_wrap
+
+  set_attribute!(Y, :is_irreducible=> true)
+  set_attribute!(Y, :is_reduced=>true)
+  set_attribute!(Y, :is_integral=>true)
+  return piY
+end
+
+     
+function weierstrass_contraction_iterative(Y::EllipticSurface)
   if isdefined(Y, :blowup)
     return Y.blowup
   end
@@ -630,24 +772,14 @@ function weierstrass_contraction(X::EllipticSurface)
       cov = Crefined
     else
       # the following leads to difficult bugs
-      #=
-      cov0 = simplified_covering(X0)
-      cov1 = _separate_disjoint_components(I_sing_X0, covering=cov0)
-      cov = _one_patch_per_component(cov1, I_sing_X0)
-      push!(X0.coverings, cov)
-      @assert has_decomposition_info(default_covering(X0))
-      inherit_decomposition_info!(X0, cov)
-      @assert has_decomposition_info(cov)
-      =#
       cov = simplified_covering(X0)
-      #inherit_decomposition_info!(cov, X0)
     end
     # take the first singular point and blow it up
     J = SimplifiedIdealSheaf(I_sing_X0[1])
     pr_X1 = blow_up(J, covering=cov, var_name=varnames[1+mod(count, length(varnames))])
 
     # Set the attribute so that the strict_transform does some extra work
-    isomorphism_on_open_subset(pr_X1)
+    #isomorphism_on_open_subset(pr_X1)
 
     X1 = domain(pr_X1)
     @vprint :EllipticSurface 1 "$(X1)\n"
@@ -701,7 +833,7 @@ function weierstrass_contraction(X::EllipticSurface)
   # We need to rewrap the last maps so that the domain is really Y
   last_pr = pop!(projectionsY)
   last_pr_wrap = CoveredSchemeMorphism(Y, codomain(last_pr), covering_morphism(last_pr))
-  set_attribute!(last_pr_wrap, :isomorphism_on_open_subset, get_attribute(last_pr, :isomorphism_on_open_subset))
+  #set_attribute!(last_pr_wrap, :isomorphism_on_open_subset, get_attribute(last_pr, :isomorphism_on_open_subset))
 
   push!(projectionsY, last_pr_wrap)
   Y.ambient_blowups = projectionsX
@@ -847,7 +979,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,18 +1104,19 @@ 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
   
-@attr Any function exceptional_divisors(S::EllipticSurface)
+@attr Vector{<:AbsIdealSheaf} function exceptional_divisors(S::EllipticSurface)
   PP = AbsIdealSheaf[]
   @vprintln :EllipticSurface 2 "computing exceptional divisors"
-  for E in S.ambient_exceptionals
+  # If we have resolution_strategy=:simultaneous, then the following is a non-trivial preprocessing step.
+  ambient_pts = reduce(vcat, maximal_associated_points.(ideal_sheaf.(S.ambient_exceptionals)))
+  for I in ambient_pts
     @vprintln :EllipticSurface 4 "decomposing divisor "
-    mp = maximal_associated_points(ideal_sheaf(pullback(S.inc_Y,E));
-                                   use_decomposition_info=true)
+    mp = maximal_associated_points(pullback(S.inc_Y, I); use_decomposition_info=true)
     append!(PP, mp)
   end 
   @vprintln :EllipticSurface 3 "done"
@@ -1629,7 +1762,10 @@ degree at most ``4`` to Weierstrass form, apply Tate's algorithm and
 return the corresponding relatively minimal elliptic surface 
 as well as the coordinate transformation.
 """
-function elliptic_surface(g::MPolyRingElem, P::Vector{<:RingElem}; minimize=true)
+function elliptic_surface(
+    g::MPolyRingElem, P::Vector{<:RingElem}; 
+    minimize::Bool=true, resolution_strategy::Symbol=:iterative
+  )
   R = parent(g)
   (x, y) = gens(R)
   P = base_ring(R).(P)
diff --git a/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Attributes.jl b/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Attributes.jl
index 5bc7f57aea79..a68601a16d10 100644
--- a/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Attributes.jl
+++ b/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Attributes.jl
@@ -303,18 +303,4 @@ end
 
 morphism_type(X::AbsAffineScheme, Y::AbsAffineScheme) = morphism_type(typeof(X), typeof(Y))
 
-@doc raw"""
-    isomorphism_on_open_subsets(f::AbsAffineSchemeMor)
-
-For a birational morphism ``f : X → Y`` of `AbsAffineScheme`s this 
-returns an isomorphism of affine schemes ``f' : U → V`` which is 
-the restriction of ``f`` to two dense open subsets ``U ⊂ X`` and 
-``V ⊂ Y``.
-"""
-function isomorphism_on_open_subsets(f::AbsAffineSchemeMor)
-  if !has_attribute(f, :iso_on_open_subset)
-    is_birational(f) # Should compute and store the attribute
-  end
-  return get_attribute(f, :iso_on_open_subset)::AbsAffineSchemeMor
-end
 
diff --git a/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Methods.jl b/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Methods.jl
index 3d1eeac5fe2a..43b76180c519 100644
--- a/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Methods.jl
+++ b/src/AlgebraicGeometry/Schemes/AffineSchemes/Morphisms/Methods.jl
@@ -138,7 +138,7 @@ function _induced_map_to_fiber_product(
   W = domain(a)
   @check gens(OO(XxY)) == vcat(pullback(gg).(gens(OO(X))), pullback(ff).(gens(OO(Y)))) "variables must be pullbacks of variables on the factors"
 
-  img_gens = vcat(pullback(a).(gens(OO(X))), pullback(b).(gens(OO(Y))))
+  img_gens = vcat(_images_of_generators(pullback(a)), _images_of_generators(pullback(b)))
   return morphism(W, XxY, img_gens, check=check)
 end
 
@@ -155,8 +155,7 @@ function _induced_map_to_fiber_product(
   XxY = fiber_product[1]
   W = domain(a)
   Y = codomain(b)
-  img_gens = pullback(b).(gens(OO(Y)))
-  return morphism(W, XxY, img_gens, check=check)
+  return morphism(W, XxY, _images_of_generators(pullback(b)), check=check)
 end
 
 function _induced_map_to_fiber_product(
@@ -169,8 +168,7 @@ function _induced_map_to_fiber_product(
   XxY = fiber_product[1]
   X = domain(f)
   W = domain(a)
-  img_gens = pullback(a).(gens(OO(X)))
-  return morphism(W, XxY, img_gens, check=check)
+  return morphism(W, XxY, _images_of_generators(pullback(a)), check=check)
 end
 
 # additional method to remove ambiguity
@@ -184,18 +182,23 @@ function _induced_map_to_fiber_product(
   # XxY is a principal open subset of Y.
   XxY = fiber_product[1]
   Y = domain(g)
-  img_gens = pullback(b).(gens(OO(Y)))
-  return morphism(W, XxY, img_gens, check=check)
+  return morphism(W, XxY, _images_of_generators(pullback(b)), check=check)
 end
 
 ### Some helper functions
+
+_images_of_generators(f::Map) = f.(gens(domain(f)))
+_images_of_generators(f::MPolyAnyMap) = _images(f)
+_images_of_generators(f::MPolyLocalizedRingHom) = _images(restricted_map(f))
+_images_of_generators(f::MPolyQuoLocalizedRingHom) = _images(restricted_map(f))
+
 function _restrict_domain(f::AbsAffineSchemeMor, D::PrincipalOpenSubset; check::Bool=true)
   D === domain(f) && return f
-  !_has_coefficient_map(pullback(f)) && ambient_scheme(D) === domain(f) && return morphism(D, codomain(f), OO(D).(pullback(f).(gens(OO(codomain(f)))); check=false), check=false)
+  !_has_coefficient_map(pullback(f)) && ambient_scheme(D) === domain(f) && return morphism(D, codomain(f), [OO(D)(x; check) for x in _images_of_generators(pullback(f))]; check=check)
 
   @check is_subscheme(D, domain(f)) "domain incompatible"
-  !_has_coefficient_map(pullback(f)) && return morphism(D, codomain(f), [OO(D)(x; check) for x in pullback(f).(gens(OO(codomain(f))))], check=check)
-  return morphism(D, codomain(f), coefficient_map(pullback(f)), [OO(D)(x; check) for x in pullback(f).(gens(OO(codomain(f))))], check=check)
+  !_has_coefficient_map(pullback(f)) && return morphism(D, codomain(f), [OO(D)(x; check) for x in _images_of_generators(pullback(f))], check=check)
+  return morphism(D, codomain(f), coefficient_map(pullback(f)), [OO(D)(x; check) for x in _images_of_generators(pullback(f))], check=check)
 end
 
 function _restrict_domain(f::AbsAffineSchemeMor, D::AbsAffineScheme; check::Bool=true)
@@ -208,12 +211,12 @@ function _restrict_codomain(f::AbsAffineSchemeMor, D::PrincipalOpenSubset; check
   D === codomain(f) && return f
   if ambient_scheme(D) === codomain(f)
     @check is_unit(pullback(f)(complement_equation(D))) "complement equation does not pull back to a unit"
-    !_has_coefficient_map(pullback(f)) && return morphism(domain(f), D, OO(domain(f)).(pullback(f).(gens(OO(codomain(f)))); check=false), check=false)
-    return morphism(domain(f), D, coefficient_map(pullback(f)), [OO(domain(f))(x; check) for x in pullback(f).(gens(OO(codomain(f))))], check=false)
+    !_has_coefficient_map(pullback(f)) && return morphism(domain(f), D, [OO(domain(f))(x; check) for x in _images_of_generators(pullback(f))]; check=check)
+    return morphism(domain(f), D, coefficient_map(pullback(f)), [OO(domain(f))(x; check) for x in _images_of_generators(pullback(f))], check=check)
   end
   @check is_subscheme(D, codomain(f)) "codomain incompatible"
   @check is_subscheme(domain(f), preimage(f, D))
-  !_has_coefficient_map(pullback(f)) && return morphism(domain(f), D, OO(domain(f)).(pullback(f).(gens(OO(codomain(f)))); check=false), check=check)
+  !_has_coefficient_map(pullback(f)) && return morphism(domain(f), D, [OO(domain(f))(x; check) for x in _images_of_generators(pullback(f))]; check=check)
   return morphism(domain(f), D, coefficient_map(pullback(f)), [OO(domain(f))(x; check) for x in pullback(f).(gens(OO(codomain(f))))], check=check)
 end
 
@@ -235,7 +238,7 @@ end
 function _restrict_codomain(f::AbsAffineSchemeMor, D::AbsAffineScheme; check::Bool=true)
   @check is_subscheme(D, codomain(f)) "codomain incompatible"
   @check is_subscheme(domain(f), preimage(f, D)) "new domain is not contained in preimage of codomain"
-  !_has_coefficient_map(pullback(f)) && return morphism(domain(f), D, OO(domain(f)).(pullback(f).(gens(OO(codomain(f)))); check=false), check=check)
+  !_has_coefficient_map(pullback(f)) && return morphism(domain(f), D, [OO(domain(f))(x; check) for x in _images_of_generators(pullback(f))]; check)
   return morphism(domain(f), D, coefficient_map(pullback(f)), [OO(domain(f))(x; check) for x in pullback(f).(gens(OO(codomain(f))))], check=check)
 end
 
@@ -393,8 +396,7 @@ function ==(f::AffineSchemeMorType, g::AffineSchemeMorType) where {AffineSchemeM
   X = domain(f)
   X == domain(g) || return false
   codomain(f) == codomain(g) || return false
-  OO(X).(pullback(f).(gens(ambient_coordinate_ring(codomain(f))))) == OO(X).(pullback(f).(gens(ambient_coordinate_ring(codomain(g))))) || return false
-  return true
+  return all(OO(X)(x) == OO(X)(y) for (x, y) in zip(_images_of_generators(pullback(f)), _images_of_generators(pullback(g))))
 end
 
 
diff --git a/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Attributes.jl b/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Attributes.jl
index 2ad15abb7fac..dd775e171fad 100644
--- a/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Attributes.jl
+++ b/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Attributes.jl
@@ -653,7 +653,7 @@ julia> singular_locus(A3)
 (scheme(1), Hom: scheme(1) -> affine 3-space)
 
 julia> singular_locus(X)
-(scheme(x^2 - y^2 + z^2, z, y, x), Hom: scheme(x^2 - y^2 + z^2, z, y, x) -> scheme(x^2 - y^2 + z^2))
+(scheme(x^2 - y^2 + z^2, 2*x, -2*y, 2*z), Hom: scheme(x^2 - y^2 + z^2, 2*x, -2*y, 2*z) -> scheme(x^2 - y^2 + z^2))
 
 julia> U = complement_of_point_ideal(R, [0,0,0])
 Complement
diff --git a/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Methods.jl b/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Methods.jl
index 51d51eeb2371..9aa319b4510f 100644
--- a/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Methods.jl
+++ b/src/AlgebraicGeometry/Schemes/AffineSchemes/Objects/Methods.jl
@@ -446,11 +446,14 @@ function _change_base_ring(phi::Any,
   )
   R = base_ring(W)
   P, Phi = _change_base_ring(phi, R)
+  @assert _has_coefficient_map(Phi)
   U = inverted_set(W)
   U_red = MPolyPowersOfElement(P, Phi.(denominators(U)))
   W_red, loc_map = localization(P, U_red)
-  res_map = hom(W, W_red, compose(Phi, loc_map), check=false)
-  #@assert _has_coefficient_map(res_map)
+  comp = hom(R, W_red, phi, gens(W_red); check=false)
+  @assert _has_coefficient_map(comp)
+  res_map = hom(W, W_red, comp, check=false)
+  @assert _has_coefficient_map(res_map)
   return W_red, res_map
 end
 
@@ -465,7 +468,8 @@ function _change_base_ring(phi::Any,
   I_red = ideal(W_red, Phi_W.(gens(I)))
   L_red, pr = quo(W_red, I_red)
   res = compose(restricted_map(Phi_W), pr)
+  @assert _has_coefficient_map(res)
   res_map = hom(L, L_red, res, check=false)
-  #@assert _has_coefficient_map(res_map)
+  @assert _has_coefficient_map(res_map)
   return L_red, res_map
 end
diff --git a/src/AlgebraicGeometry/Schemes/CoveredSchemes/Morphisms/Attributes.jl b/src/AlgebraicGeometry/Schemes/CoveredSchemes/Morphisms/Attributes.jl
index 1f662e64e86d..eb6e0aca4f2f 100644
--- a/src/AlgebraicGeometry/Schemes/CoveredSchemes/Morphisms/Attributes.jl
+++ b/src/AlgebraicGeometry/Schemes/CoveredSchemes/Morphisms/Attributes.jl
@@ -72,21 +72,7 @@ domain(f::CoveredSchemeMorphism) = f.X
 codomain(f::CoveredSchemeMorphism) = f.Y
 covering_morphism(f::CoveredSchemeMorphism) = f.f
 
-@doc raw"""
-    isomorphism_on_open_subsets(f::AbsCoveredSchemeMorphism)
-
-For a birational morphism ``f : X → Y`` of `AbsCoveredScheme`s this
-returns an isomorphism of affine schemes ``fᵣₑₛ : U → V`` which is
-the restriction of ``f`` to two dense open subsets ``U ⊂ X`` and
-``V ⊂ Y``.
-"""
-function isomorphism_on_open_subsets(f::AbsCoveredSchemeMorphism)
-  if !has_attribute(f, :iso_on_open_subset)
-    is_birational(f) # Should compute and store the attribute
-  end
-  return get_attribute(f, :iso_on_open_subset)::AbsAffineSchemeMor
-end
-
 @attr AbsCoveredSchemeMorphism function inverse(f::AbsCoveredSchemeMorphism)
   error("method not implemented")
 end
+
diff --git a/src/AlgebraicGeometry/Schemes/Covering/Morphisms/Methods.jl b/src/AlgebraicGeometry/Schemes/Covering/Morphisms/Methods.jl
index a51658e992b8..daf8d107e829 100644
--- a/src/AlgebraicGeometry/Schemes/Covering/Morphisms/Methods.jl
+++ b/src/AlgebraicGeometry/Schemes/Covering/Morphisms/Methods.jl
@@ -336,8 +336,9 @@ function fiber_product(f::CoveringMorphism, g::CoveringMorphism)
       pre_glue_double_to_simple = induced_map_to_fiber_product(
            compose(f_res_UU, compose(double_to_simple_U, inc_U)),
            compose(g_res_VV, compose(double_to_simple_V, inc_V)),
-           f[U], g[V],
-           fiber_product=(UxV, to_U, to_V)
+           f[U], g[V];
+           fiber_product=(UxV, to_U, to_V),
+           check=false
           )
       double_to_simple = restrict(pre_glue_double_to_simple, domain(pre_glue_double_to_simple), 
                                   U_UUxV_VV; check=false
@@ -345,8 +346,9 @@ function fiber_product(f::CoveringMorphism, g::CoveringMorphism)
       pre_glue_simple_to_double = induced_map_to_fiber_product(
            compose(f_res_U, compose(simple_to_double_U, inc_UU)),
            compose(g_res_V, compose(simple_to_double_V, inc_VV)),
-           f[UU], g[VV],
-           fiber_product=(UUxVV, to_UU, to_VV)
+           f[UU], g[VV];
+           fiber_product=(UUxVV, to_UU, to_VV),
+           check=false
           )
       simple_to_double = restrict(pre_glue_simple_to_double, domain(pre_glue_simple_to_double), 
                                   UU_UxVV_V; check=false
diff --git a/src/AlgebraicGeometry/Schemes/Covering/Objects/Methods.jl b/src/AlgebraicGeometry/Schemes/Covering/Objects/Methods.jl
index 89800660a684..d1ab09fcf1d5 100644
--- a/src/AlgebraicGeometry/Schemes/Covering/Objects/Methods.jl
+++ b/src/AlgebraicGeometry/Schemes/Covering/Objects/Methods.jl
@@ -369,6 +369,18 @@ function has_ancestor_in(L::Vector, U::AbsAffineScheme)
 end
 
 function is_refinement(D::Covering, C::Covering)
+  # catch the case that `C` is a simplification of `D`
+  if all(x isa SimplifiedAffineScheme for x in patches(C)) && all(y->any(x->original(x)===y, patches(C)), patches(D))
+    map_dict = IdDict{AbsAffineScheme, AbsAffineSchemeMor}()
+    for U in patches(C)
+      f, g = identification_maps(U)
+      V = domain(g)
+      map_dict[V] = g
+    end
+    return true, CoveringMorphism(D, C, map_dict, check=false)
+    # TODO: Catch the more complicated case where `D` is a proper 
+    # refinement of the original covering for `C`.
+  end
   if !all(x->has_ancestor(u->any(y->(u===y), patches(C)), x), patches(D))
     return false, nothing
   end
diff --git a/src/AlgebraicGeometry/Schemes/Divisors/base_change.jl b/src/AlgebraicGeometry/Schemes/Divisors/base_change.jl
new file mode 100644
index 000000000000..b57f1184f3f1
--- /dev/null
+++ b/src/AlgebraicGeometry/Schemes/Divisors/base_change.jl
@@ -0,0 +1,19 @@
+function base_change(
+    phi::Any, C::AbsAlgebraicCycle; 
+    scheme_base_change::Map=base_change(phi, scheme(C)[2])
+  )
+  X = scheme(C)
+  @assert codomain(scheme_base_change) === X
+  Y = domain(scheme_base_change)
+  kk = coefficient_ring(C)
+  return AlgebraicCycle(Y, kk, IdDict{AbsIdealSheaf, elem_type(kk)}(pullback(scheme_base_change, I)=>c for (I, c) in coefficient_dict(C)); check=false)
+end
+
+function base_change(
+    phi::Any, D::WeilDivisor;
+    scheme_base_change::Map=base_change(phi, scheme(D))[2]
+  )
+  und = base_change(phi, underlying_cycle(D); scheme_base_change)
+  return WeilDivisor(und; check=false)
+end
+
diff --git a/src/AlgebraicGeometry/Schemes/FunctionField/FunctionFields.jl b/src/AlgebraicGeometry/Schemes/FunctionField/FunctionFields.jl
index 4c1ab2106581..91c018e0b6e0 100644
--- a/src/AlgebraicGeometry/Schemes/FunctionField/FunctionFields.jl
+++ b/src/AlgebraicGeometry/Schemes/FunctionField/FunctionFields.jl
@@ -496,22 +496,6 @@ function is_regular(f::VarietyFunctionFieldElem, W::AffineSchemeOpenSubscheme)
   return all(U->is_regular(f, U), affine_patches(W))
 end
 
-function pushforward(f::AbsCoveredSchemeMorphism, a::VarietyFunctionFieldElem)
-  X = domain(f)
-  Y = codomain(f)
-  parent(a) === function_field(X) || error("element does not belong to the correct ring")
-  has_attribute(f, :isomorphism_on_open_subset) || error("need an isomorphism on some open subset")
-  f_res = isomorphism_on_open_subset(f)
-  U = domain(f_res)
-  V = codomain(f_res)
-  aa = a[ambient_scheme(U)]
-  f_res_inv = inverse(f_res)
-  num = pullback(f_res_inv)(numerator(aa))
-  den = pullback(f_res_inv)(denominator(aa))
-  #bb = fraction(num)//fraction(den)
-  return function_field(Y)(lifted_numerator(num)*lifted_denominator(den), lifted_numerator(den)*lifted_denominator(num))
-end
-
 function pullback(f::AbsCoveredSchemeMorphism, a::VarietyFunctionFieldElem)
   fcov = covering_morphism(f)
   KK = parent(a)
diff --git a/src/AlgebraicGeometry/Schemes/Gluing/Methods.jl b/src/AlgebraicGeometry/Schemes/Gluing/Methods.jl
index df8966835392..7024fc902e00 100644
--- a/src/AlgebraicGeometry/Schemes/Gluing/Methods.jl
+++ b/src/AlgebraicGeometry/Schemes/Gluing/Methods.jl
@@ -230,6 +230,8 @@ function base_change(phi::Any, G::AbsGluing;
     patch_change1::AbsAffineSchemeMor=base_change(phi, gluing_domains(G)[1])[2], # The base change morphism for 
     patch_change2::AbsAffineSchemeMor=base_change(phi, gluing_domains(G)[2])[2]  # the two gluing patches
   )
+  @assert _has_coefficient_map(pullback(patch_change1))
+  @assert _has_coefficient_map(pullback(patch_change2))
   gd = BaseChangeGluingData{typeof(phi), typeof(G)}(phi, G, patch_change1, patch_change2)
   return LazyGluing(domain(patch_change1), domain(patch_change2), gd)
 end
diff --git a/src/AlgebraicGeometry/Schemes/PrincipalOpenSubset/Objects/Methods.jl b/src/AlgebraicGeometry/Schemes/PrincipalOpenSubset/Objects/Methods.jl
index 7668ea01164f..55c5964647a2 100644
--- a/src/AlgebraicGeometry/Schemes/PrincipalOpenSubset/Objects/Methods.jl
+++ b/src/AlgebraicGeometry/Schemes/PrincipalOpenSubset/Objects/Methods.jl
@@ -79,6 +79,7 @@ function base_change(phi::Any, U::PrincipalOpenSubset;
   pbf = pullback(ambient_map)
   h = pbf(complement_equation(U))
   UU = PrincipalOpenSubset(Y, h)
+  @assert _has_coefficient_map(pullback(ambient_map))
   res_map = restrict(ambient_map, UU, U, check=false)
   #@assert _has_coefficient_map(pullback(res_map))
   return UU, res_map
diff --git a/src/AlgebraicGeometry/Schemes/main.jl b/src/AlgebraicGeometry/Schemes/main.jl
index 6f2eeca514ac..4777e8b57e97 100644
--- a/src/AlgebraicGeometry/Schemes/main.jl
+++ b/src/AlgebraicGeometry/Schemes/main.jl
@@ -190,3 +190,4 @@ include("FunctionField/FunctionFields.jl")
 include("Divisors/AlgebraicCycles.jl")
 include("Divisors/WeilDivisor.jl")
 include("Divisors/CartierDivisor.jl")
+include("Divisors/base_change.jl")
diff --git a/src/Rings/MPolyMap/MPolyAnyMap.jl b/src/Rings/MPolyMap/MPolyAnyMap.jl
index 8394cbb8f749..aba25f188bfe 100644
--- a/src/Rings/MPolyMap/MPolyAnyMap.jl
+++ b/src/Rings/MPolyMap/MPolyAnyMap.jl
@@ -29,51 +29,7 @@
 # - MPolyQuoRing
 # - MPolyDecRing
 
-const _DomainTypes = Union{MPolyRing, MPolyQuoRing}
-
-@attributes mutable struct MPolyAnyMap{
-    D <: _DomainTypes,
-    C <: NCRing,
-    U,
-    V} <: Map{D, C, Map, MPolyAnyMap}
-
-  domain::D
-  codomain::C
-  coeff_map::U
-  img_gens::Vector{V}
-  temp_ring           # temporary ring used when evaluating maps
-  variable_indices::Vector{Int} # a table where the i-th entry contains the 
-                                # index of the variable where it is mapped to
-                                # in case the mapping takes such a particularly
-                                # simple form.
-
-  function MPolyAnyMap{D, C, U, V}(domain::D,
-                                   codomain::C,
-                                   coeff_map::U,
-                                   img_gens::Vector{V};
-                                   check_for_mapping_of_vars::Bool=true
-                                  ) where {D, C, U, V}
-    @assert V === elem_type(C)
-    for g in img_gens
-      @assert parent(g) === codomain "elements does not have the correct parent"
-    end
-    result = new{D, C, U, V}(domain, codomain, coeff_map, img_gens)
-    # If it ever turned out that doing the checks within the following if-block
-    # is a bottleneck, consider passing on the `check_for_mapping_of_vars` kw 
-    # argument to the outer constructors or make the outer constructors 
-    # call the inner one with this argument set to `false`. This way the check 
-    # can safely be disabled.
-    if check_for_mapping_of_vars && all(_is_gen, img_gens) && _allunique(img_gens)
-      gens_codomain = gens(codomain)
-      result.variable_indices = [findfirst(==(x), gens_codomain) for x in img_gens]
-    end
-    return result
-  end
-end
-
-function MPolyAnyMap(d::D, c::C, cm::U, ig::Vector{V}) where {D, C, U, V}
-  return MPolyAnyMap{D, C, U, V}(d, c, cm, ig)
-end
+### See `Types.jl` for the declaration of the type.
 
 ################################################################################
 #
@@ -255,6 +211,32 @@ function compose(F::MPolyAnyMap{D, C, S}, G::MPolyAnyMap{C, E, U}) where {D, C,
   end
 end
 
+# No coefficient maps in the second argument
+function compose(F::MPolyAnyMap{D, C, S}, G::MPolyAnyMap{C, E, Nothing}) where {D, C, E, S <: Map}
+  @req codomain(F) === domain(G) "Incompatible (co)domain in composition"
+  f = coefficient_map(F)
+  if typeof(codomain(f)) === typeof(coefficient_ring(domain(G)))
+    return hom(domain(F), codomain(G), f, G.(_images(F)), check=false)
+  elseif typeof(codomain(f)) === typeof(domain(G))
+    new_coeff_map = compose(f, G)
+    return hom(domain(F), codomain(G), new_coeff_map, G.(_images(F)), check=false)
+  else
+    return Generic.CompositeMap(F, G)
+  end
+end
+
+# No coefficient maps in the first argument
+function compose(F::MPolyAnyMap{D, C, Nothing}, G::MPolyAnyMap{C, E, S}) where {D, C, E, S <: Map}
+  @req codomain(F) === domain(G) "Incompatible (co)domain in composition"
+  g = coefficient_map(G)
+  if domain(g) === coefficient_ring(domain(F))
+    return hom(domain(F), codomain(G), g, G.(_images(F)), check=false)
+  else
+    new_coeff_map = MapFromFunc(coefficient_ring(domain(F)), codomain(g), x->g(domain(g)(x)))
+    return hom(domain(F), codomain(G), new_coeff_map, G.(_images(F)), check=false)
+  end
+end
+
 # No coefficient maps in both maps
 function compose(F::MPolyAnyMap{D, C, Nothing}, G::MPolyAnyMap{C, E, Nothing}) where {D, C, E}
   @req codomain(F) === domain(G) "Incompatible (co)domain in composition"
diff --git a/src/Rings/MPolyMap/Types.jl b/src/Rings/MPolyMap/Types.jl
new file mode 100644
index 000000000000..9f22448a96e9
--- /dev/null
+++ b/src/Rings/MPolyMap/Types.jl
@@ -0,0 +1,47 @@
+### Types for maps from polynomial rings
+const _DomainTypes = Union{MPolyRing, MPolyQuoRing}
+
+@attributes mutable struct MPolyAnyMap{
+    D <: _DomainTypes,
+    C <: NCRing,
+    U,
+    V} <: Map{D, C, Map, MPolyAnyMap}
+
+  domain::D
+  codomain::C
+  coeff_map::U
+  img_gens::Vector{V}
+  temp_ring           # temporary ring used when evaluating maps
+  variable_indices::Vector{Int} # a table where the i-th entry contains the 
+                                # index of the variable where it is mapped to
+                                # in case the mapping takes such a particularly
+                                # simple form.
+
+  function MPolyAnyMap{D, C, U, V}(domain::D,
+                                   codomain::C,
+                                   coeff_map::U,
+                                   img_gens::Vector{V};
+                                   check_for_mapping_of_vars::Bool=true
+                                  ) where {D, C, U, V}
+    @assert V === elem_type(C)
+    for g in img_gens
+      @assert parent(g) === codomain "elements does not have the correct parent"
+    end
+    result = new{D, C, U, V}(domain, codomain, coeff_map, img_gens)
+    # If it ever turned out that doing the checks within the following if-block
+    # is a bottleneck, consider passing on the `check_for_mapping_of_vars` kw 
+    # argument to the outer constructors or make the outer constructors 
+    # call the inner one with this argument set to `false`. This way the check 
+    # can safely be disabled.
+    if check_for_mapping_of_vars && all(_is_gen, img_gens) && _allunique(img_gens)
+      gens_codomain = gens(codomain)
+      result.variable_indices = [findfirst(==(x), gens_codomain) for x in img_gens]
+    end
+    return result
+  end
+end
+
+function MPolyAnyMap(d::D, c::C, cm::U, ig::Vector{V}) where {D, C, U, V}
+  return MPolyAnyMap{D, C, U, V}(d, c, cm, ig)
+end
+
diff --git a/src/Rings/MPolyQuo.jl b/src/Rings/MPolyQuo.jl
index 7c37f1c1ae2b..79f407830116 100644
--- a/src/Rings/MPolyQuo.jl
+++ b/src/Rings/MPolyQuo.jl
@@ -546,10 +546,12 @@ end
   return is_prime(saturated_ideal(I))
 end
 
-@attr MPolyQuoIdeal function radical(I::MPolyQuoIdeal)
-  R = base_ring(I)
-  J = saturated_ideal(I)
-  return ideal(R, [g for g in R.(gens(radical(J))) if !iszero(g)])
+function radical(I::MPolyQuoIdeal; eliminate_variables::Bool=true)
+  get_attribute!(I, :radical) do
+    R = base_ring(I)
+    J = saturated_ideal(I)
+    return ideal(R, [g for g in R.(gens(radical(J; eliminate_variables))) if !iszero(g)])
+  end::typeof(I)
 end
 
 # The following is to streamline the programmer's
diff --git a/src/Rings/Rings.jl b/src/Rings/Rings.jl
index 7ec21f35cf30..2a0d9ad474f8 100644
--- a/src/Rings/Rings.jl
+++ b/src/Rings/Rings.jl
@@ -12,6 +12,7 @@ include("FractionalIdeal.jl")
 
 include("special_ideals.jl")
 
+include("MPolyMap/Types.jl")
 include("MPolyMap/MPolyAnyMap.jl")
 include("MPolyMap/MPolyRing.jl")
 include("MPolyMap/MPolyQuo.jl")
diff --git a/src/Rings/mpoly-ideals.jl b/src/Rings/mpoly-ideals.jl
index 99d510bc1c6a..27b717b3d013 100644
--- a/src/Rings/mpoly-ideals.jl
+++ b/src/Rings/mpoly-ideals.jl
@@ -463,18 +463,17 @@ end
 
 #######################################################
 @doc raw"""
-    radical(I::MPolyIdeal)
+    radical(I::MPolyIdeal; eliminate_variables::Bool=true)
 
 Return the radical of `I`.
 
 # Implemented Algorithms
 
-If the base ring of `I` is a polynomial
-ring over a field, a combination of the algorithms of Krick and Logar
-(with modifications by Laplagne) and Kemper is used. For polynomial
-rings over the integers, the algorithm proceeds as suggested by
-Pfister, Sadiq, and Steidel. See [KL91](@cite),
-[Kem02](@cite), and [PSS11](@cite).
+  * If the base ring of `I` is a polynomial ring over a the field of rational numbers or a finite field, a combination of the algorithms of Krick and Logar (with modifications by Laplagne) and Kemper is used. 
+  * If the `base_ring` of `I` is a number field, then we first expand the minimal polynomials to reduce to a computation over the rationals.
+  * For polynomial rings over the integers, the algorithm proceeds as suggested by Pfister, Sadiq, and Steidel. See [KL91](@cite), [Kem02](@cite), and [PSS11](@cite).
+
+When `eliminate_variables` is set to `true` a preprocessing heuristic is applied in order to find variables which can be eliminated using `I`.
 
 # Examples
 ```jldoctest
@@ -570,9 +569,56 @@ Ideal generated by
   39*a*c
   6*a*b*d
   3*a*b*c
+  1326*a^2*d^5
+  1989*a^2*c^5
+  102*b^4*d^5
+  153*b^4*c^5
+  663*a^2*c^5*d^5
+  51*b^4*c^5*d^5
+  78*a^2*d^15
+  117*a^2*c^15
+  78*a^15*d^5
+  117*a^15*c^5
+  6*a^2*b^4*d^15
+  9*a^2*b^4*c^15
+  39*a^2*c^5*d^15
+  39*a^2*c^15*d^5
+  6*a^2*b^15*d^5
+  9*a^2*b^15*c^5
+  6*a^15*b^4*d^5
+  9*a^15*b^4*c^5
+  39*a^15*c^5*d^5
+  3*a^2*b^4*c^5*d^15
+  3*a^2*b^4*c^15*d^5
+  3*a^2*b^15*c^5*d^5
+  3*a^15*b^4*c^5*d^5
 ```
 """
-@attr T function radical(I::T) where {T <: MPolyIdeal}
+function radical(
+    I::MPolyIdeal{T}; 
+    eliminate_variables::Bool=true
+  ) where {T <: MPolyRingElem}
+  get_attribute(I, :radical) do
+    if eliminate_variables
+      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 = ideal(R, small_generating_set(pre_res + I))
+      set_attribute!(res, :is_radical=>true)
+      return res
+    end
+    return _compute_radical(I)
+  end::MPolyIdeal{T}
+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))
@@ -595,8 +641,10 @@ end
 # Rerouting via expansion of the coefficient field
 function radical(
     I::MPolyIdeal{T};
-    factor_generators::Bool=true
+    factor_generators::Bool=true, 
+    eliminate_variables::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)
@@ -619,7 +667,7 @@ function radical(
     end
     R_flat, iso, iso_inv = _expand_coefficient_field_to_QQ(R)
     I_flat = ideal(R_flat, iso_inv.(gens(J)))
-    I_flat_rad = radical(I_flat)
+    I_flat_rad = radical(I_flat; eliminate_variables)
     Irad = iso(I_flat_rad)
     set_attribute!(Irad, :is_radical => true)
     @hassert :IdealSheaves 2 !is_one(Irad)
@@ -629,7 +677,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 +690,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 +1202,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 +1269,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 +1295,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 +2309,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(base_ring(I), c)
+end
+
diff --git a/src/Rings/mpoly-localizations.jl b/src/Rings/mpoly-localizations.jl
index 8a1d06b50524..176fd337162a 100644
--- a/src/Rings/mpoly-localizations.jl
+++ b/src/Rings/mpoly-localizations.jl
@@ -143,7 +143,7 @@ function simplify_light(S::MPolyPowersOfElement)
 
   # Try to factor the denominators if this is deemed possible
   for a in new_denom
-    if total_degree(a) > 50 || length(a) > 10000
+    if total_degree(a) > 50 || length(a) > 10000 || coefficient_ring(R) isa Hecke.RelSimpleNumField{AbsSimpleNumFieldElem}
       push!(fac_denom, a)
     else
       for (b, _) in factor(a)
@@ -1146,7 +1146,7 @@ function *(a::MPolyLocRingElem{BRT, BRET, RT, RET, MST}, b::RET) where {BRT, BRE
 end
 
 function *(a::BRET, b::MPolyLocRingElem{BRT, BRET, RT, RET, MST}) where {BRT, BRET <: RingElem, RT, RET, MST}
-  return (parent(b))(a*numerator(b), denominator(b), check=false)
+  return (parent(b))(a*fraction(b), check=false)
 end
 
 function *(a::MPolyLocRingElem{BRT, BRET, RT, RET, MST}, b::BRET) where {BRT, BRET <: RingElem, RT, RET, MST}
diff --git a/test/AlgebraicGeometry/Schemes/BlowupMorphism.jl b/test/AlgebraicGeometry/Schemes/BlowupMorphism.jl
index 3a8a74cb596f..c8194d4eb5d6 100644
--- a/test/AlgebraicGeometry/Schemes/BlowupMorphism.jl
+++ b/test/AlgebraicGeometry/Schemes/BlowupMorphism.jl
@@ -51,8 +51,6 @@ end
   @test compose(f, h) == identity_map(U)
   V = codomain(f)
   @test compose(h, f) == identity_map(V)
-  g = Oscar.isomorphism_on_open_subset(p)
-  @test is_isomorphism(g)
   KY = function_field(Y)
   KX = function_field(X)
   y, z = gens(ambient_coordinate_ring(first(affine_charts(Y))))
@@ -64,7 +62,6 @@ end
   C = domain(inc_C)
   C_up, inc_C_up, p_res = strict_transform(p, inc_C)
 
-  @test is_isomorphism(Oscar.isomorphism_on_open_subset(p_res))
 
   KC_up = function_field(C_up)
   KC = function_field(C)
@@ -93,29 +90,21 @@ end
   xx0 = KX0(x)
   yy0 = KX0(y)
   zz0 = KX0(z)
-  # Need to call this once so that the attribute is set
-  Oscar.isomorphism_on_open_subset(pr1)
   X1, inc1, pr1_res = strict_transform(pr1, inc)
   xx1 = pullback(pr1_res)(xx0)
-  @test xx0 == pushforward(pr1_res, xx1)
   yy1 = pullback(pr1_res)(yy0)
-  @test yy0 == pushforward(pr1_res, yy1)
   zz1 = pullback(pr1_res)(zz0)
-  @test zz0 == pushforward(pr1_res, zz1)
 
   I_sing = radical(pushforward(inc1, Oscar.ideal_sheaf_of_singular_locus(X1)))
 
   pr2 = blow_up(I_sing)
   @test scheme(I_sing) === domain(pr1)
   @test codomain(pr2) === domain(pr1)
-  Oscar.isomorphism_on_open_subset(pr2)
   Y2 = domain(pr2)
   X2, inc2, pr2_res = strict_transform(pr2, inc1)
 
   pr_res = Oscar.composite_map(pr2_res, pr1_res)
   pr = Oscar.compose(pr2, pr1)
-  @test pushforward(pr, pullback(pr, yy)^2) == yy^2
-  @test pushforward(pr_res, pullback(pr_res, yy0)^2) == yy0^2
 
   pr_inc = compose(compose(inc2, pr2), pr1)
   @test Oscar.underlying_morphism(pr1) isa CoveredSchemeMorphism
diff --git a/test/AlgebraicGeometry/Schemes/CoveredProjectiveSchemes.jl b/test/AlgebraicGeometry/Schemes/CoveredProjectiveSchemes.jl
index 9c77ab5bcf5f..b7464fa06d69 100644
--- a/test/AlgebraicGeometry/Schemes/CoveredProjectiveSchemes.jl
+++ b/test/AlgebraicGeometry/Schemes/CoveredProjectiveSchemes.jl
@@ -57,7 +57,7 @@ end
   @test !is_empty(U[1])
   @test dim(U[1]) == 1
   @test dim(U[2]) == 0
-  B2 = blow_up(image_ideal(inc), var_name="t") # Use a different letter for the homogeneous variables in the 2nd blowup
+  B2 = blow_up(simplify(image_ideal(inc)), var_name="t") # Use a different letter for the homogeneous variables in the 2nd blowup
   Z = domain(B2)
   @test is_smooth(Z)
 end
diff --git a/test/Rings/MPolyAnyMap/MPolyRing.jl b/test/Rings/MPolyAnyMap/MPolyRing.jl
index 99473253e188..b1913963a429 100644
--- a/test/Rings/MPolyAnyMap/MPolyRing.jl
+++ b/test/Rings/MPolyAnyMap/MPolyRing.jl
@@ -124,7 +124,7 @@
     Kxz, (z1, z2) = Kx[:z1, :z2];
     f = hom(Kxz, Kxz, hom(Kx, Kxz, [z1, z2]), [z2, z1])
     g = hom(Kxz, Kxz, hom(Kx, Kx, [y, x]), [z1 + 1, z2 + 1])
-    fg = @inferred f * g
+    fg = f * g
     @test fg(x) == g(f(x))
     @test fg(y) == g(f(y))
     @test fg(z1) == g(f(z1))
@@ -142,7 +142,7 @@
   Qix, (x, y) = Qi[:x, :y]
   f = hom(Qix, Qix, hom(Qi, Qi, -i), [x^2, y^2])
   g = hom(Qix, Qix, hom(Qi, Qi, -i), [x + 1, y + 1])
-  fg = @inferred f * g
+  fg = f * g
   @test fg(i) == g(f(i))
   @test fg(x) == g(f(x))
   @test fg(y) == g(f(y))
@@ -151,7 +151,7 @@
   Qix, (x, y) = Qi[:x, :y]
   f = hom(Qix, Qix, z -> z + 1, [x^2, y^2])
   g = hom(Qix, Qix, z -> z^2, [x + 1, y + 1])
-  fg = @inferred f * g
+  fg = f * g
   @test fg(i) == g(f(i))
   @test fg(x) == g(f(x))
   @test fg(y) == g(f(y))