Fix blowups along rays and singular cones

docs/src/PolyhedralGeometry/intro.md

@@ -85,6 +85,7 @@ accessed as follows:
 
 ```@docs
 primitive_generator(r::AbstractVector{T}) where T<:RationalUnion
+primitive_generator_with_scaling_factor(r::AbstractVector{T}) where T<:RationalUnion
 ```
 
 `AbstractCollection[PointVector]` can be given as:

experimental/Schemes/src/ToricBlowups/constructors.jl

@@ -205,8 +205,15 @@ end
 @doc raw"""
     blow_up(v::NormalToricVariety, n::Int; coordinate_name::String = "e")
 
-Blow up the toric variety by subdividing the n-th cone in the list
-of *all* cones of the fan of `v`. This cone need not be maximal.
+Blow up the toric variety $v$ with polyhedral fan $\Sigma$ by star
+subdivision along the barycenter of the $n$-th cone $\sigma$ in the list
+of all the cones of $\Sigma$.
+We remind that the barycenter of a nonzero cone is the primitive
+generator of the sum of the primitive generators of the extremal rays of
+the cone (Exercise 11.1.10 in [CLS11](@cite)).
+In the case all the cones of $\Sigma$ containing $\sigma$ are smooth,
+this coincides with the star subdivision of $\Sigma$ relative to
+$\sigma$ (Definition 3.3.17 of [CLS11](@cite)).
 This function returns the corresponding morphism.
 
 By default, we pick "e" as the name of the homogeneous coordinate for
@@ -238,9 +245,28 @@ function blow_up(v::NormalToricVarietyType, n::Int; coordinate_name::Union{Strin
   gens_S = gens(cox_ring(v))
   center_unnormalized = ideal_sheaf(v, ideal([gens_S[i] for i in 1:number_of_rays(v) if cones(v)[n,i]]))
   blown_up_variety = normal_toric_variety(star_subdivision(v, n))
-  rays_of_variety = matrix(ZZ, rays(v))
-  exceptional_ray = vec(sum([rays_of_variety[i, :] for i in 1:number_of_rays(v) if cones(v)[n, i]]))
-  return ToricBlowupMorphism(v, blown_up_variety, coordinate_name, exceptional_ray, exceptional_ray, center_unnormalized)
+
+  # minimal supercone coordinates
+  coords = zeros(QQ, n_rays(v))
+  for i in 1:number_of_rays(v)
+    cones(v)[n, i] && (coords[i] = QQ(1))
+  end
+  exceptional_ray_scaled = standard_coordinates(polyhedral_fan(v), coords)
+  exceptional_ray, scaling_factor = primitive_generator_with_scaling_factor(
+    exceptional_ray_scaled
+  )
+  coords = scaling_factor * coords
+
+  phi = ToricBlowupMorphism(
+    v,
+    blown_up_variety,
+    coordinate_name,
+    exceptional_ray,
+    exceptional_ray,
+    center_unnormalized
+  )
+  set_attribute!(phi, :minimal_supercone_coordinates_of_exceptional_ray, coords)
+  return phi
 end
 
 

src/PolyhedralGeometry/PolyhedralFan/properties.jl

@@ -284,7 +284,7 @@ end
 @doc raw"""
     primitive_generator(r::AbstractVector{T}) where T<:RationalUnion -> Vector{ZZRingElem}
 
-Returns the primitive generator, also called minimal generator, of a ray.
+Given a vector `r` which is not a zero vector, returns the primitive generator of `r`, meaning the integer point on the ray generated by `r` that is closest to the origin.
 
 # Examples
 ```jldoctest
@@ -303,10 +303,38 @@ julia> primitive_generator(r)
 ```
 """
 function primitive_generator(r::AbstractVector{T}) where {T<:RationalUnion}
-  result = numerator.(lcm([denominator(i) for i in r]) * r)
+  return Vector{ZZRingElem}(first(primitive_generator_with_scaling_factor(r)))
+end
+
+@doc raw"""
+    primitive_generator_with_scaling_factor(r::AbstractVector{T}) where T<:RationalUnion -> Tuple{Vector{ZZRingElem}, QQFieldElem}
+
+Given a vector `r` which is not a zero vector, returns the primitive generator `s` of `r` together with the rational number `q` such that `qr = s`.
+
+# Examples
+```jldoctest
+julia> PF = polyhedral_fan(IncidenceMatrix([[1, 2]]), [[2, -1], [0, 1]])
+Polyhedral fan in ambient dimension 2
+
+julia> r = rays(PF)[1]
+2-element RayVector{QQFieldElem}:
+ 1
+ -1//2
+
+julia> primitive_generator_with_scaling_factor(r)
+(ZZRingElem[2, -1], 2)
+```
+"""
+function primitive_generator_with_scaling_factor(
+  r::AbstractVector{T}
+) where {T<:RationalUnion}
+  first_scaling_factor = ZZ(lcm(denominator.(r)))
+  result = ZZ.(first_scaling_factor * r)
   g = gcd(result)
+  @req g > 0 "The vector `r` cannot be a zero vector"
+  scaling_factor = QQ(first_scaling_factor, g)
   result = map(x -> div(x, g), result)
-  return Vector{ZZRingElem}(result)
+  return Tuple{Vector{ZZRingElem},QQFieldElem}((result, scaling_factor))
 end
 
 @doc raw"""

src/PolyhedralGeometry/PolyhedralFan/standard_constructions.jl

@@ -252,7 +252,6 @@ julia> ray_indices(maximal_cones(star))
 function star_subdivision(Sigma::_FanLikeType, n::Int)
   cones_Sigma = cones(Sigma)
   tau = Polymake.row(cones_Sigma, n)
-  @req length(tau) > 1 "Cannot subdivide cone $n as it is generated by a single ray"
   R = matrix(ZZ, rays(Sigma))
   exceptional_ray = vec(sum([R[i, :] for i in tau]))
   exceptional_ray = exceptional_ray ./ gcd(exceptional_ray)

src/exports.jl

@@ -1330,6 +1330,7 @@ export prime_ideal
 export prime_of_pgroup, has_prime_of_pgroup, set_prime_of_pgroup
 export primitive_collections
 export primitive_generator
+export primitive_generator_with_scaling_factor
 export primitive_group
 export primitive_group_identification, has_primitive_group_identification
 export primorial

test/AlgebraicGeometry/ToricVarieties/toric_blowups.jl

@@ -78,5 +78,41 @@
     @test bl2 isa Oscar.ToricBlowupMorphism
     @test center_unnormalized(bl2) == II
   end
-
+
+  @testset "Toric blowups along singular cones" begin
+    # 1/2(1, 1) quotient singularity, blowup along the maximal cone
+    ray_generators = [[2, -1], [0, 1]]
+    max_cones = IncidenceMatrix([[1, 2]])
+    X = normal_toric_variety(max_cones, ray_generators)
+    f = blow_up(X, 1)
+    @test ray_vector(QQFieldElem, [1, 0]) in rays(domain(f))
+    @test (
+      minimal_supercone_coordinates_of_exceptional_ray(f)
+      ==
+      QQFieldElem[1//2, 1//2]
+    )
+
+    # Now blowing up along an existing ray
+    g = blow_up(X, 2)
+    @test n_rays(domain(g)) == 2
+    @test n_cones(domain(g)) == 3
+
+    # Quadratic cone, blowup along maximal cone
+    ray_generators = [[0, 0, 1], [1, 0, 1], [1, 1, 1], [0, 1, 1]]
+    max_cones = IncidenceMatrix([[1, 2, 3, 4]])
+    PF = polyhedral_fan(max_cones, ray_generators)
+    X = normal_toric_variety(PF)
+    f = blow_up(X, 1)
+    @test ray_vector(QQFieldElem, [1, 1, 2]) in rays(domain(f))
+    @test (
+      minimal_supercone_coordinates_of_exceptional_ray(f)
+      ==
+      QQFieldElem[1//2, 1//2, 1//2, 1//2]
+    )
+
+    # Now blowing up along an existing ray
+    g = blow_up(X, 6)
+    @test n_rays(domain(g)) == 4
+    @test n_cones(domain(g)) == 11
+  end
 end