From 0c4f03b45e9bba6970dcff0c6dd81c67206b1361 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Thu, 6 Jun 2024 11:06:14 +0100 Subject: [PATCH 01/14] TropicalGeometry: beginning wrapping Singular's tropicalVariety --- src/TropicalGeometry/variety.jl | 44 +++++++++++++++++++++++++++++++++ 1 file changed, 44 insertions(+) diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl index ec41a4cd66b7..fc468b01779e 100644 --- a/src/TropicalGeometry/variety.jl +++ b/src/TropicalGeometry/variety.jl @@ -320,6 +320,7 @@ end ################################################################################ # # Zero-dimensional ideals +# (p-adic valuation only) TODO: add trivial valuation # ################################################################################ @@ -414,6 +415,49 @@ function simultaneous_diagonalization(v::Vector{<:MatElem{T}}) where T <: Hecke. end +################################################################################ +# +# positive-dimensional ideals +# (bandaid wrapeer of Singular's tropical variety) +# +################################################################################ + +function singular_fan_string_to_oscar_fan(input_string::String) + + # Extracting the RAYS, ORTH_LINEALITY_SPACE and MAXIMAL_CONES sections + stringsParsed = Vector{SubString{String}}[] + for regexp in [r"RAYS\n([\s\S]*?)\nN_RAYS", r"ORTH_LINEALITY_SPACE\n([\s\S]*?)\nF_VECTOR", r"MAXIMAL_CONES\n([\s\S]*)"] + sectionOfInterest = match(regexp, input_string).captures[1] + linesOfInterest = split(sectionOfInterest, "\n") + linesOfInterestFiltered = [split(line, r"\s+#")[1] for line in linesOfInterest if !isempty(line)] + push!(stringsParsed, linesOfInterestFiltered) + end + + # Convert Rays and ORTH_LINEALITY_SPACE to matrices + rayGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[1]]) + linealityGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[2]]) + + # Convert MAXIMAL_CONES to a IncidenceMatrix + maxConeIncidences = IncidenceMatrix([parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]]) + + return polyhedral_fan(maxConeIncidences, rayGenerators, linealityGenerators) + +end + + +function tropical_variety_singular(I::MPolyIdeal,nu::TropicalSemiringMap{QQField,Nothing,typeof(max)}) + R = base_ring(I) + singularConstructRing = "ring r=0,("*join(string.(symbols(R)),",")*"),dp; " + singularConstructIdeal = "ideal I = "*join(string.(gens(I)), ",")*"; " + singularLoadLibrary = "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; }; " + singularComputeTropicalVariety = "fan TropI = tropicalVariety(I); string TropIString = string(TropI);" + singularCommand = singularConstructRing*singularConstructIdeal*singularLoadLibrary*singularComputeTropicalVariety + Singular.call_interpreter(singularCommand) + TropIString = Singular.lookup_library_symbol("Top", "TropIString") + return singular_fan_string_to_oscar_fan(TropIString) +end + + # #======= # tropical variety of an ideal # todo: proper documentation From 450978d2374d6415a0777919c05beec710a25d05 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Sat, 8 Jun 2024 15:02:06 +0100 Subject: [PATCH 02/14] TropicalGeometry: first draft of Singular's tropicalVariety wrapper --- src/TropicalGeometry/variety.jl | 78 ++++++++++++++++++++++++++++----- 1 file changed, 66 insertions(+), 12 deletions(-) diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl index fc468b01779e..05e4a2f498b1 100644 --- a/src/TropicalGeometry/variety.jl +++ b/src/TropicalGeometry/variety.jl @@ -205,7 +205,7 @@ function tropical_variety(I::Union{MPolyIdeal,MPolyRingElem}, nu::Union{Tropical TropV = tropical_variety_linear(P,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) else # P general - error("general tropical varieties currently unsupported") + TropV = tropical_variety_singular(P,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end push!(tropicalVarieties,TropV) end @@ -422,7 +422,11 @@ end # ################################################################################ -function singular_fan_string_to_oscar_fan(input_string::String) + +# this function takes the string of a gfan fan and converts it to an OSCAR polyhedral complex +# if negateFan==true, return the negative of the fan (because Singular is max-convention only, and OSCAR is flexible) +# if dehomogenizeFan==true, the gfan fan is a homogenized polyhedral complex, and we need to dehomogenize it +function singular_fan_string_to_oscar_complex(input_string::String, negateFan::Bool=false, dehomogenizeFan::Bool=false) # Extracting the RAYS, ORTH_LINEALITY_SPACE and MAXIMAL_CONES sections stringsParsed = Vector{SubString{String}}[] @@ -434,29 +438,79 @@ function singular_fan_string_to_oscar_fan(input_string::String) end # Convert Rays and ORTH_LINEALITY_SPACE to matrices + # and negate if necessary rayGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[1]]) linealityGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[2]]) + if negateFan + rayGenerators *= -1 + end + + # Convert MAXIMAL_CONES to a Matrix (for later manipulation) + coneIncidences = Base.stack([parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]],dims=1) + + if dehomogenizeFan + # if the singular fan is a homogenized polyhedral complex, + # identify which rays have a non-zero first entry and actually represent vertices + # then strip the first coordinate of the rays and lineality generators + rayIndices = findall(iszero, rayGenerators[:,1]) + rayGenerators = rayGenerators[:,2:end] + linealityGenerators = linealityGenerators[:,2:end] + # in some cases, the first unit vector can be a lineality generator + # hence we need to filter out potential zero rows from linealityGenerators + linealityGenerators = linealityGenerators[findall(i->!iszero(linealityGenerators[i,:]),1:nrows(linealityGenerators)),:] + return polyhedral_complex(IncidenceMatrix(coneIncidences), rayGenerators, rayIndices, linealityGenerators) + else + # if the singular fan is an honest polyhedral fan + # we need to add the origin as a vertex and add it to all polyhedra + rayIndices = collect(1:nrows(rayGenerators)) + rayGenerators = vcat(rayGenerators, zero_matrix(QQ,1,ncols(rayGenerators))) + coneIncidences = hcat(coneIncidences, ones(Int,nrows(coneIncidences),1)) + return polyhedral_complex(IncidenceMatrix(coneIncidences), rayGenerators, rayIndices, linealityGenerators) + end - # Convert MAXIMAL_CONES to a IncidenceMatrix - maxConeIncidences = IncidenceMatrix([parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]]) +end - return polyhedral_fan(maxConeIncidences, rayGenerators, linealityGenerators) +# this function takes the tropical variety as a polyhedral complex, +# computes the weights of the tropical variety, and constructs the tropical variety +function compute_weights_and_construct_tropical_variety(Sigma::PolyhedralComplex,I::MPolyIdeal,nu::TropicalSemiringMap) + # TODO: compute multiplicities + mults = ones(ZZRingElem,n_maximal_polyhedra(Sigma)) + TropI = TropicalVariety{typeof(convention(nu)),true}(Sigma,mults) + set_attribute!(TropI,:algebraic_ideal,I) + set_attribute!(TropI,:tropical_semiring_map,nu) + return TropI end -function tropical_variety_singular(I::MPolyIdeal,nu::TropicalSemiringMap{QQField,Nothing,typeof(max)}) +# trivial valuation case +function tropical_variety_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) R = base_ring(I) - singularConstructRing = "ring r=0,("*join(string.(symbols(R)),",")*"),dp; " - singularConstructIdeal = "ideal I = "*join(string.(gens(I)), ",")*"; " - singularLoadLibrary = "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; }; " - singularComputeTropicalVariety = "fan TropI = tropicalVariety(I); string TropIString = string(TropI);" - singularCommand = singularConstructRing*singularConstructIdeal*singularLoadLibrary*singularComputeTropicalVariety + singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", + "ideal I = "*join(string.(gens(I)), ",")*";", + "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; };", + "fan TropI = tropicalVariety(I);", + "string TropIString = string(TropI);"]) Singular.call_interpreter(singularCommand) TropIString = Singular.lookup_library_symbol("Top", "TropIString") - return singular_fan_string_to_oscar_fan(TropIString) + Sigma = singular_fan_string_to_oscar_complex(TropIString,convention(nu)==max,false) + TropI = compute_weights_and_construct_tropical_variety(Sigma,I,nu) + return TropI end +function tropical_variety_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,ZZRingElem,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) + R = base_ring(I) + singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", + "ideal I = "*join(string.(gens(I)), ",")*";", + "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; };", + "fan TropI = tropicalVariety(I,number("*string(uniformizer(nu))*"));", + "string TropIString = string(TropI);"]) + Singular.call_interpreter(singularCommand) + TropIString = Singular.lookup_library_symbol("Top", "TropIString") + Sigma = singular_fan_string_to_oscar_complex(TropIString,convention(nu)==max,true) + TropI = compute_weights_and_construct_tropical_variety(Sigma,I,nu) + return TropI +end # #======= # tropical variety of an ideal From 9127307679ed7736c12ce864089b21152935d84e Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Sat, 17 Aug 2024 16:17:50 +0100 Subject: [PATCH 03/14] TropicalGeometry: tropical_variety overhaul --- src/TropicalGeometry/TropicalGeometry.jl | 5 + src/TropicalGeometry/variety.jl | 689 +++--------------- src/TropicalGeometry/variety_affine_linear.jl | 38 + src/TropicalGeometry/variety_binomial.jl | 42 ++ .../variety_equidimensional.jl | 107 +++ src/TropicalGeometry/variety_principal.jl | 20 + .../variety_zerodimensional.jl | 98 +++ 7 files changed, 410 insertions(+), 589 deletions(-) create mode 100644 src/TropicalGeometry/variety_affine_linear.jl create mode 100644 src/TropicalGeometry/variety_binomial.jl create mode 100644 src/TropicalGeometry/variety_equidimensional.jl create mode 100644 src/TropicalGeometry/variety_principal.jl create mode 100644 src/TropicalGeometry/variety_zerodimensional.jl diff --git a/src/TropicalGeometry/TropicalGeometry.jl b/src/TropicalGeometry/TropicalGeometry.jl index bff650ae63fa..b7fd996bdce5 100644 --- a/src/TropicalGeometry/TropicalGeometry.jl +++ b/src/TropicalGeometry/TropicalGeometry.jl @@ -17,5 +17,10 @@ include("hypersurface.jl") include("curve.jl") include("linear_space.jl") include("variety.jl") +include("variety_binomial.jl") +include("variety_principal.jl") +include("variety_affine_linear.jl") +include("variety_zerodimensional.jl") +include("variety_equidimensional.jl") include("intersection.jl") include("groebner_fan.jl") diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl index 05e4a2f498b1..20bded6c430d 100644 --- a/src/TropicalGeometry/variety.jl +++ b/src/TropicalGeometry/variety.jl @@ -134,17 +134,16 @@ end ################################################################################ @doc raw""" - tropical_variety(I::MPolyIdeal, nu::Union{TropicalSemiringMap,Nothing}=nothing; weighted_polyhedral_complex_only::Bool=false, skip_saturation::Bool=false, skip_primary_decomposition::Bool=false) + tropical_variety(I::MPolyIdeal, nu::Union{TropicalSemiringMap,Nothing}=nothing; weighted_polyhedral_complex_only::Bool=false) -Return the tropicalization of `I` with respect to `nu` as a `Vector{TropicalVariety}`. +Return the tropicalization of `I` with respect to `nu`. If `nu==nothing`, will compute with respect to the trivial valuation and min convention. If `weighted_polyhedral_complex_only==true`, will not cache any additional information. -If `skip_saturation==true`, will not saturate `I` with respect to the product of all variables. -If `skip_primary_decomposition==true`, will not decompose `I`. !!! warning - `tropical_variety` is currently under development and only works for ideals that primary decompose into - principal, linear, and binomial ideals. + Assumes that `I` is equi-dimensional. Only special cases supported: + - any valuation: `I` principal, binomial, affine linear + - trivial and p-adic valuation only: `I` general # Examples ```jldoctest @@ -160,616 +159,128 @@ julia> tropical_variety(I) ``` """ -function tropical_variety(I::Union{MPolyIdeal,MPolyRingElem}, nu::Union{TropicalSemiringMap,Nothing}=nothing; weighted_polyhedral_complex_only::Bool=false, skip_saturation::Bool=false, skip_primary_decomposition::Bool=false) - ### - # Step 0.a: convert I to ideal if poly, - # initialize nu as the trivial valuation if not specified by user - ### - if I isa MPolyRingElem - I = ideal(parent(I),[I]) - end +function tropical_variety(I::MPolyIdeal, nu::Union{TropicalSemiringMap,Nothing}=nothing; weighted_polyhedral_complex_only::Bool=false) + # initialize nu as trivial valuation if not given if isnothing(nu) nu = tropical_semiring_map(coefficient_ring(I)) end - ### - # Step 0.b: Saturate `I` and assert that `I` is not the whole ring - ### - if !skip_saturation - R = base_ring(I) - I = saturation(I,ideal([prod(gens(R))])) + # compute a reduced GB to test whether `I` is principal, binomial, or affine linear + GB = groebner_basis(I,complete_reduction=true) + if length(GB)==1 + # I is principal + return tropical_variety_principal(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) + elseif all(g->(length(g)==2), GB) + # I is binomial + return tropical_variety_binomial(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) + elseif all(g->(total_degree(g)<2), GB) + # I affine linear + return tropical_variety_affine_linear(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end - @req !isone(I) "ideal contains a monomial, tropical varieties in OSCAR cannot be empty" - ### - # Step 0.c: Compute primary decomposition and tropical varieties of all primary factors - ### - toTropicalize = [I] - if !skip_primary_decomposition - toTropicalize = [ P for (P,_) in primary_decomposition(I) ] - end - - tropicalVarieties = TropicalVariety[] - for P in toTropicalize - # compute a reduced GB to test whether `P` is principal, binomial, or linear - GB = groebner_basis(P,complete_reduction=true) - - if length(GB)==1 - # P is principal - TropV = tropical_variety_principal(P,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) - elseif max(length.(GB)...)==2 - # P is binomial - TropV = tropical_variety_binomial(P,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) - elseif max(total_degree.(GB)...)==1 - # P linear - TropV = tropical_variety_linear(P,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) - else - # P general - TropV = tropical_variety_singular(P,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) - end - push!(tropicalVarieties,TropV) - end - - return tropicalVarieties + # I general + return tropical_variety_equidimensional(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end -################################################################################ -# -# Principal ideals -# -################################################################################ - -function tropical_variety_principal(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) - ### - # Construct TropicalVariety from TropicalHypersurface - ### - g = first(gens(I)) - TropV = tropical_variety(tropical_hypersurface(g,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) - if !weighted_polyhedral_complex_only - set_attribute!(TropV,:algebraic_ideal,I) - set_attribute!(TropV,:tropical_semiring_map,nu) - end - return TropV +function tropical_variety(f::MPolyRingElem, nu::Union{TropicalSemiringMap,Nothing}=nothing; weighted_polyhedral_complex_only::Bool=false) + return tropical_variety(ideal(parent(f),[f]),nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end - -################################################################################ -# -# Binomial ideals -# -################################################################################ - -function tropical_variety_binomial(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) - ### - # Construct matrix of exponent vector differences - # and vector of coefficient valuation differences - ### - G = gens(I) - A = matrix(ZZ,[first(collect(expv))-last(collect(expv)) for expv in exponents.(G)]) - b = [ QQ(nu(last(collect(coeff))/first(collect(coeff)))) for coeff in coefficients.(G)] - - ### - # Compute tropical variety multiplicity - ### - snfAdiag = elementary_divisors(A) - weight = abs(prod([m for m in snfAdiag if !iszero(m)])) - +function homogenize_pre_tropicalization(I::MPolyIdeal) ### - # Constructing tropical variety set-theoretically + # Compute reduced Groebner basis (usually already cached), and construct homogenization ### - A = QQMatrix(A) - L = transpose(kernel(A, side = :right)) - can_solve, V = can_solve_with_solution(transpose(A),matrix(QQ,[b]),side=:left) - @req can_solve "tropical variety cannot be empty" - SigmaV = polyhedral_complex(IncidenceMatrix([[1]]), V, nothing, L) - - ### - # Assemble tropical variety - ### - TropV = tropical_variety(SigmaV,[weight],convention(nu)) - if !weighted_polyhedral_complex_only - set_attribute!(TropV,:algebraic_ideal,I) - set_attribute!(TropV,:tropical_semiring_map,nu) - end - return TropV -end - - - -################################################################################ -# -# Linear ideals -# -################################################################################ - -function tropical_variety_linear(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) - ### - # Compute reduced Groebner basis (usually already cached), - # and check whether the linear polynomials have a constant term - ### - R = base_ring(I) G = groebner_basis(I,complete_reduction=true) - if min(total_degree.(Iterators.flatten(collect.(terms.(G))))...)==1 - # input homogneeous, construct TropicalVariety via TropicalLinearSpace - TropV = tropical_variety(tropical_linear_space(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) - if !weighted_polyhedral_complex_only - set_attribute!(TropV,:algebraic_ideal,I) - set_attribute!(TropV,:tropical_semiring_map,nu) - end - return TropV - else - # input inhomogeneous, homogenise first - Ih,_,_ = homogenize_pre_tropicalization(I) - TropLh = tropical_linear_space(Ih,nu,weighted_polyhedral_complex_only=true) - Sigma = dehomogenize_post_tropicalization(polyhedral_complex(TropLh)) - - multiplicities = ones(ZZRingElem, n_maximal_polyhedra(Sigma)) - TropV = tropical_variety(Sigma,multiplicities) - if !weighted_polyhedral_complex_only - set_attribute!(TropV,:algebraic_ideal,I) - set_attribute!(TropV,:tropical_semiring_map,nu) + + Kx = base_ring(I) + K = coefficient_ring(Kx) + x = symbols(Kx) + Kxhx,_ = polynomial_ring(K,vcat([:xh],x)) + + Gh = Vector{elem_type(Kx)}(undef,length(G)) + for (i,g) in enumerate(G) + gh = MPolyBuildCtx(Kxhx) + d = total_degree(g) + for (c,alpha) in coefficients_and_exponents(g) + pushfirst!(alpha,d-sum(alpha)) # homogenize exponent vector + push_term!(gh,c,alpha) end - return TropV + Gh[i] = finish(gh) end -end - - -################################################################################ -# -# Zero-dimensional ideals -# (p-adic valuation only) TODO: add trivial valuation -# -################################################################################ - -function tropical_variety_zerodimensional(I::MPolyIdeal,nu::TropicalSemiringMap{QQField,QQFieldElem,<:Union{typeof(min),typeof(max)}}) - - k,(a,_) = number_field(I) - zk = maximal_order(k) - p = uniformizer(nu) - lp = [x[1] for x = prime_decomposition(zk,p)] - ma = representation_matrix(a) - mb = representation_matrix(k(lp[1].gen_two*lp[2].gen_two^2)) - @assert iszero(ma*mb - mb*ma) - Qp = PadicField(p, 10) - TropVDict = simultaneous_diagonalization([map_entries(Qp, ma),map_entries(Qp, mb)]) - - TropVPoints = collect(values(TropVDict)) - TropVPointsUnique = unique(TropVPointsMults) - Sigma = polyhedral_complex(IncidenceMatrix([[i] for i in 1:length(TropVPointsUnique)]), TropVPointsUnique) - TropVMults = [ZZ(length(findall(isequal(p),TropVPoints))) for p in TropVPointsUnique] - TropV = tropical_variety(Sigma,TropVMults) - set_attribute!(TropV,:algebraic_points,collect(keys(TropVDict))) - return TropV + return ideal(Gh) end -### -# Code by Claus: -### -function slope_eigenspace(M::MatElem{T}) where T <: Hecke.NonArchLocalFieldElem - f = charpoly(M) - lf = Hecke.slope_factorization(f) - # @req all(x->x==1, values(lf)) - - se = Dict{typeof(f), typeof(M)}() - k = base_ring(M) - zk = maximal_order(k) +function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) + @req lineality_dim(Sigma)>0 "dehomogenizing polyhedral complex without lineality" - for f = keys(lf) - se[f] = kernel(f(M), side = :right) #hopefully, this is in rref - end - @assert sum(ncols(x) for x = values(se)) == nrows(M) - return se -end + ### + # Construct hyperplane {first coord = 0} + ### + n = ambient_dim(Sigma) + zerothUnitRowVector = zeros(Int,1,n) + zerothUnitRowVector[1,1] = 1 + dehomogenisingHyperplane = polyhedron((zeros(Int,0,n),zeros(Int,0)), (zerothUnitRowVector,[0])) -function _intersect(M::MatElem{T}, N::MatElem{T}) where T <: Hecke.FieldElem - k = base_ring(M) - I = [M N] - PR = maximum(precision, I) - pr = minimum(precision, I) - if pr != PR - for i = eachindex(I) - I[i] = setprecision(I[i], pr) + ### + # Construct matrix and incidence matrix of vertices and rays + ### + incidenceMatrixVertices = Vector{Int}[] + dehomogenizedVertices = Vector{QQFieldElem}[] + incidenceMatrixRays = Vector{Int}[] + dehomogenizedRays = Vector{QQFieldElem}[] + for sigma in maximal_polyhedra(Sigma) + sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) + incidenceVectorVertices = Int[] + V,_ = minimal_faces(sigmaDehomogenized) + for vertex in V + vertex = vertex[2:end] + i = findfirst(isequal(vertex),dehomogenizedVertices) + if i === nothing + push!(dehomogenizedVertices,vertex) + push!(incidenceVectorVertices,length(dehomogenizedVertices)) + else + push!(incidenceVectorVertices,i) + end end - end - - v = kernel(I, side = :right) #precision issues... - l = M*v[1:ncols(M), 1:ncols(v)] - return transpose(rref(transpose(l))[2]) -end - -function valuation_of_roots(f::PolyRingElem{<:Hecke.NonArchLocalFieldElem}) - iszero(f) && error("polynomial must not be zero") - return (valuation(constant_coefficient(f)) - valuation(leading_coefficient(f)))//degree(f) -end - -function simultaneous_diagonalization(v::Vector{<:MatElem{T}}) where T <: Hecke.NonArchLocalFieldElem - - k = base_ring(v[1]) - @assert all(x->base_ring(x) == k, v) - n = nrows(v[1]) - @assert all(x->ncols(x) == nrows(x) == n, v) - - vv = map(slope_eigenspace, v) - - d = Dict(v => [valuation_of_roots(k)] for (k,v) = vv[1]) - @assert sum(ncols(x) for x = keys(d)) == n - for i=2:length(vv) - dd = typeof(d)() - for (mat, pol_vec) = d - for (p, m) = vv[i] - j = _intersect(mat, m) - if ncols(j) > 0 - dd[j] = push!(copy(pol_vec), valuation_of_roots(p)) - end + push!(incidenceMatrixVertices,incidenceVectorVertices) + + incidenceVectorRays = Int[] + R,_ = rays_modulo_lineality(sigmaDehomogenized) + for ray in R + ray = ray[2:end] + i = findfirst(isequal(ray),dehomogenizedRays) + if i === nothing + push!(dehomogenizedRays,ray) + push!(incidenceVectorRays,length(dehomogenizedRays)) + else + push!(incidenceVectorRays,i) end end - d = dd - @assert sum(ncols(x) for x = keys(d)) == n + push!(incidenceMatrixRays,incidenceVectorRays) end - return d -end - - -################################################################################ -# -# positive-dimensional ideals -# (bandaid wrapeer of Singular's tropical variety) -# -################################################################################ - - -# this function takes the string of a gfan fan and converts it to an OSCAR polyhedral complex -# if negateFan==true, return the negative of the fan (because Singular is max-convention only, and OSCAR is flexible) -# if dehomogenizeFan==true, the gfan fan is a homogenized polyhedral complex, and we need to dehomogenize it -function singular_fan_string_to_oscar_complex(input_string::String, negateFan::Bool=false, dehomogenizeFan::Bool=false) - - # Extracting the RAYS, ORTH_LINEALITY_SPACE and MAXIMAL_CONES sections - stringsParsed = Vector{SubString{String}}[] - for regexp in [r"RAYS\n([\s\S]*?)\nN_RAYS", r"ORTH_LINEALITY_SPACE\n([\s\S]*?)\nF_VECTOR", r"MAXIMAL_CONES\n([\s\S]*)"] - sectionOfInterest = match(regexp, input_string).captures[1] - linesOfInterest = split(sectionOfInterest, "\n") - linesOfInterestFiltered = [split(line, r"\s+#")[1] for line in linesOfInterest if !isempty(line)] - push!(stringsParsed, linesOfInterestFiltered) - end - - # Convert Rays and ORTH_LINEALITY_SPACE to matrices - # and negate if necessary - rayGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[1]]) - linealityGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[2]]) - if negateFan - rayGenerators *= -1 - end - - # Convert MAXIMAL_CONES to a Matrix (for later manipulation) - coneIncidences = Base.stack([parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]],dims=1) - - if dehomogenizeFan - # if the singular fan is a homogenized polyhedral complex, - # identify which rays have a non-zero first entry and actually represent vertices - # then strip the first coordinate of the rays and lineality generators - rayIndices = findall(iszero, rayGenerators[:,1]) - rayGenerators = rayGenerators[:,2:end] - linealityGenerators = linealityGenerators[:,2:end] - # in some cases, the first unit vector can be a lineality generator - # hence we need to filter out potential zero rows from linealityGenerators - linealityGenerators = linealityGenerators[findall(i->!iszero(linealityGenerators[i,:]),1:nrows(linealityGenerators)),:] - return polyhedral_complex(IncidenceMatrix(coneIncidences), rayGenerators, rayIndices, linealityGenerators) - else - # if the singular fan is an honest polyhedral fan - # we need to add the origin as a vertex and add it to all polyhedra - rayIndices = collect(1:nrows(rayGenerators)) - rayGenerators = vcat(rayGenerators, zero_matrix(QQ,1,ncols(rayGenerators))) - coneIncidences = hcat(coneIncidences, ones(Int,nrows(coneIncidences),1)) - return polyhedral_complex(IncidenceMatrix(coneIncidences), rayGenerators, rayIndices, linealityGenerators) - end - -end - - -# this function takes the tropical variety as a polyhedral complex, -# computes the weights of the tropical variety, and constructs the tropical variety -function compute_weights_and_construct_tropical_variety(Sigma::PolyhedralComplex,I::MPolyIdeal,nu::TropicalSemiringMap) - # TODO: compute multiplicities - mults = ones(ZZRingElem,n_maximal_polyhedra(Sigma)) - TropI = TropicalVariety{typeof(convention(nu)),true}(Sigma,mults) - set_attribute!(TropI,:algebraic_ideal,I) - set_attribute!(TropI,:tropical_semiring_map,nu) - return TropI -end - - -# trivial valuation case -function tropical_variety_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) - R = base_ring(I) - singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", - "ideal I = "*join(string.(gens(I)), ",")*";", - "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; };", - "fan TropI = tropicalVariety(I);", - "string TropIString = string(TropI);"]) - Singular.call_interpreter(singularCommand) - TropIString = Singular.lookup_library_symbol("Top", "TropIString") - Sigma = singular_fan_string_to_oscar_complex(TropIString,convention(nu)==max,false) - TropI = compute_weights_and_construct_tropical_variety(Sigma,I,nu) - return TropI -end + ### + # Concatenate vertically matrixes of vertices and rays, + # shift incidence matrix of rays and concatenate it horizontally to incicende matrix of vertices, + # dehomogenize generators of lineality space + ### + dehomogenizedVerticesAndRays = matrix(QQ,vcat(dehomogenizedVertices,dehomogenizedRays)) + incidenceMatrixRaysShifted = (x -> x .+length(dehomogenizedVertices)).(incidenceMatrixRays) + incidenceMatrixVerticesAndRays = IncidenceMatrix([vcat(iv,ir) for (iv,ir) in zip(incidenceMatrixVertices,incidenceMatrixRaysShifted)]) -function tropical_variety_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,ZZRingElem,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) - R = base_ring(I) - singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", - "ideal I = "*join(string.(gens(I)), ",")*";", - "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; };", - "fan TropI = tropicalVariety(I,number("*string(uniformizer(nu))*"));", - "string TropIString = string(TropI);"]) - Singular.call_interpreter(singularCommand) - TropIString = Singular.lookup_library_symbol("Top", "TropIString") - Sigma = singular_fan_string_to_oscar_complex(TropIString,convention(nu)==max,true) - TropI = compute_weights_and_construct_tropical_variety(Sigma,I,nu) - return TropI + ### + # Dehomogenize lineality space + ### + sigma = first(maximal_polyhedra(Sigma)) + sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) + dehomogenizedLineality = [linealityVector[2:end] for linealityVector in lineality_space(sigmaDehomogenized)] + + return polyhedral_complex(incidenceMatrixVerticesAndRays, + dehomogenizedVerticesAndRays, + collect(length(dehomogenizedVertices)+1:length(dehomogenizedVertices)+length(dehomogenizedRays)), + dehomogenizedLineality) end - -# #======= -# tropical variety of an ideal -# todo: proper documentation -# Example: -# import Random -# K,s = RationalFunctionField(QQ,"s"); -# Kx,(x1,x2,x3,x4) = polynomial_ring(K,4); -# val = TropicalSemiringMap(K,s); -# I = ideal([x1-s*x2+(s+1)*x3,3*x2-s^2*x3+(s^2+1)*x4]); -# Random.seed!(3847598273423); -# TropI = tropical_variety(I,val) -# =======# -# function tropical_variety(I::MPolyIdeal, val::TropicalSemiringMap, convention::Union{typeof(min),typeof(max)}=min) - -# ### -# # Part 0: Preprocessing -# # Check whether valuation is on the coefficient ring of input polynomials, -# # homogenize input ideal if not homogeneous -# ### -# if coefficient_ring(base_ring(I))!=val.valued_field -# error("input valuation not on coefficient ring of input ideal") -# end -# was_input_homogeneous = true -# for g in groebner_basis(I,complete_reduction=true) # todo: replace GB computation with interreduction -# if !sloppy_is_homogeneous(g) -# was_input_homogeneous = false -# I = homogenize(I) -# break -# end -# end - - - -# ### -# # Part 1: Tropical starting polyhedra (todo: avoid recomputation) -# # - compute and recompute starting points until they lie in the relative interior of maximal cells -# # - initialize working lists -# ### - -# # Note: a working list entry consists of a triple (w,C,G) where -# # * G is a (tropical) Groebner basis -# # * C is the Groebner polyhedron -# # * w is the sum of vertices and rays of C -# # In particular: -# # * w is a weight vector with respect to which G is a Groebner basis, -# # * w is compatible with coordinate permutations if symmetries exist, -# # * instead of comparing C or G it suffices to compare w. -# working_list_todo = [] # list of groebner polyhedra with potentially unknown neighbors -# working_list_done = [] # list of groebner polyhedra with known neighbors -# facet_points_done = [] # list of facet points whose tropical links were computed and traversed - -# compute_starting_points = true -# while compute_starting_points -# print("computing random starting points... ") -# starting_points = tropical_points(I,val) -# println("done") - -# working_list_todo = [] -# working_list_done = [] -# facet_points_done = [] -# compute_starting_points = false -# for starting_point in starting_points -# print("computing groebner_basis for starting point ",starting_point,"... ") -# G = groebner_basis(I,val,starting_point) -# println("done") -# C = groebner_polyhedron(G,val,starting_point) -# w = anchor_point(C) - -# # if C is lower-dimensional, recompute all starting points -# if dim(C)!=dim(I) -# println("starting point on lower-dimensional cell, recomputing...") -# compute_starting_points = true -# break -# end - -# # if (w,C,G) already is in working_list_todo, skip -# i = searchsortedfirst(working_list_todo,(w,C,G),by=x->x[1]) -# if i<=length(working_list_todo) && working_list_todo[i][1]==w -# continue -# end -# # otherwise, add (w,C,G) to todo list -# insert!(working_list_todo, i, (w,C,G)) -# end -# end - - - -# ### -# # Part 2: Tropical traversal -# ### -# while !isempty(working_list_todo) -# print("#working_list_todo: ",length(working_list_todo)," ") -# println("#working_list_done: ",length(working_list_done)) - -# # pick a groebner polyhedron from todo list, add it to the done list, and compute its facet points -# (w,C,G) = popfirst!(working_list_todo) -# i = searchsortedfirst(working_list_done,(w,C,G),by=x->x[1]) -# insert!(working_list_done, i, (w,C,G)) - -# points_to_traverse = facet_points(C) -# for point_to_traverse in points_to_traverse -# # if point was traversed before, skip -# i = searchsortedfirst(facet_points_done,point_to_traverse) -# if i<=length(facet_points_done) && facet_points_done[i]==point_to_traverse -# continue -# end -# # otherwise add point_to_traverse to facet_points_done -# insert!(facet_points_done, i, point_to_traverse) - -# directions_to_traverse = tropical_link(ideal(G),val,point_to_traverse) # todo, this output can be wrong -# for direction_to_traverse in directions_to_traverse -# # compute neighbor -# print("computing groebner_basis for ",point_to_traverse,direction_to_traverse,"... ") -# G_neighbor = groebner_flip(G,val,w,point_to_traverse,direction_to_traverse) -# println("done") -# C_neighbor = groebner_polyhedron(G_neighbor,val,point_to_traverse,perturbation=direction_to_traverse) -# w_neighbor = anchor_point(C_neighbor) - -# # if neighbor is already in done list, skip -# i = searchsortedfirst(working_list_done, -# (w_neighbor,C_neighbor,G_neighbor), -# by=x->x[1]) -# if i<=length(working_list_done) && working_list_done[i][1]==w_neighbor -# continue -# end -# # if neighbor is already in todo list, skip -# i = searchsortedfirst(working_list_todo, -# (w_neighbor,C_neighbor,G_neighbor), -# by=x->x[1]) -# if i<=length(working_list_todo) && working_list_todo[i][1]==w_neighbor -# continue -# end -# # otherwise, add neighbor to todo list -# insert!(working_list_todo, i, (w_neighbor,C_neighbor,G_neighbor)) -# end -# end -# end - - - -# ### -# # Part 3: Postprocessing -# ### -# # 3.0: dehomogenize data if input was homogenized -# if !was_input_homogeneous -# n = length(gens(base_ring(I))) - -# # 3.0.1: dehomogenize Groebner polyhedra -# zeroth_unit_vector_as_row_vector = zeros(Int,1,n) -# zeroth_unit_vector_as_row_vector[1,1] = 1 -# dehomogenising_hyperplane = polyhedron((zeros(Int,0,n),zeros(Int,0)), -# (zeroth_unit_vector_as_row_vector,[1])) -# for wCG in working_list_done -# wCG[2] = intersect(wCG[2],dehomogenising_hyperplane) -# end -# # 3.0.2: check that initial ideals are distinct (todo) -# end - -# # 3.1: construct PolyhedralComplex -# # 3.1.1: construct incidence_matrix, vertices_and_rays, and far_vertices -# incidence_matrix = Vector{Vector{Int}}() -# verts_rays = Vector{Vector{Polymake.Rational}}() -# far_vertices = Vector{Int}() -# for (w,C,G) in working_list_done -# incidence_vector = Vector{Int}() -# for vert in vertices(C) -# i = findfirst(isequal(vert),verts_rays) -# if i === nothing -# # if vert does not occur in verts_rays -# # add it to verts_rays -# push!(verts_rays,vert) -# push!(incidence_vector,length(verts_rays)) -# else -# push!(incidence_vector,i) -# end -# end -# for ray in rays(C) -# i = findfirst(isequal(ray),verts_rays) -# if i === nothing || !(i in far_vertices) -# # if ray does not occur in verts_rays or if it occurs but not as a ray, -# # add it to verts_rays -# push!(verts_rays,ray) -# push!(far_vertices,length(verts_rays)) -# push!(incidence_vector,length(verts_rays)) -# else -# push!(incidence_vector,i) -# end -# end -# push!(incidence_matrix,incidence_vector) - -# end -# verts_rays_matrix = permutedims(reduce(hcat, verts_rays)) # convert Vector{Vector} to Matrix - -# # 3.1.2: construct lineality space -# (w,C,G) = first(working_list_done) -# lineality_space_gens = matrix(QQ,lineality_space(C)) - - - -# # 3.2: Construct lists for weight_vectors, initial_ideals and multiplicities -# weight_vectors = [w for (w,C,G) in working_list_done] -# initial_ideals = [ideal(initial(G,val,w)) for (w,C,G) in working_list_done] -# multiplicities = [multiplicity(inI) for inI in initial_ideals] - - - -# PC = polyhedral_complex(IncidenceMatrix(incidence_matrix), -# verts_rays_matrix, -# far_vertices, -# lineality_space_gens) - -# verts_rays_perm = collect(vertices_and_rays(PC)) -# verts_rays_perm = Vector{Int64}.(verts_rays_perm) -# permutation = [findfirst(isequal(l), verts_rays_perm) for l in verts_rays] -# # inc = [findall(c -> c, incidence_matrix[i, :]) for i in 1:size(incidence_matrix, 1)] -# new_incidence = [permutation[incidence] for incidence in incidence_matrix] -# mults = Dict(new_incidence[i] => multiplicities[i] for i in 1:length(working_list_done)) - -# TropI = TropicalVariety{typeof(max),true}(PC, mults) - - -# set_attribute!(TropI,:weight_vectors,weight_vectors) -# set_attribute!(TropI,:initial_ideals,initial_ideals) - -# return TropI -# end - - - -# #======= -# Example: -# P = cube(4) -# anchor_point(P) -# facet_points(P) -# =======# -# function anchor_point(P::Polyhedron) -# # compute the sum of vertices and rays in homogenized coordinates -# pt = convert(Vector{QQFieldElem},sum([vertices(P)...,rays(P)...])) -# pushfirst!(pt,n_vertices(P)) - -# # project to orthogonal complement of lineality space if necessary -# if lineality_dim(P)>0 -# pt = Polymake.Matrix{Polymake.Rational}(vcat(transpose(pt))) -# Polymake.common.project_to_orthogonal_complement(pt, P.pm_polytope.LINEALITY_SPACE) -# pt = convert(Matrix{QQFieldElem}, pt)[1,:] -# end - -# # rescale until first entry is 1 and remove it -# pt = [pt[i]//pt[1] for i in 2:length(pt)] -# return pt -# end - -# function facet_points(P::Polyhedron) -# points = [] -# for facet in faces(P,dim(P)-1) -# if length(vertices(facet))>0 # skipping facets at infinity -# push!(points,anchor_point(facet)) -# end -# end -# return points -# end diff --git a/src/TropicalGeometry/variety_affine_linear.jl b/src/TropicalGeometry/variety_affine_linear.jl new file mode 100644 index 000000000000..7a84aacc1a5b --- /dev/null +++ b/src/TropicalGeometry/variety_affine_linear.jl @@ -0,0 +1,38 @@ +################################################################################ +# +# Tropicalization of affine linear ideals +# +# WARNING: assumes without test that `gens(I)` has degree one +# +################################################################################ + +function tropical_variety_affine_linear(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) + ### + # Compute reduced Groebner basis (usually already cached), + # and check whether the linear polynomials have a constant term + ### + R = base_ring(I) + G = groebner_basis(I,complete_reduction=true) + if min(total_degree.(Iterators.flatten(collect.(terms.(G))))...)==1 + # input homogneeous, construct TropicalVariety via TropicalLinearSpace + TropV = tropical_variety(tropical_linear_space(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) + if !weighted_polyhedral_complex_only + set_attribute!(TropV,:algebraic_ideal,I) + set_attribute!(TropV,:tropical_semiring_map,nu) + end + return TropV + else + # input inhomogeneous, homogenise first + Ih = homogenize_pre_tropicalization(I) + TropLh = tropical_linear_space(Ih,nu,weighted_polyhedral_complex_only=true) + Sigma = dehomogenize_post_tropicalization(polyhedral_complex(TropLh)) + + multiplicities = ones(ZZRingElem, n_maximal_polyhedra(Sigma)) + TropV = tropical_variety(Sigma,multiplicities) + if !weighted_polyhedral_complex_only + set_attribute!(TropV,:algebraic_ideal,I) + set_attribute!(TropV,:tropical_semiring_map,nu) + end + return TropV + end +end diff --git a/src/TropicalGeometry/variety_binomial.jl b/src/TropicalGeometry/variety_binomial.jl new file mode 100644 index 000000000000..3a014886e284 --- /dev/null +++ b/src/TropicalGeometry/variety_binomial.jl @@ -0,0 +1,42 @@ +################################################################################ +# +# Tropicalization of binomial ideals +# +# WARNING: assumes without test that `gens(I)` is binomial +# +################################################################################ + +function tropical_variety_binomial(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) + ### + # Construct matrix of exponent vector differences + # and vector of coefficient valuation differences + ### + G = gens(I) + A = matrix(ZZ,[first(collect(expv))-last(collect(expv)) for expv in exponents.(G)]) + b = [ QQ(nu(last(collect(coeff))/first(collect(coeff)))) for coeff in coefficients.(G)] + + ### + # Compute tropical variety multiplicity + ### + snfAdiag = elementary_divisors(A) + weight = abs(prod([m for m in snfAdiag if !iszero(m)])) + + ### + # Constructing tropical variety set-theoretically + ### + A = QQMatrix(A) + L = transpose(kernel(A, side = :right)) + can_solve, V = can_solve_with_solution(transpose(A),matrix(QQ,[b]),side=:left) + @req can_solve "tropical variety cannot be empty" + SigmaV = polyhedral_complex(IncidenceMatrix([[1]]), V, nothing, L) + + ### + # Assemble tropical variety + ### + TropV = tropical_variety(SigmaV,[weight],convention(nu)) + if !weighted_polyhedral_complex_only + set_attribute!(TropV,:algebraic_ideal,I) + set_attribute!(TropV,:tropical_semiring_map,nu) + end + return TropV +end diff --git a/src/TropicalGeometry/variety_equidimensional.jl b/src/TropicalGeometry/variety_equidimensional.jl new file mode 100644 index 000000000000..12a1a2842f51 --- /dev/null +++ b/src/TropicalGeometry/variety_equidimensional.jl @@ -0,0 +1,107 @@ +################################################################################ +# +# Tropicalization of equi-dimensional ideals +# +# WARNING: assumes without test that `I` is equi-dimensional +# +################################################################################ + +function tropical_variety_equidimensional(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) + return tropical_variety_equidimensional_singular(I,nu; weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) +end + +# trivial valuation case +function tropical_variety_equidimensional_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) + R = base_ring(I) + singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", + "ideal I = "*join(string.(gens(I)), ",")*";", + "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; };", + "fan TropI = tropicalVariety(I);", + "string TropIString = string(TropI);"]) + Singular.call_interpreter(singularCommand) + TropIString = Singular.lookup_library_symbol("Top", "TropIString") + Sigma = gfan_fan_string_to_oscar_complex(TropIString,convention(nu)==max,false) + TropI = compute_weights_and_construct_tropical_variety(Sigma,I,nu) + if !weighted_polyhedral_complex_only + set_attribute!(TropI,:algebraic_ideal,I) + set_attribute!(TropI,:tropical_semiring_map,nu) + end + return TropI +end + +function tropical_variety_equidimensional_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,ZZRingElem,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) + R = base_ring(I) + singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", + "ideal I = "*join(string.(gens(I)), ",")*";", + "if (!defined(tropicalVariety)) { LIB \"tropical.lib\"; };", + "fan TropI = tropicalVariety(I,number("*string(uniformizer(nu))*"));", + "string TropIString = string(TropI);"]) + Singular.call_interpreter(singularCommand) + TropIString = Singular.lookup_library_symbol("Top", "TropIString") + Sigma = gfan_fan_string_to_oscar_complex(TropIString,convention(nu)==max,true) + TropI = compute_weights_and_construct_tropical_variety(Sigma,I,nu) + if !weighted_polyhedral_complex_only + set_attribute!(TropI,:algebraic_ideal,I) + set_attribute!(TropI,:tropical_semiring_map,nu) + end + return TropI +end + + +# Takes the string of a gfan fan (= Singular fan) and converts it to an OSCAR polyhedral complex +# if negateFan==true, return the negative of the fan (because Singular is max-convention only, and OSCAR is flexible) +# if dehomogenizeFan==true, the gfan fan is a homogenized polyhedral complex, and we need to dehomogenize it +function gfan_fan_string_to_oscar_complex(input_string::String, negateFan::Bool=false, dehomogenizeFan::Bool=false) + + # Extracting the RAYS, LINEALITY_SPACE and MAXIMAL_CONES sections + stringsParsed = Vector{SubString{String}}[] + for regexp in [r"RAYS\n([\s\S]*?)\nN_RAYS", r"LINEALITY_SPACE\n([\s\S]*?)\nORTH_LINEALITY_SPACE", r"MAXIMAL_CONES\n([\s\S]*)"] + sectionOfInterest = match(regexp, input_string).captures[1] + linesOfInterest = split(sectionOfInterest, "\n") + linesOfInterestFiltered = [split(line, r"\s+#")[1] for line in linesOfInterest if !isempty(line)] + push!(stringsParsed, linesOfInterestFiltered) + end + + # Convert Rays and ORTH_LINEALITY_SPACE to matrices + # and negate if necessary + rayGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[1]]) + linealityGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[2]]) + if negateFan + rayGenerators *= -1 + end + + # Convert MAXIMAL_CONES to a Vector{Vector{Int}} + coneIncidences = [parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]] + + if dehomogenizeFan + # if the singular fan is a homogenized polyhedral complex, + # identify which rays have a non-zero first entry and actually represent vertices + # then strip the first coordinate of the rays and lineality generators + rayIndices = findall(iszero, rayGenerators[:,1]) + rayGenerators = rayGenerators[:,2:end] + linealityGenerators = linealityGenerators[:,2:end] + # in some cases, the first unit vector can be a lineality generator + # hence we need to filter out potential zero rows from linealityGenerators + linealityGenerators = linealityGenerators[findall(i->!iszero(linealityGenerators[i,:]),1:nrows(linealityGenerators)),:] + return polyhedral_complex(IncidenceMatrix(coneIncidences), rayGenerators, rayIndices, linealityGenerators) + else + # if the singular fan is an honest polyhedral fan + # we need to add the origin as a vertex and add it to all polyhedra + rayIndices = collect(1:nrows(rayGenerators)) + rayGenerators = vcat(rayGenerators, zero_matrix(QQ,1,ncols(rayGenerators))) + originIndex = nrows(rayGenerators) + coneIncidences = [ vcat(incidence,[originIndex]) for incidence in coneIncidences ] + return polyhedral_complex(IncidenceMatrix(coneIncidences), rayGenerators, rayIndices, linealityGenerators) + end + +end + + +# this function takes the tropical variety as a polyhedral complex, +# computes the weights of the tropical variety, and constructs the tropical variety +function compute_weights_and_construct_tropical_variety(Sigma::PolyhedralComplex,::MPolyIdeal,nu::TropicalSemiringMap) + # TODO: compute multiplicities + mults = ones(ZZRingElem,n_maximal_polyhedra(Sigma)) + TropI = TropicalVariety{typeof(convention(nu)),true}(Sigma,mults) + return TropI +end diff --git a/src/TropicalGeometry/variety_principal.jl b/src/TropicalGeometry/variety_principal.jl new file mode 100644 index 000000000000..46a1a837df71 --- /dev/null +++ b/src/TropicalGeometry/variety_principal.jl @@ -0,0 +1,20 @@ +################################################################################ +# +# Tropicalization of principal ideals +# +# WARNING: assumes without test that `gens(I)` has cardinality one +# +################################################################################ + +function tropical_variety_principal(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) + ### + # Construct TropicalVariety from TropicalHypersurface + ### + g = first(gens(I)) + TropV = tropical_variety(tropical_hypersurface(g,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) + if !weighted_polyhedral_complex_only + set_attribute!(TropV,:algebraic_ideal,I) + set_attribute!(TropV,:tropical_semiring_map,nu) + end + return TropV +end diff --git a/src/TropicalGeometry/variety_zerodimensional.jl b/src/TropicalGeometry/variety_zerodimensional.jl new file mode 100644 index 000000000000..66443912d7fa --- /dev/null +++ b/src/TropicalGeometry/variety_zerodimensional.jl @@ -0,0 +1,98 @@ +################################################################################ +# +# Tropicalization of zero-dimensional ideals +# +# WARNING: assumes without test that `I` is zero-dimensional +# +# WARNING: currently p-adic valuation only +# +################################################################################ + +function tropical_variety_zerodimensional_eigenvalue(I::MPolyIdeal,nu::TropicalSemiringMap{QQField,QQFieldElem,<:Union{typeof(min),typeof(max)}}) + k,(a,_) = number_field(I) + zk = maximal_order(k) + p = uniformizer(nu) + lp = [x[1] for x = prime_decomposition(zk,p)] + ma = representation_matrix(a) + mb = representation_matrix(k(lp[1].gen_two*lp[2].gen_two^2)) + @assert iszero(ma*mb - mb*ma) + Qp = PadicField(p, 10) + TropVDict = simultaneous_diagonalization([map_entries(Qp, ma),map_entries(Qp, mb)]) + + TropVPoints = collect(values(TropVDict)) + TropVPointsUnique = unique(TropVPointsMults) + Sigma = polyhedral_complex(IncidenceMatrix([[i] for i in 1:length(TropVPointsUnique)]), TropVPointsUnique) + TropVMults = [ZZ(length(findall(isequal(p),TropVPoints))) for p in TropVPointsUnique] + TropV = tropical_variety(Sigma,TropVMults) + set_attribute!(TropV,:algebraic_points,collect(keys(TropVDict))) + return TropV +end + + +### +# Code by Claus: +### +function slope_eigenspace(M::MatElem{T}) where T <: Hecke.NonArchLocalFieldElem + f = charpoly(M) + lf = Hecke.slope_factorization(f) + # @req all(x->x==1, values(lf)) + + se = Dict{typeof(f), typeof(M)}() + k = base_ring(M) + zk = maximal_order(k) + + for f = keys(lf) + se[f] = kernel(f(M), side = :right) #hopefully, this is in rref + end + @assert sum(ncols(x) for x = values(se)) == nrows(M) + return se +end + +function _intersect(M::MatElem{T}, N::MatElem{T}) where T <: Hecke.FieldElem + k = base_ring(M) + I = [M N] + PR = maximum(precision, I) + pr = minimum(precision, I) + if pr != PR + for i = eachindex(I) + I[i] = setprecision(I[i], pr) + end + end + + v = kernel(I, side = :right) #precision issues... + l = M*v[1:ncols(M), 1:ncols(v)] + return transpose(rref(transpose(l))[2]) +end + +function valuation_of_roots(f::PolyRingElem{<:Hecke.NonArchLocalFieldElem}) + iszero(f) && error("polynomial must not be zero") + return (valuation(constant_coefficient(f)) - valuation(leading_coefficient(f)))//degree(f) +end + +function simultaneous_diagonalization(v::Vector{<:MatElem{T}}) where T <: Hecke.NonArchLocalFieldElem + + k = base_ring(v[1]) + @assert all(x->base_ring(x) == k, v) + n = nrows(v[1]) + @assert all(x->ncols(x) == nrows(x) == n, v) + + vv = map(slope_eigenspace, v) + + d = Dict(v => [valuation_of_roots(k)] for (k,v) = vv[1]) + @assert sum(ncols(x) for x = keys(d)) == n + for i=2:length(vv) + dd = typeof(d)() + for (mat, pol_vec) = d + for (p, m) = vv[i] + j = _intersect(mat, m) + if ncols(j) > 0 + dd[j] = push!(copy(pol_vec), valuation_of_roots(p)) + end + end + end + d = dd + @assert sum(ncols(x) for x = keys(d)) == n + end + + return d +end From 2a6bcbed9b6029ace86566568fefe35611d86ba6 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Fri, 30 Aug 2024 14:06:41 +0100 Subject: [PATCH 04/14] TropicalGeometry: shorten attributes copying --- src/TropicalGeometry/variety.jl | 16 ++++------------ 1 file changed, 4 insertions(+), 12 deletions(-) diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl index 20bded6c430d..ca7b0a9471bf 100644 --- a/src/TropicalGeometry/variety.jl +++ b/src/TropicalGeometry/variety.jl @@ -81,12 +81,8 @@ end function tropical_variety(TropH::TropicalHypersurface{minOrMax,true}) where {minOrMax<:Union{typeof(min),typeof(max)}} TropV = tropical_variety(polyhedral_complex(TropH),multiplicities(TropH), convention(TropH)) - if has_attribute(TropH,:algebraic_polynomial) - set_attribute!(TropV,:algebraic_ideal,ideal([get_attribute(TropH,:algebraic_polynomial)])) - end - if has_attribute(TropH,:tropical_semiring_map) - set_attribute!(TropV,:tropical_semiring_map,get_attribute(TropH,:tropical_semiring_map)) - end + # copy all attributes + TropV.__attrs = deepcopy(TropH.__attrs) return TropV end @@ -95,12 +91,8 @@ end function tropical_variety(TropL::TropicalLinearSpace{minOrMax,true}) where {minOrMax<:Union{typeof(min),typeof(max)}} TropV = tropical_variety(polyhedral_complex(TropL),multiplicities(TropL), convention(TropL)) - if has_attribute(TropL,:algebraic_ideal) - set_attribute!(TropV,:algebraic_ideal,get_attribute(TropL,:algebraic_ideal)) - end - if has_attribute(TropL,:tropical_semiring_map) - set_attribute!(TropV,:tropical_semiring_map,get_attribute(TropL,:tropical_semiring_map)) - end + # copy all attributes + TropV.__attrs = deepcopy(TropL.__attrs) return TropV end From dbee8e4da00834eeee80b6670037ba2a80823352 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Fri, 30 Aug 2024 14:07:56 +0100 Subject: [PATCH 05/14] TropicalGeometry: tropical_variety overhaul --- src/TropicalGeometry/TropicalGeometry.jl | 2 +- src/TropicalGeometry/variety.jl | 68 ++++++++++++------- src/TropicalGeometry/variety_affine_linear.jl | 34 +++------- ...ty_equidimensional.jl => variety_prime.jl} | 48 ++++++++++--- test/TropicalGeometry/variety.jl | 32 +++++++-- 5 files changed, 121 insertions(+), 63 deletions(-) rename src/TropicalGeometry/{variety_equidimensional.jl => variety_prime.jl} (71%) diff --git a/src/TropicalGeometry/TropicalGeometry.jl b/src/TropicalGeometry/TropicalGeometry.jl index b7fd996bdce5..15b0c05804c3 100644 --- a/src/TropicalGeometry/TropicalGeometry.jl +++ b/src/TropicalGeometry/TropicalGeometry.jl @@ -21,6 +21,6 @@ include("variety_binomial.jl") include("variety_principal.jl") include("variety_affine_linear.jl") include("variety_zerodimensional.jl") -include("variety_equidimensional.jl") +include("variety_prime.jl") include("intersection.jl") include("groebner_fan.jl") diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl index ca7b0a9471bf..f4f78b0837bc 100644 --- a/src/TropicalGeometry/variety.jl +++ b/src/TropicalGeometry/variety.jl @@ -24,10 +24,10 @@ end # ################################################################################ -function Base.show(io::IO, tv::TropicalVariety{typeof(min), true}) +function Base.show(io::IO, ::TropicalVariety{typeof(min), true}) print(io, "Min tropical variety") end -function Base.show(io::IO, tv::TropicalVariety{typeof(max), true}) +function Base.show(io::IO, ::TropicalVariety{typeof(max), true}) print(io, "Max tropical variety") end @@ -40,7 +40,7 @@ end ################################################################################ @doc raw""" - tropical_variety(Sigma::PolyhedralComplex, mult, minOrMax::Union{typeof(min),typeof(max)}=min) + tropical_variety(Sigma::PolyhedralComplex, mult::Vector{ZZRingElem}, minOrMax::Union{typeof(min),typeof(max)}=min) Return the `TropicalVariety` whose polyhedral complex is `Sigma` with multiplicities `mult` and convention `minOrMax`. Here, `mult` is optional can be specified as a `Vector{ZZRingElem}` which represents a list of multiplicities on the maximal polyhedra in the order of `maximal_polyhedra(Sigma)`. If `mult` is unspecified, then all multiplicities are set to one. @@ -133,21 +133,37 @@ If `nu==nothing`, will compute with respect to the trivial valuation and min con If `weighted_polyhedral_complex_only==true`, will not cache any additional information. !!! warning - Assumes that `I` is equi-dimensional. Only special cases supported: - - any valuation: `I` principal, binomial, affine linear - - trivial and p-adic valuation only: `I` general + Experimental feature, only special cases supported: + - any coefficient field and any valuation: `I` principal, binomial, or affine linear + - QQ and trivial / p-adic valuation only: `I` prime # Examples ```jldoctest -julia> R,(x,y) = QQ["x","y"]; +julia> R,(x,y,z) = QQ["x","y","z"]; -julia> I = ideal([(x^2+y)*(x+y^2)*(x+y)]); +julia> nu_2 = tropical_semiring_map(QQ,2) +Map into Min tropical semiring encoding the 2-adic valuation on Rational field -julia> tropical_variety(I) -3-element Vector{TropicalVariety}: - Min tropical variety - Min tropical variety - Min tropical variety +julia> f1 = 8*x^2 + x*y + x*z + x + 8*y^2 + y*z + y + 8*z^2 + z + 8; + +julia> f2 = x + 1; + +julia> I = ideal([f1,f2]); + +julia> TropI_0 = tropical_variety(I) +Min tropical variety + +julia> vertices(TropI_0) +1-element SubObjectIterator{PointVector{QQFieldElem}}: + [0, 0, 0] + +julia> TropI_2 = tropical_variety(I,nu_2) +Min tropical variety + +julia> vertices(TropI_2) +2-element SubObjectIterator{PointVector{QQFieldElem}}: + [0, -3, 3] + [0, 3, -3] ``` """ @@ -171,7 +187,7 @@ function tropical_variety(I::MPolyIdeal, nu::Union{TropicalSemiringMap,Nothing}= end # I general - return tropical_variety_equidimensional(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) + return tropical_variety_prime(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end @@ -206,13 +222,14 @@ function homogenize_pre_tropicalization(I::MPolyIdeal) end -function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) - @req lineality_dim(Sigma)>0 "dehomogenizing polyhedral complex without lineality" +function dehomogenize_post_tropicalization(TropV::TropicalVarietySupertype) + + @req lineality_dim(TropV)>0 "dehomogenizing polyhedral complex without lineality" ### # Construct hyperplane {first coord = 0} ### - n = ambient_dim(Sigma) + n = ambient_dim(TropV) zerothUnitRowVector = zeros(Int,1,n) zerothUnitRowVector[1,1] = 1 dehomogenisingHyperplane = polyhedron((zeros(Int,0,n),zeros(Int,0)), (zerothUnitRowVector,[0])) @@ -224,7 +241,7 @@ function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) dehomogenizedVertices = Vector{QQFieldElem}[] incidenceMatrixRays = Vector{Int}[] dehomogenizedRays = Vector{QQFieldElem}[] - for sigma in maximal_polyhedra(Sigma) + for sigma in maximal_polyhedra(TropV) sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) incidenceVectorVertices = Int[] V,_ = minimal_faces(sigmaDehomogenized) @@ -267,12 +284,17 @@ function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) ### # Dehomogenize lineality space ### - sigma = first(maximal_polyhedra(Sigma)) + sigma = first(maximal_polyhedra(TropV)) sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) dehomogenizedLineality = [linealityVector[2:end] for linealityVector in lineality_space(sigmaDehomogenized)] - return polyhedral_complex(incidenceMatrixVerticesAndRays, - dehomogenizedVerticesAndRays, - collect(length(dehomogenizedVertices)+1:length(dehomogenizedVertices)+length(dehomogenizedRays)), - dehomogenizedLineality) + SigmaDehom = polyhedral_complex(incidenceMatrixVerticesAndRays, + dehomogenizedVerticesAndRays, + collect(length(dehomogenizedVertices)+1:length(dehomogenizedVertices)+length(dehomogenizedRays)), + dehomogenizedLineality) + + TropVDehom = tropical_variety(SigmaDehom,multiplicities(TropV),convention(TropV)) + # TropVDehom.__attrs = deepcopy(TropV.__attrs) # TODO: not working, how to copy all attributes? + return TropVDehom + end diff --git a/src/TropicalGeometry/variety_affine_linear.jl b/src/TropicalGeometry/variety_affine_linear.jl index 7a84aacc1a5b..f8b6516cd842 100644 --- a/src/TropicalGeometry/variety_affine_linear.jl +++ b/src/TropicalGeometry/variety_affine_linear.jl @@ -7,32 +7,16 @@ ################################################################################ function tropical_variety_affine_linear(I::MPolyIdeal,nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) - ### - # Compute reduced Groebner basis (usually already cached), - # and check whether the linear polynomials have a constant term - ### - R = base_ring(I) + # compute a reduced GB to check whether I is homogeneous G = groebner_basis(I,complete_reduction=true) - if min(total_degree.(Iterators.flatten(collect.(terms.(G))))...)==1 - # input homogneeous, construct TropicalVariety via TropicalLinearSpace - TropV = tropical_variety(tropical_linear_space(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) - if !weighted_polyhedral_complex_only - set_attribute!(TropV,:algebraic_ideal,I) - set_attribute!(TropV,:tropical_semiring_map,nu) - end - return TropV - else - # input inhomogeneous, homogenise first - Ih = homogenize_pre_tropicalization(I) - TropLh = tropical_linear_space(Ih,nu,weighted_polyhedral_complex_only=true) - Sigma = dehomogenize_post_tropicalization(polyhedral_complex(TropLh)) - multiplicities = ones(ZZRingElem, n_maximal_polyhedra(Sigma)) - TropV = tropical_variety(Sigma,multiplicities) - if !weighted_polyhedral_complex_only - set_attribute!(TropV,:algebraic_ideal,I) - set_attribute!(TropV,:tropical_semiring_map,nu) - end - return TropV + if all(Oscar._is_homogeneous.(G)) + # input homogneeous, construct TropicalVariety via TropicalLinearSpace + return tropical_variety(tropical_linear_space(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) end + + # input inhomogeneous, homogenise first + Ih = homogenize_pre_tropicalization(I) + TropLh = tropical_linear_space(Ih,nu,weighted_polyhedral_complex_only=true) + return dehomogenize_post_tropicalization(TropLh) end diff --git a/src/TropicalGeometry/variety_equidimensional.jl b/src/TropicalGeometry/variety_prime.jl similarity index 71% rename from src/TropicalGeometry/variety_equidimensional.jl rename to src/TropicalGeometry/variety_prime.jl index 12a1a2842f51..541501e0672d 100644 --- a/src/TropicalGeometry/variety_equidimensional.jl +++ b/src/TropicalGeometry/variety_prime.jl @@ -1,17 +1,26 @@ ################################################################################ # -# Tropicalization of equi-dimensional ideals +# Tropicalization of prime ideals # -# WARNING: assumes without test that `I` is equi-dimensional +# WARNING: assumes without test that `I` is prime # ################################################################################ -function tropical_variety_equidimensional(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) - return tropical_variety_equidimensional_singular(I,nu; weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) +function tropical_variety_prime(I::MPolyIdeal, nu::TropicalSemiringMap; weighted_polyhedral_complex_only::Bool=false) + # compute a reduced GB to check whether I is homogeneous + G = groebner_basis(I,complete_reduction=true) + + if all(Oscar._is_homogeneous.(G)) + return tropical_variety_prime_singular(I,nu; weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) + end + + Ih = homogenize_pre_tropicalization(I) + TropIh = tropical_variety_prime_singular(Ih,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) + return dehomogenize_post_tropicalization(TropIh) end -# trivial valuation case -function tropical_variety_equidimensional_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) +# trivial valuation +function tropical_variety_prime_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,Nothing,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) R = base_ring(I) singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", "ideal I = "*join(string.(gens(I)), ",")*";", @@ -29,7 +38,8 @@ function tropical_variety_equidimensional_singular(I::MPolyIdeal, nu::TropicalSe return TropI end -function tropical_variety_equidimensional_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,ZZRingElem,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) +# p-adic valuation +function tropical_variety_prime_singular(I::MPolyIdeal, nu::TropicalSemiringMap{QQField,ZZRingElem,<:Union{typeof(min),typeof(max)}}; weighted_polyhedral_complex_only::Bool=false) R = base_ring(I) singularCommand = join(["ring r=0,("*join(string.(symbols(R)),",")*"),dp;", "ideal I = "*join(string.(gens(I)), ",")*";", @@ -53,6 +63,9 @@ end # if dehomogenizeFan==true, the gfan fan is a homogenized polyhedral complex, and we need to dehomogenize it function gfan_fan_string_to_oscar_complex(input_string::String, negateFan::Bool=false, dehomogenizeFan::Bool=false) + # Extracting AMBIENT_DIM + ambientDim = parse(Int, match(r"AMBIENT_DIM\n(\d+)", input_string).captures[1]) + # Extracting the RAYS, LINEALITY_SPACE and MAXIMAL_CONES sections stringsParsed = Vector{SubString{String}}[] for regexp in [r"RAYS\n([\s\S]*?)\nN_RAYS", r"LINEALITY_SPACE\n([\s\S]*?)\nORTH_LINEALITY_SPACE", r"MAXIMAL_CONES\n([\s\S]*)"] @@ -64,14 +77,29 @@ function gfan_fan_string_to_oscar_complex(input_string::String, negateFan::Bool= # Convert Rays and ORTH_LINEALITY_SPACE to matrices # and negate if necessary - rayGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[1]]) - linealityGenerators = matrix(QQ,[parse.(Int, split(line)) for line in stringsParsed[2]]) + rayGenerators = [parse.(Int, split(line)) for line in stringsParsed[1]] + if isempty(rayGenerators) + rayGenerators = zero_matrix(QQ,0,ambientDim) + else + rayGenerators = matrix(QQ,rayGenerators) + end if negateFan rayGenerators *= -1 end + linealityGenerators = [parse.(Int, split(line)) for line in stringsParsed[2]] + if isempty(linealityGenerators) + linealityGenerators = zero_matrix(QQ,0,ambientDim) + else + linealityGenerators = matrix(QQ,linealityGenerators) + end + # Convert MAXIMAL_CONES to a Vector{Vector{Int}} - coneIncidences = [parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]] + if length(stringsParsed[3])==1 && first(stringsParsed[3])=="{}" + coneIncidences = [Int[]] + else + coneIncidences = [parse.(Int, split(replace(line, r"[{}]" => ""), r"\s+")) .+ 1 for line in stringsParsed[3]] + end if dehomogenizeFan # if the singular fan is a homogenized polyhedral complex, diff --git a/test/TropicalGeometry/variety.jl b/test/TropicalGeometry/variety.jl index a7d3a5776d4a..bfe139340aff 100644 --- a/test/TropicalGeometry/variety.jl +++ b/test/TropicalGeometry/variety.jl @@ -19,11 +19,11 @@ R,(x,y,z,w) = QQ["x","y","z","w"] I = ideal(R,[x+2*y,z+4*w]) nu = tropical_semiring_map(QQ,2) - TropV = first(tropical_variety(I,nu)) + TropV = tropical_variety(I,nu) TropL = tropical_linear_space(I,nu) @test issetequal(maximal_polyhedra(TropV),maximal_polyhedra(TropL)) nu = tropical_semiring_map(QQ,2,max) - TropV = first(tropical_variety(I,nu)) + TropV = tropical_variety(I,nu) TropL = tropical_linear_space(I,nu) @test issetequal(maximal_polyhedra(TropV),maximal_polyhedra(TropL)) end @@ -33,12 +33,36 @@ f = x*y*z+2 nu = tropical_semiring_map(QQ,2) TropH = tropical_hypersurface(f,nu) - TropV = first(tropical_variety(ideal(R,f),nu)) + TropV = tropical_variety(ideal(R,f),nu) @test issetequal(maximal_polyhedra(TropH),maximal_polyhedra(TropV)) nu = tropical_semiring_map(QQ,2,max) - TropV = first(tropical_variety(ideal(R,f),nu)) + TropV = tropical_variety(ideal(R,f),nu) TropH = tropical_hypersurface(f,nu) @test issetequal(maximal_polyhedra(TropH),maximal_polyhedra(TropV)) end + # running tropical_variety and all its subroutines + @testset "testing tropical_variety" begin + + # principal ideals + R,(x,y,z) = QQ["x","y","z"] + f = x^2+y^2+z^2+1 + TropV = tropical_variety(ideal(R,f)) + @test f_vector(TropV) == [1,4,6] + + # binomial ideals + f = x^2+1 + g = y^2+1 + TropV = tropical_variety(ideal(R,[f,g])) + @test f_vector(TropV) == [0,1] + + # affine linear ideals + f = x+z+1 + g = y+z+1 + TropV = tropical_variety(ideal(R,[f,g])) + @test f_vector(TropV) == [1,3] + + # general ideals, see doctests + + end end From 8cc7c6827bb729304ebf06ab11c400af6f11d563 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Tue, 10 Sep 2024 20:29:44 +0100 Subject: [PATCH 06/14] Update src/TropicalGeometry/variety_binomial.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_binomial.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/TropicalGeometry/variety_binomial.jl b/src/TropicalGeometry/variety_binomial.jl index 3a014886e284..70e029c891ea 100644 --- a/src/TropicalGeometry/variety_binomial.jl +++ b/src/TropicalGeometry/variety_binomial.jl @@ -12,8 +12,8 @@ function tropical_variety_binomial(I::MPolyIdeal,nu::TropicalSemiringMap; weight # and vector of coefficient valuation differences ### G = gens(I) - A = matrix(ZZ,[first(collect(expv))-last(collect(expv)) for expv in exponents.(G)]) - b = [ QQ(nu(last(collect(coeff))/first(collect(coeff)))) for coeff in coefficients.(G)] + A = matrix(ZZ,[first(expv)-last(expv) for expv in exponents.(G)]) + b = [ QQ(nu(last(coeff)/first(coeff))) for coeff in coefficients.(G)] ### # Compute tropical variety multiplicity From cdb1fa2f33a58579e40e3098030c70f5993482aa Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Tue, 10 Sep 2024 20:30:07 +0100 Subject: [PATCH 07/14] Update src/TropicalGeometry/variety_binomial.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_binomial.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/TropicalGeometry/variety_binomial.jl b/src/TropicalGeometry/variety_binomial.jl index 70e029c891ea..f3cbca373279 100644 --- a/src/TropicalGeometry/variety_binomial.jl +++ b/src/TropicalGeometry/variety_binomial.jl @@ -24,7 +24,7 @@ function tropical_variety_binomial(I::MPolyIdeal,nu::TropicalSemiringMap; weight ### # Constructing tropical variety set-theoretically ### - A = QQMatrix(A) + A = matrix(QQ, A) L = transpose(kernel(A, side = :right)) can_solve, V = can_solve_with_solution(transpose(A),matrix(QQ,[b]),side=:left) @req can_solve "tropical variety cannot be empty" From 323f905af6c73de7ccdc0cdc562cd0aaf3da87e1 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Tue, 10 Sep 2024 20:30:33 +0100 Subject: [PATCH 08/14] Update src/TropicalGeometry/variety_binomial.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_binomial.jl | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/src/TropicalGeometry/variety_binomial.jl b/src/TropicalGeometry/variety_binomial.jl index f3cbca373279..89ef25b4e35a 100644 --- a/src/TropicalGeometry/variety_binomial.jl +++ b/src/TropicalGeometry/variety_binomial.jl @@ -25,8 +25,7 @@ function tropical_variety_binomial(I::MPolyIdeal,nu::TropicalSemiringMap; weight # Constructing tropical variety set-theoretically ### A = matrix(QQ, A) - L = transpose(kernel(A, side = :right)) - can_solve, V = can_solve_with_solution(transpose(A),matrix(QQ,[b]),side=:left) + can_solve, V, L = can_solve_with_solution_and_kernel(transpose(A), matrix(QQ,[b]); side=:left) @req can_solve "tropical variety cannot be empty" SigmaV = polyhedral_complex(IncidenceMatrix([[1]]), V, nothing, L) From 218f5413f87d3a869a7c593c09ba676188efb242 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Sun, 15 Sep 2024 15:24:52 +0100 Subject: [PATCH 09/14] Update src/TropicalGeometry/variety_prime.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_prime.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/TropicalGeometry/variety_prime.jl b/src/TropicalGeometry/variety_prime.jl index 541501e0672d..5670492db3d5 100644 --- a/src/TropicalGeometry/variety_prime.jl +++ b/src/TropicalGeometry/variety_prime.jl @@ -10,7 +10,7 @@ function tropical_variety_prime(I::MPolyIdeal, nu::TropicalSemiringMap; weighted # compute a reduced GB to check whether I is homogeneous G = groebner_basis(I,complete_reduction=true) - if all(Oscar._is_homogeneous.(G)) + if all(Oscar._is_homogeneous, G) return tropical_variety_prime_singular(I,nu; weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end From 46b4e04d0eaea7e43621decd7b68f2885564b606 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Mon, 16 Sep 2024 09:54:33 +0100 Subject: [PATCH 10/14] Update src/TropicalGeometry/variety_zerodimensional.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_zerodimensional.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/TropicalGeometry/variety_zerodimensional.jl b/src/TropicalGeometry/variety_zerodimensional.jl index 66443912d7fa..f1c0ae3a5d75 100644 --- a/src/TropicalGeometry/variety_zerodimensional.jl +++ b/src/TropicalGeometry/variety_zerodimensional.jl @@ -16,7 +16,7 @@ function tropical_variety_zerodimensional_eigenvalue(I::MPolyIdeal,nu::TropicalS ma = representation_matrix(a) mb = representation_matrix(k(lp[1].gen_two*lp[2].gen_two^2)) @assert iszero(ma*mb - mb*ma) - Qp = PadicField(p, 10) + Qp = padic_field(p; precision=10) TropVDict = simultaneous_diagonalization([map_entries(Qp, ma),map_entries(Qp, mb)]) TropVPoints = collect(values(TropVDict)) From 4ff4584041950b9e25cc42b747836cb5de970bde Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Mon, 16 Sep 2024 09:55:11 +0100 Subject: [PATCH 11/14] Update src/TropicalGeometry/variety_affine_linear.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_affine_linear.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/TropicalGeometry/variety_affine_linear.jl b/src/TropicalGeometry/variety_affine_linear.jl index f8b6516cd842..dc78caafd5f9 100644 --- a/src/TropicalGeometry/variety_affine_linear.jl +++ b/src/TropicalGeometry/variety_affine_linear.jl @@ -10,7 +10,7 @@ function tropical_variety_affine_linear(I::MPolyIdeal,nu::TropicalSemiringMap; w # compute a reduced GB to check whether I is homogeneous G = groebner_basis(I,complete_reduction=true) - if all(Oscar._is_homogeneous.(G)) + if all(Oscar._is_homogeneous,G) # input homogneeous, construct TropicalVariety via TropicalLinearSpace return tropical_variety(tropical_linear_space(I,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) end From 05818f967f264af9f62900ff5309f9136a185306 Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Mon, 16 Sep 2024 09:55:23 +0100 Subject: [PATCH 12/14] Update src/TropicalGeometry/variety_principal.jl Co-authored-by: Max Horn --- src/TropicalGeometry/variety_principal.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/TropicalGeometry/variety_principal.jl b/src/TropicalGeometry/variety_principal.jl index 46a1a837df71..451aec10ceb6 100644 --- a/src/TropicalGeometry/variety_principal.jl +++ b/src/TropicalGeometry/variety_principal.jl @@ -10,7 +10,7 @@ function tropical_variety_principal(I::MPolyIdeal,nu::TropicalSemiringMap; weigh ### # Construct TropicalVariety from TropicalHypersurface ### - g = first(gens(I)) + g = gen(I,1) TropV = tropical_variety(tropical_hypersurface(g,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only)) if !weighted_polyhedral_complex_only set_attribute!(TropV,:algebraic_ideal,I) From 44cd831861f8e600b106a3091b9c6c57ce7aacbb Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Wed, 18 Sep 2024 07:16:04 +0100 Subject: [PATCH 13/14] Udate src/TropicalGeometry/variety_binomial.jl --- src/TropicalGeometry/variety_binomial.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/TropicalGeometry/variety_binomial.jl b/src/TropicalGeometry/variety_binomial.jl index 89ef25b4e35a..81b2dc997ec7 100644 --- a/src/TropicalGeometry/variety_binomial.jl +++ b/src/TropicalGeometry/variety_binomial.jl @@ -12,8 +12,8 @@ function tropical_variety_binomial(I::MPolyIdeal,nu::TropicalSemiringMap; weight # and vector of coefficient valuation differences ### G = gens(I) - A = matrix(ZZ,[first(expv)-last(expv) for expv in exponents.(G)]) - b = [ QQ(nu(last(coeff)/first(coeff))) for coeff in coefficients.(G)] + A = matrix(ZZ,[first(expv)-last(expv) for expv in collect.(exponents.(G))]) + b = [ QQ(nu(last(coeff)/first(coeff))) for coeff in collect.(coefficients.(G))] ### # Compute tropical variety multiplicity From b582d0bdf111bdf71cb5ff9aafe17336b9a299ab Mon Sep 17 00:00:00 2001 From: Yue Ren Date: Wed, 18 Sep 2024 07:16:30 +0100 Subject: [PATCH 14/14] TropicalGeometry: resolved merge conflict --- src/TropicalGeometry/TropicalGeometry.jl | 2 +- src/TropicalGeometry/homogenization.jl | 24 ++-- src/TropicalGeometry/variety.jl | 104 ------------------ src/TropicalGeometry/variety_affine_linear.jl | 2 +- src/TropicalGeometry/variety_prime.jl | 2 +- 5 files changed, 18 insertions(+), 116 deletions(-) diff --git a/src/TropicalGeometry/TropicalGeometry.jl b/src/TropicalGeometry/TropicalGeometry.jl index 15b0c05804c3..c9cdc655bbea 100644 --- a/src/TropicalGeometry/TropicalGeometry.jl +++ b/src/TropicalGeometry/TropicalGeometry.jl @@ -7,12 +7,12 @@ include("semiring.jl") include("semiring_map.jl") include("matrix.jl") include("poly.jl") +include("variety_supertype.jl") include("homogenization.jl") include("groebner_basis.jl") include("initial.jl") # include("groebner_polyhedron.jl") # include("points.jl") -include("variety_supertype.jl") include("hypersurface.jl") include("curve.jl") include("linear_space.jl") diff --git a/src/TropicalGeometry/homogenization.jl b/src/TropicalGeometry/homogenization.jl index d58d78093319..a99d56f21f54 100644 --- a/src/TropicalGeometry/homogenization.jl +++ b/src/TropicalGeometry/homogenization.jl @@ -34,13 +34,14 @@ function homogenize_pre_tropicalization(I::MPolyIdeal) end -function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) - @req lineality_dim(Sigma)>0 "dehomogenizing polyhedral complex without lineality" +function dehomogenize_post_tropicalization(TropV::TropicalVarietySupertype) + + @req lineality_dim(TropV)>0 "dehomogenizing polyhedral complex without lineality" ### # Construct hyperplane {first coord = 0} ### - n = ambient_dim(Sigma) + n = ambient_dim(TropV) zerothUnitRowVector = zeros(Int,1,n) zerothUnitRowVector[1,1] = 1 dehomogenisingHyperplane = polyhedron((zeros(Int,0,n),zeros(Int,0)), (zerothUnitRowVector,[0])) @@ -52,7 +53,7 @@ function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) dehomogenizedVertices = Vector{QQFieldElem}[] incidenceMatrixRays = Vector{Int}[] dehomogenizedRays = Vector{QQFieldElem}[] - for sigma in maximal_polyhedra(Sigma) + for sigma in maximal_polyhedra(TropV) sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) incidenceVectorVertices = Int[] V,_ = minimal_faces(sigmaDehomogenized) @@ -95,12 +96,17 @@ function dehomogenize_post_tropicalization(Sigma::PolyhedralComplex) ### # Dehomogenize lineality space ### - sigma = first(maximal_polyhedra(Sigma)) + sigma = first(maximal_polyhedra(TropV)) sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) dehomogenizedLineality = [linealityVector[2:end] for linealityVector in lineality_space(sigmaDehomogenized)] - return polyhedral_complex(incidenceMatrixVerticesAndRays, - dehomogenizedVerticesAndRays, - collect(length(dehomogenizedVertices)+1:length(dehomogenizedVertices)+length(dehomogenizedRays)), - dehomogenizedLineality) + SigmaDehom = polyhedral_complex(incidenceMatrixVerticesAndRays, + dehomogenizedVerticesAndRays, + collect(length(dehomogenizedVertices)+1:length(dehomogenizedVertices)+length(dehomogenizedRays)), + dehomogenizedLineality) + + TropVDehom = tropical_variety(SigmaDehom,multiplicities(TropV),convention(TropV)) + # TropVDehom.__attrs = deepcopy(TropV.__attrs) # TODO: not working, how to copy all attributes? + return TropVDehom + end diff --git a/src/TropicalGeometry/variety.jl b/src/TropicalGeometry/variety.jl index f4f78b0837bc..f4e209134edd 100644 --- a/src/TropicalGeometry/variety.jl +++ b/src/TropicalGeometry/variety.jl @@ -194,107 +194,3 @@ end function tropical_variety(f::MPolyRingElem, nu::Union{TropicalSemiringMap,Nothing}=nothing; weighted_polyhedral_complex_only::Bool=false) return tropical_variety(ideal(parent(f),[f]),nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end - - -function homogenize_pre_tropicalization(I::MPolyIdeal) - ### - # Compute reduced Groebner basis (usually already cached), and construct homogenization - ### - G = groebner_basis(I,complete_reduction=true) - - Kx = base_ring(I) - K = coefficient_ring(Kx) - x = symbols(Kx) - Kxhx,_ = polynomial_ring(K,vcat([:xh],x)) - - Gh = Vector{elem_type(Kx)}(undef,length(G)) - for (i,g) in enumerate(G) - gh = MPolyBuildCtx(Kxhx) - d = total_degree(g) - for (c,alpha) in coefficients_and_exponents(g) - pushfirst!(alpha,d-sum(alpha)) # homogenize exponent vector - push_term!(gh,c,alpha) - end - Gh[i] = finish(gh) - end - - return ideal(Gh) -end - - -function dehomogenize_post_tropicalization(TropV::TropicalVarietySupertype) - - @req lineality_dim(TropV)>0 "dehomogenizing polyhedral complex without lineality" - - ### - # Construct hyperplane {first coord = 0} - ### - n = ambient_dim(TropV) - zerothUnitRowVector = zeros(Int,1,n) - zerothUnitRowVector[1,1] = 1 - dehomogenisingHyperplane = polyhedron((zeros(Int,0,n),zeros(Int,0)), (zerothUnitRowVector,[0])) - - ### - # Construct matrix and incidence matrix of vertices and rays - ### - incidenceMatrixVertices = Vector{Int}[] - dehomogenizedVertices = Vector{QQFieldElem}[] - incidenceMatrixRays = Vector{Int}[] - dehomogenizedRays = Vector{QQFieldElem}[] - for sigma in maximal_polyhedra(TropV) - sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) - incidenceVectorVertices = Int[] - V,_ = minimal_faces(sigmaDehomogenized) - for vertex in V - vertex = vertex[2:end] - i = findfirst(isequal(vertex),dehomogenizedVertices) - if i === nothing - push!(dehomogenizedVertices,vertex) - push!(incidenceVectorVertices,length(dehomogenizedVertices)) - else - push!(incidenceVectorVertices,i) - end - end - push!(incidenceMatrixVertices,incidenceVectorVertices) - - incidenceVectorRays = Int[] - R,_ = rays_modulo_lineality(sigmaDehomogenized) - for ray in R - ray = ray[2:end] - i = findfirst(isequal(ray),dehomogenizedRays) - if i === nothing - push!(dehomogenizedRays,ray) - push!(incidenceVectorRays,length(dehomogenizedRays)) - else - push!(incidenceVectorRays,i) - end - end - push!(incidenceMatrixRays,incidenceVectorRays) - end - - ### - # Concatenate vertically matrixes of vertices and rays, - # shift incidence matrix of rays and concatenate it horizontally to incicende matrix of vertices, - # dehomogenize generators of lineality space - ### - dehomogenizedVerticesAndRays = matrix(QQ,vcat(dehomogenizedVertices,dehomogenizedRays)) - incidenceMatrixRaysShifted = (x -> x .+length(dehomogenizedVertices)).(incidenceMatrixRays) - incidenceMatrixVerticesAndRays = IncidenceMatrix([vcat(iv,ir) for (iv,ir) in zip(incidenceMatrixVertices,incidenceMatrixRaysShifted)]) - - ### - # Dehomogenize lineality space - ### - sigma = first(maximal_polyhedra(TropV)) - sigmaDehomogenized = intersect(sigma,dehomogenisingHyperplane) - dehomogenizedLineality = [linealityVector[2:end] for linealityVector in lineality_space(sigmaDehomogenized)] - - SigmaDehom = polyhedral_complex(incidenceMatrixVerticesAndRays, - dehomogenizedVerticesAndRays, - collect(length(dehomogenizedVertices)+1:length(dehomogenizedVertices)+length(dehomogenizedRays)), - dehomogenizedLineality) - - TropVDehom = tropical_variety(SigmaDehom,multiplicities(TropV),convention(TropV)) - # TropVDehom.__attrs = deepcopy(TropV.__attrs) # TODO: not working, how to copy all attributes? - return TropVDehom - -end diff --git a/src/TropicalGeometry/variety_affine_linear.jl b/src/TropicalGeometry/variety_affine_linear.jl index dc78caafd5f9..add0a3f48ea4 100644 --- a/src/TropicalGeometry/variety_affine_linear.jl +++ b/src/TropicalGeometry/variety_affine_linear.jl @@ -16,7 +16,7 @@ function tropical_variety_affine_linear(I::MPolyIdeal,nu::TropicalSemiringMap; w end # input inhomogeneous, homogenise first - Ih = homogenize_pre_tropicalization(I) + Ih,_,_ = homogenize_pre_tropicalization(I) TropLh = tropical_linear_space(Ih,nu,weighted_polyhedral_complex_only=true) return dehomogenize_post_tropicalization(TropLh) end diff --git a/src/TropicalGeometry/variety_prime.jl b/src/TropicalGeometry/variety_prime.jl index 5670492db3d5..8ac648c2f100 100644 --- a/src/TropicalGeometry/variety_prime.jl +++ b/src/TropicalGeometry/variety_prime.jl @@ -14,7 +14,7 @@ function tropical_variety_prime(I::MPolyIdeal, nu::TropicalSemiringMap; weighted return tropical_variety_prime_singular(I,nu; weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) end - Ih = homogenize_pre_tropicalization(I) + Ih,_,_ = homogenize_pre_tropicalization(I) TropIh = tropical_variety_prime_singular(Ih,nu,weighted_polyhedral_complex_only=weighted_polyhedral_complex_only) return dehomogenize_post_tropicalization(TropIh) end