From a827b6403b1ae5ca72e12ec4ccefa2e8ec83e4ca Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Fri, 14 Jul 2023 22:02:17 +0200 Subject: [PATCH 01/21] introduced scalar_type_or_field --- src/PolyhedralGeometry/Cone/constructors.jl | 58 ++++++++-------- .../PolyhedralComplex/constructors.jl | 60 ++++++++-------- .../PolyhedralFan/constructors.jl | 68 +++++++++---------- .../Polyhedron/constructors.jl | 62 ++++++++--------- .../Polyhedron/standard_constructions.jl | 12 ++-- .../SubdivisionOfPoints/constructors.jl | 38 +++++------ src/PolyhedralGeometry/helpers.jl | 2 + src/PolyhedralGeometry/iterators.jl | 32 ++++----- src/PolyhedralGeometry/linear_program.jl | 4 +- 9 files changed, 169 insertions(+), 167 deletions(-) diff --git a/src/PolyhedralGeometry/Cone/constructors.jl b/src/PolyhedralGeometry/Cone/constructors.jl index ade7eae79831..31b0d06d989a 100644 --- a/src/PolyhedralGeometry/Cone/constructors.jl +++ b/src/PolyhedralGeometry/Cone/constructors.jl @@ -62,25 +62,25 @@ julia> HS = positive_hull(R, L) Polyhedral cone in ambient dimension 2 ``` """ -function positive_hull(f::Union{Type{T}, Field}, R::AbstractCollection[RayVector], L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, R, L) - inputrays = remove_zero_rows(unhomogenized_matrix(R)) - if isnothing(L) || isempty(L) - L = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, _ambient_dim(R)) - end - - if non_redundant - return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(RAYS = inputrays, LINEALITY_SPACE = unhomogenized_matrix(L),), parent_field) - else - return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(INPUT_RAYS = inputrays, INPUT_LINEALITY = unhomogenized_matrix(L),), parent_field) - end +function positive_hull(f::scalar_type_or_field, R::AbstractCollection[RayVector], L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) + parent_field, scalar_type = _determine_parent_and_scalar(f, R, L) + inputrays = remove_zero_rows(unhomogenized_matrix(R)) + if isnothing(L) || isempty(L) + L = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, _ambient_dim(R)) + end + + if non_redundant + return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(RAYS = inputrays, LINEALITY_SPACE = unhomogenized_matrix(L),), parent_field) + else + return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(INPUT_RAYS = inputrays, INPUT_LINEALITY = unhomogenized_matrix(L),), parent_field) + end end # Redirect everything to the above constructor, use QQFieldElem as default for the # scalar type T. positive_hull(R::AbstractCollection[RayVector], L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) = positive_hull(QQFieldElem, R, L; non_redundant=non_redundant) cone(R::AbstractCollection[RayVector], L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) = positive_hull(QQFieldElem, R, L; non_redundant=non_redundant) -cone(f::Union{Type{T}, Field}, R::AbstractCollection[RayVector], L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) where T<:scalar_types = positive_hull(f, R, L; non_redundant=non_redundant) -cone(f::Union{Type{T}, Field}, x...) where T<:scalar_types = positive_hull(f, x...) +cone(f::scalar_type_or_field, R::AbstractCollection[RayVector], L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) = positive_hull(f, R, L; non_redundant=non_redundant) +cone(f::scalar_type_or_field, x...) = positive_hull(f, x...) function ==(C0::Cone{T}, C1::Cone{T}) where T<:scalar_types @@ -115,16 +115,16 @@ julia> rays(C) [1, 1] ``` """ -function cone_from_inequalities(f::Union{Type{T}, Field}, I::AbstractCollection[LinearHalfspace], E::Union{Nothing, AbstractCollection[LinearHyperplane]} = nothing; non_redundant::Bool = false) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, I, E) - IM = -linear_matrix_for_polymake(I) - EM = isnothing(E) || isempty(E) ? Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(IM, 2)) : linear_matrix_for_polymake(E) - - if non_redundant - return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(FACETS = IM, LINEAR_SPAN = EM), parent_field) - else - return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(INEQUALITIES = IM, EQUATIONS = EM), parent_field) - end +function cone_from_inequalities(f::scalar_type_or_field, I::AbstractCollection[LinearHalfspace], E::Union{Nothing, AbstractCollection[LinearHyperplane]} = nothing; non_redundant::Bool = false) + parent_field, scalar_type = _determine_parent_and_scalar(f, I, E) + IM = -linear_matrix_for_polymake(I) + EM = isnothing(E) || isempty(E) ? Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(IM, 2)) : linear_matrix_for_polymake(E) + + if non_redundant + return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(FACETS = IM, LINEAR_SPAN = EM), parent_field) + else + return Cone{scalar_type}(Polymake.polytope.Cone{_scalar_type_to_polymake(scalar_type)}(INEQUALITIES = IM, EQUATIONS = EM), parent_field) + end end @doc raw""" @@ -152,11 +152,11 @@ julia> dim(C) 1 ``` """ -function cone_from_equations(f::Union{Type{T}, Field}, E::AbstractCollection[LinearHyperplane]; non_redundant::Bool = false) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, E) - EM = linear_matrix_for_polymake(E) - IM = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(EM, 2)) - return cone_from_inequalities(f, IM, EM; non_redundant = non_redundant) +function cone_from_equations(f::scalar_type_or_field, E::AbstractCollection[LinearHyperplane]; non_redundant::Bool = false) + parent_field, scalar_type = _determine_parent_and_scalar(f, E) + EM = linear_matrix_for_polymake(E) + IM = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(EM, 2)) + return cone_from_inequalities(f, IM, EM; non_redundant = non_redundant) end cone_from_inequalities(x...) = cone_from_inequalities(QQFieldElem, x...) diff --git a/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl b/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl index a0bd21e680e8..e508b354c753 100644 --- a/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl +++ b/src/PolyhedralGeometry/PolyhedralComplex/constructors.jl @@ -75,36 +75,36 @@ julia> lineality_dim(PC) 1 ``` """ -function polyhedral_complex(f::Union{Type{T}, Field}, +function polyhedral_complex(f::scalar_type_or_field, polyhedra::IncidenceMatrix, vr::AbstractCollection[PointVector], far_vertices::Union{Vector{Int}, Nothing} = nothing, L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false - ) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, vr, L) - points = homogenized_matrix(vr, 1) - LM = isnothing(L) || isempty(L) ? zero_matrix(QQ, 0, size(points, 2)) : homogenized_matrix(L, 0) - - # Rays and Points are homogenized and combined and - # If some vertices are far vertices, give them a leading 0 - if !isnothing(far_vertices) - points[far_vertices,1] .= 0 - end - - if non_redundant - return PolyhedralComplex{scalar_type}(Polymake.fan.PolyhedralComplex{_scalar_type_to_polymake(scalar_type)}( - VERTICES = points, - LINEALITY_SPACE = LM, - MAXIMAL_CONES = polyhedra, - ), parent_field) - else - return PolyhedralComplex{scalar_type}(Polymake.fan.PolyhedralComplex{_scalar_type_to_polymake(scalar_type)}( - POINTS = points, - INPUT_LINEALITY = LM, - INPUT_CONES = polyhedra, - ), parent_field) - end + ) + parent_field, scalar_type = _determine_parent_and_scalar(f, vr, L) + points = homogenized_matrix(vr, 1) + LM = isnothing(L) || isempty(L) ? zero_matrix(QQ, 0, size(points, 2)) : homogenized_matrix(L, 0) + + # Rays and Points are homogenized and combined and + # If some vertices are far vertices, give them a leading 0 + if !isnothing(far_vertices) + points[far_vertices,1] .= 0 + end + + if non_redundant + return PolyhedralComplex{scalar_type}(Polymake.fan.PolyhedralComplex{_scalar_type_to_polymake(scalar_type)}( + VERTICES = points, + LINEALITY_SPACE = LM, + MAXIMAL_CONES = polyhedra, + ), parent_field) + else + return PolyhedralComplex{scalar_type}(Polymake.fan.PolyhedralComplex{_scalar_type_to_polymake(scalar_type)}( + POINTS = points, + INPUT_LINEALITY = LM, + INPUT_CONES = polyhedra, + ), parent_field) + end end # default scalar type: `QQFieldElem` polyhedral_complex(polyhedra::IncidenceMatrix, @@ -114,11 +114,11 @@ polyhedral_complex(polyhedra::IncidenceMatrix, non_redundant::Bool = false) = polyhedral_complex(QQFieldElem, polyhedra, vr, far_vertices, L; non_redundant=non_redundant) -function polyhedral_complex(f::Union{Type{T}, Field}, v::AbstractCollection[PointVector], vi::IncidenceMatrix, r::AbstractCollection[RayVector], ri::IncidenceMatrix, L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) where T<:scalar_types - vr = [unhomogenized_matrix(v); unhomogenized_matrix(r)] - far_vertices = collect((size(v, 1) + 1):size(vr, 1)) - polyhedra = hcat(vi, ri) - return polyhedral_complex(f, polyhedra, vr, far_vertices, L; non_redundant = non_redundant) +function polyhedral_complex(f::scalar_type_or_field, v::AbstractCollection[PointVector], vi::IncidenceMatrix, r::AbstractCollection[RayVector], ri::IncidenceMatrix, L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) + vr = [unhomogenized_matrix(v); unhomogenized_matrix(r)] + far_vertices = collect((size(v, 1) + 1):size(vr, 1)) + polyhedra = hcat(vi, ri) + return polyhedral_complex(f, polyhedra, vr, far_vertices, L; non_redundant = non_redundant) end # TODO: Only works for this specific case; implement generalization using `iter.Acc` diff --git a/src/PolyhedralGeometry/PolyhedralFan/constructors.jl b/src/PolyhedralGeometry/PolyhedralFan/constructors.jl index eaa018608ea9..0b6a16df06d5 100644 --- a/src/PolyhedralGeometry/PolyhedralFan/constructors.jl +++ b/src/PolyhedralGeometry/PolyhedralFan/constructors.jl @@ -64,33 +64,33 @@ julia> lineality_dim(PF) 1 ``` """ -function polyhedral_fan(f::Union{Type{T}, Field}, +function polyhedral_fan(f::scalar_type_or_field, Rays::AbstractCollection[RayVector], LS::Union{AbstractCollection[RayVector], Nothing}, Incidence::IncidenceMatrix; - non_redundant::Bool = false) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, Rays, LS) - RM = unhomogenized_matrix(Rays) - if isnothing(LS) - LM = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(RM, 2)) - else - LM = unhomogenized_matrix(LS) - end - if non_redundant - return PolyhedralFan{scalar_type}(Polymake.fan.PolyhedralFan{_scalar_type_to_polymake(scalar_type)}( - RAYS = RM, - LINEALITY_SPACE = LM, - MAXIMAL_CONES = Incidence, - ), parent_field) - else - return PolyhedralFan{scalar_type}(Polymake.fan.PolyhedralFan{_scalar_type_to_polymake(scalar_type)}( - INPUT_RAYS = RM, - INPUT_LINEALITY = LM, - INPUT_CONES = Incidence, - ), parent_field) - end + non_redundant::Bool = false) + parent_field, scalar_type = _determine_parent_and_scalar(f, Rays, LS) + RM = unhomogenized_matrix(Rays) + if isnothing(LS) + LM = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(RM, 2)) + else + LM = unhomogenized_matrix(LS) + end + if non_redundant + return PolyhedralFan{scalar_type}(Polymake.fan.PolyhedralFan{_scalar_type_to_polymake(scalar_type)}( + RAYS = RM, + LINEALITY_SPACE = LM, + MAXIMAL_CONES = Incidence, + ), parent_field) + else + return PolyhedralFan{scalar_type}(Polymake.fan.PolyhedralFan{_scalar_type_to_polymake(scalar_type)}( + INPUT_RAYS = RM, + INPUT_LINEALITY = LM, + INPUT_CONES = Incidence, + ), parent_field) + end end -polyhedral_fan(f::Union{Type{T}, Field}, Rays::AbstractCollection[RayVector], Incidence::IncidenceMatrix; non_redundant::Bool = false) where T<:scalar_types = polyhedral_fan(f, Rays, nothing, Incidence; non_redundant=non_redundant) +polyhedral_fan(f::scalar_type_or_field, Rays::AbstractCollection[RayVector], Incidence::IncidenceMatrix; non_redundant::Bool = false) = polyhedral_fan(f, Rays, nothing, Incidence; non_redundant=non_redundant) polyhedral_fan(Rays::AbstractCollection[RayVector], LS::Union{AbstractCollection[RayVector], Nothing}, Incidence::IncidenceMatrix; non_redundant::Bool = false) = polyhedral_fan(QQFieldElem, Rays, LS, Incidence; non_redundant = non_redundant) polyhedral_fan(Rays::AbstractCollection[RayVector], Incidence::IncidenceMatrix; non_redundant::Bool = false) = polyhedral_fan(QQFieldElem, Rays, Incidence; non_redundant = non_redundant) @@ -104,9 +104,9 @@ pm_object(PF::PolyhedralFan) = PF.pm_fan polyhedral_fan(itr::AbstractVector{Cone{T}}) where T<:scalar_types = PolyhedralFan{T}(Polymake.fan.check_fan_objects(pm_object.(itr)...), coefficient_field(iterate(itr)[1])) #Same construction for when the user gives Matrix{Bool} as incidence matrix -polyhedral_fan(f::Union{Type{T}, Field}, Rays::AbstractCollection[RayVector], LS::AbstractCollection[RayVector], Incidence::Matrix{Bool}) where T<:scalar_types = +polyhedral_fan(f::scalar_type_or_field, Rays::AbstractCollection[RayVector], LS::AbstractCollection[RayVector], Incidence::Matrix{Bool}) = polyhedral_fan(f, Rays, LS, IncidenceMatrix(Polymake.IncidenceMatrix(Incidence))) -polyhedral_fan(f::Union{Type{T}, Field}, Rays::AbstractCollection[RayVector], Incidence::Matrix{Bool}) where T<:scalar_types = +polyhedral_fan(f::scalar_type_or_field, Rays::AbstractCollection[RayVector], Incidence::Matrix{Bool}) = polyhedral_fan(f, Rays, IncidenceMatrix(Polymake.IncidenceMatrix(Incidence))) @@ -135,15 +135,15 @@ parent `Field`. group acting on the rays of the fan. """ -function polyhedral_fan_from_rays_action(f::Union{Type{T}, Field}, Rays::AbstractCollection[RayVector], MC_reps::IncidenceMatrix, perms::AbstractVector{PermGroupElem}) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, Rays) - pf = Polymake.fan.PolyhedralFan{_scalar_type_to_polymake(scalar_type)}() - Polymake.take(pf, "RAYS", Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(unhomogenized_matrix(Rays))) - d = length(Rays) - gp = _group_generators_to_pm_arr_arr(perms, d) - Polymake.take(pf, "GROUP.RAYS_ACTION.GENERATORS", gp) - Polymake.take(pf, "GROUP.MAXIMAL_CONES_ACTION.MAXIMAL_CONES_GENERATORS", MC_reps) - return PolyhedralFan{scalar_type}(pf, parent_field) +function polyhedral_fan_from_rays_action(f::scalar_type_or_field, Rays::AbstractCollection[RayVector], MC_reps::IncidenceMatrix, perms::AbstractVector{PermGroupElem}) + parent_field, scalar_type = _determine_parent_and_scalar(f, Rays) + pf = Polymake.fan.PolyhedralFan{_scalar_type_to_polymake(scalar_type)}() + Polymake.take(pf, "RAYS", Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(unhomogenized_matrix(Rays))) + d = length(Rays) + gp = _group_generators_to_pm_arr_arr(perms, d) + Polymake.take(pf, "GROUP.RAYS_ACTION.GENERATORS", gp) + Polymake.take(pf, "GROUP.MAXIMAL_CONES_ACTION.MAXIMAL_CONES_GENERATORS", MC_reps) + return PolyhedralFan{scalar_type}(pf, parent_field) end polyhedral_fan_from_rays_action(Rays::AbstractCollection[RayVector], MC_reps::IncidenceMatrix, perms::AbstractVector{PermGroupElem}) = polyhedral_fan_from_rays_action(QQFieldElem, Rays, MC_reps, perms) diff --git a/src/PolyhedralGeometry/Polyhedron/constructors.jl b/src/PolyhedralGeometry/Polyhedron/constructors.jl index 45250c090245..0f2bf83b2a90 100644 --- a/src/PolyhedralGeometry/Polyhedron/constructors.jl +++ b/src/PolyhedralGeometry/Polyhedron/constructors.jl @@ -57,17 +57,17 @@ julia> polyhedron(A,b) Polyhedron in ambient dimension 2 ``` """ -polyhedron(f::Union{Type{T}, Field}, A::AnyVecOrMat, b::AbstractVector) where T<:scalar_types = polyhedron(f, (A, b)) +polyhedron(f::scalar_type_or_field, A::AnyVecOrMat, b::AbstractVector) = polyhedron(f, (A, b)) -polyhedron(f::Union{Type{T}, Field}, A::AbstractVector, b::Any) where T<:scalar_types = polyhedron(f, ([A], [b])) +polyhedron(f::scalar_type_or_field, A::AbstractVector, b::Any) = polyhedron(f, ([A], [b])) -polyhedron(f::Union{Type{T}, Field}, A::AbstractVector, b::AbstractVector) where T<:scalar_types = polyhedron(f, ([A], b)) +polyhedron(f::scalar_type_or_field, A::AbstractVector, b::AbstractVector) = polyhedron(f, ([A], b)) -polyhedron(f::Union{Type{T}, Field}, A::AbstractVector{<:AbstractVector}, b::Any) where T<:scalar_types = polyhedron(f, (A, [b])) +polyhedron(f::scalar_type_or_field, A::AbstractVector{<:AbstractVector}, b::Any) = polyhedron(f, (A, [b])) -polyhedron(f::Union{Type{T}, Field}, A::AbstractVector{<:AbstractVector}, b::AbstractVector) where T<:scalar_types = polyhedron(f, (A, b)) +polyhedron(f::scalar_type_or_field, A::AbstractVector{<:AbstractVector}, b::AbstractVector) = polyhedron(f, (A, b)) -polyhedron(f::Union{Type{T}, Field}, A::AnyVecOrMat, b::Any) where T<:scalar_types = polyhedron(f, A, [b]) +polyhedron(f::scalar_type_or_field, A::AnyVecOrMat, b::Any) = polyhedron(f, A, [b]) @doc raw""" polyhedron(::Union{Type{T}, Field}, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) where T<:scalar_types @@ -106,17 +106,17 @@ julia> vertices(P) [0, 0] ``` """ -function polyhedron(f::Union{Type{T}, Field}, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing; parent_field::Union{Nothing, Field} = nothing) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, I, E) - if isnothing(I) || _isempty_halfspace(I) - EM = affine_matrix_for_polymake(E) - IM = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(EM, 2)) - else - IM = -affine_matrix_for_polymake(I) - EM = isnothing(E) || _isempty_halfspace(E) ? Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(IM, 2)) : affine_matrix_for_polymake(E) - end - - return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(INEQUALITIES = remove_zero_rows(IM), EQUATIONS = remove_zero_rows(EM)), parent_field) +function polyhedron(f::scalar_type_or_field, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing; parent_field::Union{Nothing, Field} = nothing) + parent_field, scalar_type = _determine_parent_and_scalar(f, I, E) + if isnothing(I) || _isempty_halfspace(I) + EM = affine_matrix_for_polymake(E) + IM = Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(EM, 2)) + else + IM = -affine_matrix_for_polymake(I) + EM = isnothing(E) || _isempty_halfspace(E) ? Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(IM, 2)) : affine_matrix_for_polymake(E) + end + + return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(INEQUALITIES = remove_zero_rows(IM), EQUATIONS = remove_zero_rows(EM)), parent_field) end """ @@ -198,20 +198,20 @@ julia> XA = convex_hull(V, R, L) Polyhedron in ambient dimension 2 ``` """ -function convex_hull(f::Union{Type{T}, Field}, V::AbstractCollection[PointVector], R::Union{AbstractCollection[RayVector], Nothing} = nothing, L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, V, R, L) - # Rays and Points are homogenized and combined and - # Lineality is homogenized - points = stack(homogenized_matrix(V, 1), homogenized_matrix(R, 0)) - lineality = isnothing(L) || isempty(L) ? zero_matrix(QQ, 0, size(points,2)) : homogenized_matrix(L, 0) - - # These matrices are in the right format for polymake. - # given non_redundant can avoid unnecessary redundancy checks - if non_redundant - return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(VERTICES = points, LINEALITY_SPACE = lineality), parent_field) - else - return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(POINTS = remove_zero_rows(points), INPUT_LINEALITY = remove_zero_rows(lineality)), parent_field) - end +function convex_hull(f::scalar_type_or_field, V::AbstractCollection[PointVector], R::Union{AbstractCollection[RayVector], Nothing} = nothing, L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) + parent_field, scalar_type = _determine_parent_and_scalar(f, V, R, L) + # Rays and Points are homogenized and combined and + # Lineality is homogenized + points = stack(homogenized_matrix(V, 1), homogenized_matrix(R, 0)) + lineality = isnothing(L) || isempty(L) ? zero_matrix(QQ, 0, size(points,2)) : homogenized_matrix(L, 0) + + # These matrices are in the right format for polymake. + # given non_redundant can avoid unnecessary redundancy checks + if non_redundant + return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(VERTICES = points, LINEALITY_SPACE = lineality), parent_field) + else + return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(POINTS = remove_zero_rows(points), INPUT_LINEALITY = remove_zero_rows(lineality)), parent_field) + end end convex_hull(V::AbstractCollection[PointVector], R::Union{AbstractCollection[RayVector], Nothing} = nothing, L::Union{AbstractCollection[RayVector], Nothing} = nothing; non_redundant::Bool = false) = convex_hull(QQFieldElem, V, R, L; non_redundant=non_redundant) diff --git a/src/PolyhedralGeometry/Polyhedron/standard_constructions.jl b/src/PolyhedralGeometry/Polyhedron/standard_constructions.jl index 3dbab25fb135..e5edf65c7279 100644 --- a/src/PolyhedralGeometry/Polyhedron/standard_constructions.jl +++ b/src/PolyhedralGeometry/Polyhedron/standard_constructions.jl @@ -206,14 +206,14 @@ julia> normalized_volume(C) 120 ``` """ -function cube(f::Union{Type{T},Field}, d::Int) where {T<:scalar_types} +function cube(f::scalar_type_or_field, d::Int) parent_field, scalar_type = _determine_parent_and_scalar(f) return Polyhedron{scalar_type}( Polymake.polytope.cube{_scalar_type_to_polymake(scalar_type)}(d), parent_field ) end cube(d::Int) = cube(QQFieldElem, d) -function cube(f::Union{Type{T},Field}, d::Int, l, u) where {T<:scalar_types} +function cube(f::scalar_type_or_field, d::Int, l, u) parent_field, scalar_type = _determine_parent_and_scalar(f, l, u) return Polyhedron{scalar_type}( Polymake.polytope.cube{_scalar_type_to_polymake(scalar_type)}(d, u, l), parent_field @@ -622,14 +622,14 @@ julia> facets(t) x₁ + x₂ + x₃ + x₄ + x₅ + x₆ + x₇ ≦ 5 ``` """ -function simplex(f::Union{Type{T},Field}, d::Int, n) where {T<:scalar_types} +function simplex(f::scalar_type_or_field, d::Int, n) parent_field, scalar_type = _determine_parent_and_scalar(f, n) return Polyhedron{scalar_type}( Polymake.polytope.simplex{_scalar_type_to_polymake(scalar_type)}(d, n), parent_field ) end simplex(d::Int, n) = simplex(QQFieldElem, d, n) -function simplex(f::Union{Type{T},Field}, d::Int) where {T<:scalar_types} +function simplex(f::scalar_type_or_field, d::Int) parent_field, scalar_type = _determine_parent_and_scalar(f) return Polyhedron{scalar_type}( Polymake.polytope.simplex{_scalar_type_to_polymake(scalar_type)}(d), parent_field @@ -678,14 +678,14 @@ x₁ - x₂ - x₃ ≦ 2 -x₁ - x₂ - x₃ ≦ 2 ``` """ -function cross_polytope(f::Union{Type{T},Field}, d::Int64, n) where {T<:scalar_types} +function cross_polytope(f::scalar_type_or_field, d::Int64, n) parent_field, scalar_type = _determine_parent_and_scalar(f, n) return Polyhedron{scalar_type}( Polymake.polytope.cross{_scalar_type_to_polymake(scalar_type)}(d, n), parent_field ) end cross_polytope(d::Int64, n) = cross_polytope(QQFieldElem, d, n) -function cross_polytope(f::Union{Type{T},Field}, d::Int64) where {T<:scalar_types} +function cross_polytope(f::scalar_type_or_field, d::Int64) parent_field, scalar_type = _determine_parent_and_scalar(f) return Polyhedron{scalar_type}( Polymake.polytope.cross{_scalar_type_to_polymake(scalar_type)}(d), parent_field diff --git a/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl b/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl index 9a314dd42c03..ab6e22d7c9db 100644 --- a/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl +++ b/src/PolyhedralGeometry/SubdivisionOfPoints/constructors.jl @@ -53,14 +53,14 @@ julia> MOAE = subdivision_of_points(moaepts, moaeimnonreg0) Subdivision of points in ambient dimension 3 ``` """ -function subdivision_of_points(f::Union{Type{T}, Field}, points::AbstractCollection[PointVector], cells::IncidenceMatrix) where T<:scalar_types - @req size(points)[1] == ncols(cells) "Number of points must be the same as columns of IncidenceMatrix" - parent_field, scalar_type = _determine_parent_and_scalar(f, points) - arr = @Polymake.convert_to Array{Set{Int}} Polymake.common.rows(cells) - SubdivisionOfPoints{scalar_type}(Polymake.fan.SubdivisionOfPoints{_scalar_type_to_polymake(scalar_type)}( - POINTS = homogenize(points,1), - MAXIMAL_CELLS = arr, - ), parent_field) +function subdivision_of_points(f::scalar_type_or_field, points::AbstractCollection[PointVector], cells::IncidenceMatrix) + @req size(points)[1] == ncols(cells) "Number of points must be the same as columns of IncidenceMatrix" + parent_field, scalar_type = _determine_parent_and_scalar(f, points) + arr = @Polymake.convert_to Array{Set{Int}} Polymake.common.rows(cells) + SubdivisionOfPoints{scalar_type}(Polymake.fan.SubdivisionOfPoints{_scalar_type_to_polymake(scalar_type)}( + POINTS = homogenize(points,1), + MAXIMAL_CELLS = arr, + ), parent_field) end @@ -91,13 +91,13 @@ julia> n_maximal_cells(SOP) 1 ``` """ -function subdivision_of_points(f::Union{Type{T}, Field}, points::AbstractCollection[PointVector], weights::AbstractVector) where T<:scalar_types - @req size(points)[1] == length(weights) "Number of points must equal number of weights" - parent_field, scalar_type = _determine_parent_and_scalar(f, points, weights) - SubdivisionOfPoints{scalar_type}(Polymake.fan.SubdivisionOfPoints{_scalar_type_to_polymake(scalar_type)}( - POINTS = homogenize(points,1), - WEIGHTS = weights, - ), parent_field) +function subdivision_of_points(f::scalar_type_or_field, points::AbstractCollection[PointVector], weights::AbstractVector) + @req size(points)[1] == length(weights) "Number of points must equal number of weights" + parent_field, scalar_type = _determine_parent_and_scalar(f, points, weights) + SubdivisionOfPoints{scalar_type}(Polymake.fan.SubdivisionOfPoints{_scalar_type_to_polymake(scalar_type)}( + POINTS = homogenize(points,1), + WEIGHTS = weights, + ), parent_field) end """ @@ -109,13 +109,13 @@ pm_object(SOP::SubdivisionOfPoints) = SOP.pm_subdivision #Same construction for when the user provides maximal cells -function subdivision_of_points(::Union{Type{T}, Field}, points::AbstractCollection[PointVector], cells::Vector{Vector{Int64}}) where T<:scalar_types - subdivision_of_points(T, points, IncidenceMatrix(cells)) +function subdivision_of_points(f::scalar_type_or_field, points::AbstractCollection[PointVector], cells::Vector{Vector{Int64}}) + subdivision_of_points(f, points, IncidenceMatrix(cells)) end #Same construction for when the user gives Matrix{Bool} as incidence matrix -function subdivision_of_points(::Union{Type{T}, Field}, Points::Union{Oscar.MatElem,AbstractMatrix}, cells::Matrix{Bool}) where T<:scalar_types - subdivision_of_points(T, points, IncidenceMatrix(Polymake.IncidenceMatrix(cells))) +function subdivision_of_points(f::scalar_type_or_field, Points::Union{Oscar.MatElem,AbstractMatrix}, cells::Matrix{Bool}) + subdivision_of_points(f, points, IncidenceMatrix(Polymake.IncidenceMatrix(cells))) end ############################################################################### diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index 2ffedf4c83ac..cbc05391ef42 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -386,6 +386,8 @@ const scalar_type_to_oscar = Dict{String, Type}([("Rational", QQFieldElem), const scalar_types_extended = Union{scalar_types, ZZRingElem} +const scalar_type_or_field = Union{Type{<:scalar_types}, Field} + _scalar_type_to_polymake(::Type{QQFieldElem}) = Polymake.Rational _scalar_type_to_polymake(::Type{<:FieldElem}) = Polymake.OscarNumber _scalar_type_to_polymake(::Type{Float64}) = Float64 diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 5779e161ee60..126e1f8abb4f 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -115,13 +115,13 @@ struct AffineHalfspace{T} <: Halfspace{T} end halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_halfspace(a, b) -halfspace(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) where T<:scalar_types = affine_halfspace(f, a, b) +halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_halfspace(f, a, b) invert(H::AffineHalfspace{T}) where T<:scalar_types = AffineHalfspace{T}(coefficient_field(H), -normal_vector(H), -negbias(H)) -function affine_halfspace(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) - return AffineHalfspace{scalar_type}(parent_field, a, b) +function affine_halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) + parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) + return AffineHalfspace{scalar_type}(parent_field, a, b) end affine_halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) = affine_halfspace(QQ, a, b) @@ -136,13 +136,13 @@ struct LinearHalfspace{T} <: Halfspace{T} end halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_halfspace(a) -halfspace(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}) where T<:scalar_types = linear_halfspace(f, a) +halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_halfspace(f, a) invert(H::LinearHalfspace{T}) where T<:scalar_types = LinearHalfspace{T}(coefficient_field(H), -normal_vector(H)) -function linear_halfspace(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, a) - return LinearHalfspace{scalar_type}(parent_field, a) +function linear_halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) + parent_field, scalar_type = _determine_parent_and_scalar(f, a) + return LinearHalfspace{scalar_type}(parent_field, a) end linear_halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_halfspace(QQ, a) @@ -170,11 +170,11 @@ struct AffineHyperplane{T} <: Hyperplane{T} end hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_hyperplane(a, b) -hyperplane(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) where T<:scalar_types = affine_hyperplane(f, a, b) +hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_hyperplane(f, a, b) -function affine_hyperplane(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) - return AffineHyperplane{scalar_type}(parent_field, a, b) +function affine_hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) + parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) + return AffineHyperplane{scalar_type}(parent_field, a, b) end affine_hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) = affine_hyperplane(QQ, a, b) @@ -189,11 +189,11 @@ struct LinearHyperplane{T} <: Hyperplane{T} end hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_hyperplane(a) -hyperplane(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}) where T<:scalar_types = linear_hyperplane(f, a) +hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_hyperplane(f, a) -function linear_hyperplane(f::Union{Type{T}, Field}, a::Union{MatElem, AbstractMatrix, AbstractVector}) where T<:scalar_types - parent_field, scalar_type = _determine_parent_and_scalar(f, a) - return LinearHyperplane{scalar_type}(parent_field, a) +function linear_hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) + parent_field, scalar_type = _determine_parent_and_scalar(f, a) + return LinearHyperplane{scalar_type}(parent_field, a) end linear_hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_hyperplane(QQ, a) diff --git a/src/PolyhedralGeometry/linear_program.jl b/src/PolyhedralGeometry/linear_program.jl index dc0ee8882d4d..7497c6658229 100644 --- a/src/PolyhedralGeometry/linear_program.jl +++ b/src/PolyhedralGeometry/linear_program.jl @@ -44,13 +44,13 @@ function linear_program( end linear_program( - f::Union{Type{T},Field}, + f::scalar_type_or_field, A::AbstractCollection[AffineHalfspace], b, c::AbstractVector; k=0, convention=:max, -) where {T<:scalar_types} = +) = linear_program(polyhedron(f, A, b), c; k=k, convention=convention) pm_object(lp::LinearProgram) = lp.polymake_lp From 8b0f0bcc1bd22149be39ea0025fffca4329317b5 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Mon, 17 Jul 2023 18:24:03 +0200 Subject: [PATCH 02/21] re-factor for PointVector --- src/PolyhedralGeometry/Cone/properties.jl | 4 +- .../PolyhedralComplex/properties.jl | 18 +- .../Polyhedron/properties.jl | 10 +- .../SubdivisionOfPoints/properties.jl | 2 +- src/PolyhedralGeometry/helpers.jl | 18 +- src/PolyhedralGeometry/iterators.jl | 44 ++-- src/PolyhedralGeometry/linear_program.jl | 2 +- .../mixed_integer_linear_program.jl | 2 +- src/exports.jl | 1 + test/PolyhedralGeometry/types.jl | 200 +++++++++--------- 10 files changed, 156 insertions(+), 145 deletions(-) diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index c432b40358d5..06230cf55750 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -530,7 +530,7 @@ function hilbert_basis(C::Cone{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(C, _hilbert_generator, size(pm_object(C).HILBERT_BASIS_GENERATORS[1], 1)) end -_hilbert_generator(::Type{PointVector{ZZRingElem}}, C::Cone, i::Base.Integer) = PointVector{ZZRingElem}(view(pm_object(C).HILBERT_BASIS_GENERATORS[1], i, :)) +_hilbert_generator(T::Type{PointVector{ZZRingElem}}, C::Cone, i::Base.Integer) = point_vector(ZZ, view(pm_object(C).HILBERT_BASIS_GENERATORS[1], i, :))::T _generator_matrix(::Val{_hilbert_generator}, C::Cone; homogenized=false) = homogenized ? homogenize(pm_object(C).HILBERT_BASIS_GENERATORS[1], 0) : pm_object(C).HILBERT_BASIS_GENERATORS[1] @@ -586,4 +586,4 @@ Base.in(v::AbstractVector, C::Cone) = Polymake.polytope.contains(pm_object(C), v Compute a point in the relative interior point of `C`, i.e. a point in `C` not contained in any facet. """ -relative_interior_point(C::Cone{T}) where T<:scalar_types = PointVector{T}(coefficient_field(C), view(Polymake.common.dense(pm_object(C).REL_INT_POINT), :)) # broadcast_view +relative_interior_point(C::Cone{T}) where T<:scalar_types = point_vector(coefficient_field(C), view(Polymake.common.dense(pm_object(C).REL_INT_POINT), :))::PointVector{T} # broadcast_view diff --git a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl index 5d0c39d2ea3b..79835daa803c 100644 --- a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl @@ -63,12 +63,10 @@ vertices(as::Type{PointVector{T}}, PC::PolyhedralComplex{T}) where {T<:scalar_ty _vertices(as::Type{PointVector{T}}, PC::PolyhedralComplex) where {T<:scalar_types} = SubObjectIterator{as}(PC, _vertex_complex, length(_vertex_indices(pm_object(PC)))) -_vertex_complex( - ::Type{PointVector{T}}, PC::PolyhedralComplex, i::Base.Integer -) where {T<:scalar_types} = PointVector{T}( +_vertex_complex(U::Type{PointVector{T}}, PC::PolyhedralComplex{T}, i::Base.Integer) where {T<:scalar_types} = point_vector( coefficient_field(PC), @view pm_object(PC).VERTICES[_vertex_indices(pm_object(PC))[i], 2:end] -) +)::U _point_matrix(::Val{_vertex_complex}, PC::PolyhedralComplex; homogenized=false) = @view pm_object(PC).VERTICES[_vertex_indices(pm_object(PC)), (homogenized ? 1 : 2):end] @@ -85,20 +83,18 @@ function _all_vertex_indices(P::Polymake.BigObject) return vi end -function _vertex_or_ray_complex( - ::Type{Union{PointVector{T},RayVector{T}}}, PC::PolyhedralComplex, i::Base.Integer -) where {T<:scalar_types} +function _vertex_or_ray_complex(::Type{Union{PointVector{T}, RayVector{T}}}, PC::PolyhedralComplex{T}, i::Base.Integer) where T<:scalar_types A = pm_object(PC).VERTICES if iszero(A[_all_vertex_indices(pm_object(PC))[i], 1]) - return RayVector{T}( + return ray_vector( coefficient_field(PC), @view pm_object(PC).VERTICES[_all_vertex_indices(pm_object(PC))[i], 2:end] - ) + )::RayVector{T} else - return PointVector{T}( + return point_vector( coefficient_field(PC), @view pm_object(PC).VERTICES[_all_vertex_indices(pm_object(PC))[i], 2:end] - ) + )::PointVector{T} end end diff --git a/src/PolyhedralGeometry/Polyhedron/properties.jl b/src/PolyhedralGeometry/Polyhedron/properties.jl index 33330dc1ddca..396ac8185f76 100644 --- a/src/PolyhedralGeometry/Polyhedron/properties.jl +++ b/src/PolyhedralGeometry/Polyhedron/properties.jl @@ -200,7 +200,7 @@ julia> vertices(PointVector, P) vertices(as::Type{PointVector{T}}, P::Polyhedron) where T<:scalar_types = lineality_dim(P) == 0 ? _vertices(as, P) : _empty_subobjectiterator(as, P) _vertices(as::Type{PointVector{T}}, P::Polyhedron) where T<:scalar_types = SubObjectIterator{as}(P, _vertex_polyhedron, length(_vertex_indices(pm_object(P)))) -_vertex_polyhedron(::Type{PointVector{T}}, P::Polyhedron, i::Base.Integer) where T<:scalar_types = PointVector{T}(coefficient_field(P).(@view pm_object(P).VERTICES[_vertex_indices(pm_object(P))[i], 2:end])) +_vertex_polyhedron(U::Type{PointVector{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types = point_vector(coefficient_field(P), @view pm_object(P).VERTICES[_vertex_indices(pm_object(P))[i], 2:end])::U _point_matrix(::Val{_vertex_polyhedron}, P::Polyhedron; homogenized=false) = @view pm_object(P).VERTICES[_vertex_indices(pm_object(P)), (homogenized ? 1 : 2):end] @@ -615,7 +615,7 @@ function lattice_points(P::Polyhedron{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(P, _lattice_point, size(pm_object(P).LATTICE_POINTS_GENERATORS[1], 1)) end -_lattice_point(::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = PointVector{ZZRingElem}(@view pm_object(P).LATTICE_POINTS_GENERATORS[1][i, 2:end]) +_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).LATTICE_POINTS_GENERATORS[1][i, 2:end])::T _point_matrix(::Val{_lattice_point}, P::Polyhedron; homogenized=false) = @view pm_object(P).LATTICE_POINTS_GENERATORS[1][:, (homogenized ? 1 : 2):end] @@ -646,7 +646,7 @@ function interior_lattice_points(P::Polyhedron{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(P, _interior_lattice_point, size(pm_object(P).INTERIOR_LATTICE_POINTS, 1)) end -_interior_lattice_point(::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = PointVector{ZZRingElem}(@view pm_object(P).INTERIOR_LATTICE_POINTS[i, 2:end]) +_interior_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).INTERIOR_LATTICE_POINTS[i, 2:end])::T _point_matrix(::Val{_interior_lattice_point}, P::Polyhedron; homogenized=false) = homogenized ? pm_object(P).INTERIOR_LATTICE_POINTS : @view pm_object(P).INTERIOR_LATTICE_POINTS[:, 2:end] @@ -686,7 +686,7 @@ function boundary_lattice_points(P::Polyhedron{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(P, _boundary_lattice_point, size(pm_object(P).BOUNDARY_LATTICE_POINTS, 1)) end -_boundary_lattice_point(::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = PointVector{ZZRingElem}(@view pm_object(P).BOUNDARY_LATTICE_POINTS[i, 2:end]) +_boundary_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).BOUNDARY_LATTICE_POINTS[i, 2:end])::T _point_matrix(::Val{_boundary_lattice_point}, P::Polyhedron; homogenized=false) = homogenized ? pm_object(P).BOUNDARY_LATTICE_POINTS : @view pm_object(P).BOUNDARY_LATTICE_POINTS[:, 2:end] @@ -1213,7 +1213,7 @@ julia> matrix(QQ, vertices(square)) [ 1 1] ``` """ -relative_interior_point(P::Polyhedron{T}) where T<:scalar_types = PointVector{T}(coefficient_field(P).(view(dehomogenize(Polymake.common.dense(pm_object(P).REL_INT_POINT)), :))) #view_broadcast +relative_interior_point(P::Polyhedron{T}) where T<:scalar_types = point_vector(coefficient_field(P), view(dehomogenize(Polymake.common.dense(pm_object(P).REL_INT_POINT)), :))::PointVector{T} #view_broadcast @doc raw""" diff --git a/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl b/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl index be7bc2d2e111..76d2e8f176e8 100644 --- a/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl +++ b/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl @@ -33,7 +33,7 @@ function points(SOP::SubdivisionOfPoints) return SubObjectIterator{PointVector{QQFieldElem}}(SOP, _point, size(pm_object(SOP).POINTS, 1)) end -_point(::Type{PointVector{QQFieldElem}}, SOP::SubdivisionOfPoints, i::Base.Integer) = PointVector{QQFieldElem}(pm_object(SOP).POINTS[i, 2:end]) +_point(T::Type{PointVector{QQFieldElem}}, SOP::SubdivisionOfPoints, i::Base.Integer) = point_vector(QQ, pm_object(SOP).POINTS[i, 2:end])::T _point_matrix(::Val{_point}, SOP::SubdivisionOfPoints; homogenized=false) = pm_object(SOP).POINTS[:, (homogenized ? 1 : 2):end] diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index cbc05391ef42..c94abdf20fe7 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -435,7 +435,7 @@ function _find_parent_field(::Type{T}, x::AbstractArray) where T <: scalar_types return QQ end -_determine_parent_and_scalar(f::Field, x...) = (f, elem_type(f)) +_determine_parent_and_scalar(f::Union{Field, ZZRing}, x...) = (f, elem_type(f)) # isempty(x) => standard/trivial field? function _determine_parent_and_scalar(::Type{T}, x...) where T <: scalar_types if T == QQFieldElem @@ -504,13 +504,13 @@ function _detect_scalar_and_field(::Type{U}, p::Polymake.BigObject) where U<:Pol end # promotion helpers -function _promote_scalar_field(f::Field...) - try - x = sum([g(0) for g in f]) - p = parent(x) - return (elem_type(p), p) - catch e - throw(ArgumentError("Can not find a mutual parent field for $f.")) +function _promote_scalar_field(f::Union{Field, ZZRing}...) + try + x = sum([g(0) for g in f]) + p = parent(x) + return (elem_type(p), p) + catch e + throw(ArgumentError("Can not find a mutual parent field for $f.")) end end @@ -519,6 +519,8 @@ function _promote_scalar_field(a::AbstractArray{<:FieldElem}) return _promote_scalar_field(parent.(a)...) end +_parent_or_coefficient_field(r::Base.RefValue{<:Union{FieldElem, ZZRingElem}}) = parent(r.x) + function _promoted_bigobject(::Type{T}, obj::PolyhedralObject{U}) where {T <: scalar_types, U <: scalar_types} T == U ? pm_object(obj) : Polymake.common.convert_to{_scalar_type_to_polymake(T)}(pm_object(obj)) end diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 126e1f8abb4f..231764d4a5fc 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -3,45 +3,55 @@ ################################################################################ struct PointVector{U} <: AbstractVector{U} - p::Vector{U} - - PointVector{U}(p::AbstractVector) where U<:scalar_types_extended = new{U}(p) - PointVector(p::AbstractVector) = new{QQFieldElem}(p) + p::MatElem{U} + + PointVector{U}(p::MatElem{U}) where U<:scalar_types_extended = new{U}(p) end Base.IndexStyle(::Type{<:PointVector}) = IndexLinear() -Base.getindex(po::PointVector{T}, i::Base.Integer) where T<:scalar_types_extended = po.p[i] +Base.getindex(po::PointVector{T}, i::Base.Integer) where T<:scalar_types_extended = po.p[1, i] function Base.setindex!(po::PointVector, val, i::Base.Integer) - @boundscheck checkbounds(po.p, i) - po.p[i] = val + @boundscheck 1 <= length(po) <= i + po.p[1, i] = val return val end Base.firstindex(::PointVector) = 1 Base.lastindex(iter::PointVector) = length(iter) -Base.size(po::PointVector) = size(po.p) +Base.size(po::PointVector) = (size(po.p, 2),) -# Forward multiplication with oscar matrices. -Base.:*(A::MatElem, v::PointVector) = A*v.p +coefficient_field(po::PointVector) = base_ring(po.p) -PointVector{U}(p::Union{Field, ZZRing}, v::AbstractVector) where U<:scalar_types_extended = PointVector{U}(p.(v)) +function point_vector(p::Union{scalar_type_or_field, ZZRing}, v::AbstractVector) + parent_field, scalar_type = _determine_parent_and_scalar(p, v) + n = length(v) + mat = matrix(parent_field, 1, n, collect(v)) # collect: workaround for constructor typing + return PointVector{scalar_type}(mat) +end -PointVector{U}(p::Union{Field, ZZRing}, n::Base.Integer) where U<:scalar_types_extended = PointVector{U}(p.(zeros(Int, n))) +function point_vector(p::Union{scalar_type_or_field, ZZRing}, n::Base.Integer) + parent_field, scalar_type = _determine_parent_and_scalar(p) + mat = zero_matrix(parent_field, 1, n) + return PointVector{scalar_type}(mat) +end -PointVector{U}(::UndefInitializer, n::Base.Integer) where U<:scalar_types_extended = PointVector{U}(Vector{U}(undef, n)) +point_vector(x::Union{AbstractVector, Base.Integer}) = point_vector(QQ, x) -PointVector{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = PointVector{U}(undef, length(n)) +# PointVector{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = PointVector{U}(undef, length(n)) function Base.similar(X::PointVector, ::Type{S}, dims::Dims{1}) where S <: scalar_types_extended - return PointVector{S}(undef, dims...) + return point_vector(coefficient_field(X), dims...) end Base.BroadcastStyle(::Type{<:PointVector}) = Broadcast.ArrayStyle{PointVector}() +_parent_or_coefficient_field(po::PointVector) = coefficient_field(po) + function Base.similar(bc::Broadcast.Broadcasted{Broadcast.ArrayStyle{PointVector}}, ::Type{ElType}) where ElType - return PointVector{ElType}(axes(bc)...) + T, f = _promote_scalar_field(_parent_or_coefficient_field.(bc.args)...) + return point_vector(f, axes(bc)...) end ################################################################################ @@ -74,6 +84,8 @@ RayVector{U}(p::Union{Field, ZZRing}, v::AbstractVector) where U<:scalar_types_e RayVector{U}(p::Union{Field, ZZRing}, n::Base.Integer) where U<:scalar_types_extended = RayVector{U}(p.(zeros(Int, n))) +RayVector(p::Union{Field, ZZRing}, x::Union{<:AbstractVector, <:Base.Integer}) = RayVector{elem_type(p)}(p, x) + RayVector{U}(::UndefInitializer, n::Base.Integer) where U<:scalar_types_extended = RayVector{U}(Vector{U}(undef, n)) RayVector{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = RayVector{U}(undef, length(n)) diff --git a/src/PolyhedralGeometry/linear_program.jl b/src/PolyhedralGeometry/linear_program.jl index 7497c6658229..8f825528e149 100644 --- a/src/PolyhedralGeometry/linear_program.jl +++ b/src/PolyhedralGeometry/linear_program.jl @@ -153,7 +153,7 @@ function optimal_vertex(lp::LinearProgram{T}) where {T<:scalar_types} opt_vert = lp.polymake_lp.MINIMAL_VERTEX end if !isnothing(opt_vert) - return PointVector{T}(coefficient_field(lp), view(dehomogenize(opt_vert), :)) + return point_vector(coefficient_field(lp), view(dehomogenize(opt_vert), :))::PointVector{T} else return nothing end diff --git a/src/PolyhedralGeometry/mixed_integer_linear_program.jl b/src/PolyhedralGeometry/mixed_integer_linear_program.jl index 614dfd14c34e..08809e82c9e6 100644 --- a/src/PolyhedralGeometry/mixed_integer_linear_program.jl +++ b/src/PolyhedralGeometry/mixed_integer_linear_program.jl @@ -176,7 +176,7 @@ function optimal_solution(milp::MixedIntegerLinearProgram{T}) where {T<:scalar_t opt_vert = milp.polymake_milp.MINIMAL_SOLUTION end if !isnothing(opt_vert) - return PointVector{T}(dehomogenize(opt_vert)) + return point_vector(coefficient_field(milp), dehomogenize(opt_vert))::PointVector{T} else return nothing end diff --git a/src/exports.jl b/src/exports.jl index 782fee7393b6..2ea8cd877866 100644 --- a/src/exports.jl +++ b/src/exports.jl @@ -1120,6 +1120,7 @@ export pitman_stanley_polytope export platonic_solid export point_coordinates export point_matrix +export point_vector export points export pol_elementary_divisors export polarize diff --git a/test/PolyhedralGeometry/types.jl b/test/PolyhedralGeometry/types.jl index 756d43e26bca..b932d72d273c 100644 --- a/test/PolyhedralGeometry/types.jl +++ b/test/PolyhedralGeometry/types.jl @@ -1,157 +1,157 @@ @testset "types" begin - @testset "IncidenceMatrix" begin + @testset "IncidenceMatrix" begin - im = IncidenceMatrix([[1,2,3],[4,5,6]]) - @test nrows(im) == 2 - @test ncols(im) == 6 - @test row(im, 1) isa Set{Int} - @test row(im, 1) == Set{Int}([1, 2, 3]) - @test column(im, 2) isa Set{Int} - @test column(im, 2) == Set{Int}([1]) + im = IncidenceMatrix([[1,2,3],[4,5,6]]) + @test nrows(im) == 2 + @test ncols(im) == 6 + @test row(im, 1) isa Set{Int} + @test row(im, 1) == Set{Int}([1, 2, 3]) + @test column(im, 2) isa Set{Int} + @test column(im, 2) == Set{Int}([1]) - end - - a = [1, 2, 3] - b = [8, 6, 4] + end - NF, sr2 = quadratic_field(2) + a = [1, 2, 3] + b = [8, 6, 4] - #TODO: RayVector - @testset "$T" for T in (PointVector, RayVector) + NF, sr2 = quadratic_field(2) - @test T(a) isa T{QQFieldElem} + #TODO: RayVector + @testset "$T" for (T, fun) in ((PointVector, point_vector),) - for f in (ZZ, QQ, NF) + @test fun(a) isa T{QQFieldElem} - U = elem_type(f) - @testset "$T{$U}" begin + for f in (ZZ, QQ, NF) - @test T{U} <: AbstractVector - @test T{U} <: AbstractVector{U} + U = elem_type(f) + @testset "$T{$U}" begin - @test T{U}(f.(a)) isa T - @test T{U}(f.(a)) isa T{U} + @test T{U} <: AbstractVector + @test T{U} <: AbstractVector{U} - @test T{U}(f, 7) == f.(zeros(Int, 7)) + @test fun(f, a) isa T + @test fun(f, a) isa T{U} - A = T{U}(f.(a)) + @test fun(f, 7) == f.(zeros(Int, 7)) - @test A[1] isa U - @test A[2] == 2 - @test A[begin] == 1 - @test A[end] == 3 - A[3] = f(7) - @test A[end] == 7 - A[3] = f(3) + A = fun(f, a) - @test size(A) == (3,) + @test A[1] isa U + @test A[2] == 2 + @test A[begin] == 1 + @test A[end] == 3 + A[3] = f(7) + @test A[end] == 7 + A[3] = f(3) - @test_throws BoundsError A[0] - @test_throws BoundsError A[4] + @test size(A) == (3,) - for g in (ZZ, QQ, NF) - V = elem_type(g) + @test_throws BoundsError A[0] + @test_throws BoundsError A[4] - B = T{V}(g.(b)) + for g in (ZZ, QQ, NF) + V = elem_type(g) - for op in [+, -] - @test op(A, B) isa T + B = fun(g, b) - @test op(A, B) == op(a, b) - end + for op in [+, -] + @test op(A, B) isa T - @test *(g(3), A) isa T + @test op(A, B) == op(a, b) + end - @test *(g(3), A) == 3 * a + @test *(g(3), A) isa T - @test [A; B] == [a; b] + @test *(g(3), A) == 3 * a - end + # @test [A; B] == [a; b] - end end + end end + end - @testset "$T" for (T, f) in ((AffineHalfspace, affine_halfspace), (AffineHyperplane, affine_hyperplane)) + + @testset "$T" for (T, f) in ((AffineHalfspace, affine_halfspace), (AffineHyperplane, affine_hyperplane)) - for p in [QQ, NF] - U = elem_type(p) - @test f(p, a, 0) isa T{U} - @test f(p, permutedims(a), 0) isa T{U} + for p in [QQ, NF] + U = elem_type(p) + @test f(p, a, 0) isa T{U} + @test f(p, permutedims(a), 0) isa T{U} - @test f(p, a, 0) == f(p, permutedims(a), 0) == f(p, a) == f(p, permutedims(a)) + @test f(p, a, 0) == f(p, permutedims(a), 0) == f(p, a) == f(p, permutedims(a)) - A = f(p, a, 0) - B = f(p, b, 2) + A = f(p, a, 0) + B = f(p, b, 2) - @test A != B + @test A != B - @test normal_vector(A) isa Vector{U} - @test normal_vector(A) == a - @test negbias(A) isa U - @test negbias(A) == 0 + @test normal_vector(A) isa Vector{U} + @test normal_vector(A) == a + @test negbias(A) isa U + @test negbias(A) == 0 - @test normal_vector(B) == b - @test negbias(B) == 2 - end + @test normal_vector(B) == b + @test negbias(B) == 2 + end - @test f(a, 0) isa T{QQFieldElem} + @test f(a, 0) isa T{QQFieldElem} - end + end - @testset "$T" for (T, f) in ((LinearHalfspace, linear_halfspace), (LinearHyperplane, linear_hyperplane)) + @testset "$T" for (T, f) in ((LinearHalfspace, linear_halfspace), (LinearHyperplane, linear_hyperplane)) - for p in [QQ, NF] - U = elem_type(p) - @test f(p, a) isa T{U} - @test f(p, permutedims(a)) isa T{U} + for p in [QQ, NF] + U = elem_type(p) + @test f(p, a) isa T{U} + @test f(p, permutedims(a)) isa T{U} - @test f(p, a) == f(p, permutedims(a)) + @test f(p, a) == f(p, permutedims(a)) - A = f(p, a) - B = f(p, b) + A = f(p, a) + B = f(p, b) - @test A != B + @test A != B - @test normal_vector(A) isa Vector{U} - @test normal_vector(A) == a - @test negbias(A) isa U - @test negbias(A) == 0 + @test normal_vector(A) isa Vector{U} + @test normal_vector(A) == a + @test negbias(A) isa U + @test negbias(A) == 0 - @test normal_vector(B) == b - @test negbias(B) == 0 - end + @test normal_vector(B) == b + @test negbias(B) == 0 + end - @test f(a) isa T{QQFieldElem} + @test f(a) isa T{QQFieldElem} - end + end - for p in [QQ, NF] - U = elem_type(p) - let A = linear_halfspace(p, a) - @test invert(A) isa LinearHalfspace{U} - Ai = invert(A) - @test normal_vector(Ai) == -a - end - let A = affine_halfspace(p, a, 8) - @test invert(A) isa AffineHalfspace{U} - Ai = invert(A) - @test normal_vector(Ai) == -a - @test negbias(Ai) == -8 - end + for p in [QQ, NF] + U = elem_type(p) + let A = linear_halfspace(p, a) + @test invert(A) isa LinearHalfspace{U} + Ai = invert(A) + @test normal_vector(Ai) == -a + end + let A = affine_halfspace(p, a, 8) + @test invert(A) isa AffineHalfspace{U} + Ai = invert(A) + @test normal_vector(Ai) == -a + @test negbias(Ai) == -8 end + end - @test halfspace(a) isa LinearHalfspace{QQFieldElem} - @test hyperplane(a) isa LinearHyperplane{QQFieldElem} - @test halfspace(a, 0) isa AffineHalfspace{QQFieldElem} - @test hyperplane(a, 0) isa AffineHyperplane{QQFieldElem} + @test halfspace(a) isa LinearHalfspace{QQFieldElem} + @test hyperplane(a) isa LinearHyperplane{QQFieldElem} + @test halfspace(a, 0) isa AffineHalfspace{QQFieldElem} + @test hyperplane(a, 0) isa AffineHyperplane{QQFieldElem} end From 0ccae44333bf4ecdc42b5f4e507b15a27b64f7fa Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Mon, 17 Jul 2023 19:33:37 +0200 Subject: [PATCH 03/21] re-factor for RayVector and PointVector --- src/PolyhedralGeometry/Cone/properties.jl | 4 +- .../PolyhedralComplex/properties.jl | 14 +- .../PolyhedralFan/properties.jl | 4 +- .../Polyhedron/properties.jl | 4 +- src/PolyhedralGeometry/iterators.jl | 127 ++++++------------ src/exports.jl | 1 + test/PolyhedralGeometry/types.jl | 3 +- 7 files changed, 58 insertions(+), 99 deletions(-) diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index 06230cf55750..9492a22b64a3 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -7,7 +7,7 @@ rays(as::Type{RayVector{T}}, C::Cone) where T<:scalar_types = lineality_dim(C) == 0 ? _rays(as, C) : _empty_subobjectiterator(as, C) _rays(as::Type{RayVector{T}}, C::Cone) where T<:scalar_types = SubObjectIterator{as}(C, _ray_cone, _nrays(C)) -_ray_cone(::Type{RayVector{T}}, C::Cone, i::Base.Integer) where T<:scalar_types = RayVector{T}(coefficient_field(C).(view(pm_object(C).RAYS, i, :))) +_ray_cone(U::Type{RayVector{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(C), view(pm_object(C).RAYS, i, :))::U _vector_matrix(::Val{_ray_cone}, C::Cone; homogenized=false) = homogenized ? homogenize(pm_object(C).RAYS, 0) : pm_object(C).RAYS @@ -477,7 +477,7 @@ julia> lineality_space(UH) """ lineality_space(C::Cone{T}) where T<:scalar_types = SubObjectIterator{RayVector{T}}(C, _lineality_cone, lineality_dim(C)) -_lineality_cone(::Type{RayVector{T}}, C::Cone, i::Base.Integer) where T<:scalar_types = RayVector{T}(coefficient_field(C).(view(pm_object(C).LINEALITY_SPACE, i, :))) +_lineality_cone(U::Type{RayVector{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(C), view(pm_object(C).LINEALITY_SPACE, i, :))::U _generator_matrix(::Val{_lineality_cone}, C::Cone; homogenized=false) = homogenized ? homogenize(pm_object(C).LINEALITY_SPACE, 0) : pm_object(C).LINEALITY_SPACE diff --git a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl index 79835daa803c..ad32f051816c 100644 --- a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl @@ -268,12 +268,10 @@ _rays(PC::PolyhedralComplex{T}) where {T<:scalar_types} = _rays(RayVector{T}, PC _ray_indices_polyhedral_complex(PC::Polymake.BigObject) = collect(Polymake.to_one_based_indexing(PC.FAR_VERTICES)) -_ray_polyhedral_complex( - ::Type{RayVector{T}}, PC::PolyhedralComplex, i::Base.Integer -) where {T<:scalar_types} = RayVector{T}( +_ray_polyhedral_complex(U::Type{RayVector{T}}, PC::PolyhedralComplex{T}, i::Base.Integer) where {T<:scalar_types} = ray_vector( coefficient_field(PC), @view pm_object(PC).VERTICES[_ray_indices_polyhedral_complex(pm_object(PC))[i], 2:end] -) +)::U _matrix_for_polymake(::Val{_ray_polyhedral_complex}) = _vector_matrix @@ -493,10 +491,10 @@ end lineality_space(PC::PolyhedralComplex{T}) where {T<:scalar_types} = SubObjectIterator{RayVector{T}}(PC, _lineality_complex, lineality_dim(PC)) -_lineality_complex( - ::Type{RayVector{T}}, PC::PolyhedralComplex, i::Base.Integer -) where {T<:scalar_types} = - RayVector{T}(coefficient_field(PC), @view pm_object(PC).LINEALITY_SPACE[i, 2:end]) +_lineality_complex(U::Type{RayVector{T}}, PC::PolyhedralComplex{T}, i::Base.Integer) where {T<:scalar_types} = ray_vector( + coefficient_field(PC), + @view pm_object(PC).LINEALITY_SPACE[i, 2:end] +)::U _generator_matrix(::Val{_lineality_complex}, PC::PolyhedralComplex; homogenized=false) = if homogenized diff --git a/src/PolyhedralGeometry/PolyhedralFan/properties.jl b/src/PolyhedralGeometry/PolyhedralFan/properties.jl index 05bf4beb8e19..6db0afa94141 100644 --- a/src/PolyhedralGeometry/PolyhedralFan/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralFan/properties.jl @@ -46,7 +46,7 @@ julia> matrix(QQ, rays(NF)) rays(PF::_FanLikeType) = lineality_dim(PF) == 0 ? _rays(PF) : _empty_subobjectiterator(RayVector{_get_scalar_type(PF)}, PF) _rays(PF::_FanLikeType) = SubObjectIterator{RayVector{_get_scalar_type(PF)}}(PF, _ray_fan, _nrays(PF)) -_ray_fan(::Type{RayVector{T}}, PF::_FanLikeType, i::Base.Integer) where T<:scalar_types = RayVector{T}(coefficient_field(PF), view(pm_object(PF).RAYS, i, :)) +_ray_fan(U::Type{RayVector{T}}, PF::_FanLikeType{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(PF), view(pm_object(PF).RAYS, i, :))::U _vector_matrix(::Val{_ray_fan}, PF::_FanLikeType; homogenized=false) = homogenized ? homogenize(pm_object(PF).RAYS, 0) : pm_object(PF).RAYS @@ -385,7 +385,7 @@ julia> lineality_space(PF) """ lineality_space(PF::_FanLikeType) = SubObjectIterator{RayVector{_get_scalar_type(PF)}}(PF, _lineality_fan, lineality_dim(PF)) -_lineality_fan(::Type{RayVector{T}}, PF::_FanLikeType, i::Base.Integer) where T<:scalar_types = RayVector{T}(coefficient_field(PF), view(pm_object(PF).LINEALITY_SPACE, i, :)) +_lineality_fan(U::Type{RayVector{T}}, PF::_FanLikeType{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(PF), view(pm_object(PF).LINEALITY_SPACE, i, :))::U _generator_matrix(::Val{_lineality_fan}, PF::_FanLikeType; homogenized=false) = homogenized ? homogenize(pm_object(PF).LINEALITY_SPACE, 0) : pm_object(PF).LINEALITY_SPACE diff --git a/src/PolyhedralGeometry/Polyhedron/properties.jl b/src/PolyhedralGeometry/Polyhedron/properties.jl index 396ac8185f76..71c949d07b64 100644 --- a/src/PolyhedralGeometry/Polyhedron/properties.jl +++ b/src/PolyhedralGeometry/Polyhedron/properties.jl @@ -314,7 +314,7 @@ julia> rays(RayVector, PO) rays(as::Type{RayVector{T}}, P::Polyhedron) where T<:scalar_types = lineality_dim(P) == 0 ? _rays(as, P) : _empty_subobjectiterator(as, P) _rays(as::Type{RayVector{T}}, P::Polyhedron) where T<:scalar_types = SubObjectIterator{as}(P, _ray_polyhedron, length(_ray_indices(pm_object(P)))) -_ray_polyhedron(::Type{RayVector{T}}, P::Polyhedron, i::Base.Integer) where T<:scalar_types = RayVector{T}(coefficient_field(P).(@view pm_object(P).VERTICES[_ray_indices(pm_object(P))[i], 2:end])) +_ray_polyhedron(U::Type{RayVector{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(P), @view pm_object(P).VERTICES[_ray_indices(pm_object(P))[i], 2:end])::U _facet_indices(::Val{_ray_polyhedron}, P::Polyhedron) = pm_object(P).FACETS_THRU_RAYS[_ray_indices(pm_object(P)), _facet_indices(pm_object(P))] @@ -753,7 +753,7 @@ julia> lineality_space(UH) """ lineality_space(P::Polyhedron{T}) where T<:scalar_types = SubObjectIterator{RayVector{T}}(P, _lineality_polyhedron, lineality_dim(P)) -_lineality_polyhedron(::Type{RayVector{T}}, P::Polyhedron, i::Base.Integer) where T<:scalar_types = RayVector{T}(coefficient_field(P).(@view pm_object(P).LINEALITY_SPACE[i, 2:end])) +_lineality_polyhedron(U::Type{RayVector{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(P), @view pm_object(P).LINEALITY_SPACE[i, 2:end])::U _generator_matrix(::Val{_lineality_polyhedron}, P::Polyhedron; homogenized=false) = homogenized ? pm_object(P).LINEALITY_SPACE : @view pm_object(P).LINEALITY_SPACE[:, 2:end] diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 231764d4a5fc..e672b42b342c 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -2,108 +2,69 @@ ######## Vector types ################################################################################ -struct PointVector{U} <: AbstractVector{U} - p::MatElem{U} +for (T, _t) in ((:PointVector, :point_vector), (:RayVector, :ray_vector)) + + @eval begin - PointVector{U}(p::MatElem{U}) where U<:scalar_types_extended = new{U}(p) -end - -Base.IndexStyle(::Type{<:PointVector}) = IndexLinear() - -Base.getindex(po::PointVector{T}, i::Base.Integer) where T<:scalar_types_extended = po.p[1, i] - -function Base.setindex!(po::PointVector, val, i::Base.Integer) - @boundscheck 1 <= length(po) <= i - po.p[1, i] = val - return val -end - -Base.firstindex(::PointVector) = 1 -Base.lastindex(iter::PointVector) = length(iter) -Base.size(po::PointVector) = (size(po.p, 2),) - -coefficient_field(po::PointVector) = base_ring(po.p) - -function point_vector(p::Union{scalar_type_or_field, ZZRing}, v::AbstractVector) - parent_field, scalar_type = _determine_parent_and_scalar(p, v) - n = length(v) - mat = matrix(parent_field, 1, n, collect(v)) # collect: workaround for constructor typing - return PointVector{scalar_type}(mat) -end - -function point_vector(p::Union{scalar_type_or_field, ZZRing}, n::Base.Integer) - parent_field, scalar_type = _determine_parent_and_scalar(p) - mat = zero_matrix(parent_field, 1, n) - return PointVector{scalar_type}(mat) -end - -point_vector(x::Union{AbstractVector, Base.Integer}) = point_vector(QQ, x) - -# PointVector{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = PointVector{U}(undef, length(n)) - -function Base.similar(X::PointVector, ::Type{S}, dims::Dims{1}) where S <: scalar_types_extended - return point_vector(coefficient_field(X), dims...) -end - -Base.BroadcastStyle(::Type{<:PointVector}) = Broadcast.ArrayStyle{PointVector}() - -_parent_or_coefficient_field(po::PointVector) = coefficient_field(po) + struct $T{U} <: AbstractVector{U} + p::MatElem{U} -function Base.similar(bc::Broadcast.Broadcasted{Broadcast.ArrayStyle{PointVector}}, ::Type{ElType}) where ElType - T, f = _promote_scalar_field(_parent_or_coefficient_field.(bc.args)...) - return point_vector(f, axes(bc)...) -end + $T{U}(p::MatElem{U}) where U<:scalar_types_extended = new{U}(p) + end -################################################################################ + Base.IndexStyle(::Type{<:$T}) = IndexLinear() -struct RayVector{U} <: AbstractVector{U} - p::Vector{U} - - RayVector{U}(p::AbstractVector) where U<:scalar_types_extended = new{U}(p) - RayVector(p::AbstractVector) = new{QQFieldElem}(p) -end + Base.getindex(po::$T, i::Base.Integer) = po.p[1, i] -Base.IndexStyle(::Type{<:RayVector}) = IndexLinear() + function Base.setindex!(po::$T, val, i::Base.Integer) + @boundscheck 1 <= length(po) <= i + po.p[1, i] = val + return val + end -Base.getindex(po::RayVector{T}, i::Base.Integer) where T<:scalar_types_extended = po.p[i] + Base.firstindex(::$T) = 1 + Base.lastindex(iter::$T) = length(iter) + Base.size(po::$T) = (size(po.p, 2),) -function Base.setindex!(po::RayVector, val, i::Base.Integer) - @boundscheck checkbounds(po.p, i) - po.p[i] = val - return val -end + coefficient_field(po::$T) = base_ring(po.p) -Base.firstindex(::RayVector) = 1 -Base.lastindex(iter::RayVector) = length(iter) -Base.size(po::RayVector) = size(po.p) + function $_t(p::Union{scalar_type_or_field, ZZRing}, v::AbstractVector) + parent_field, scalar_type = _determine_parent_and_scalar(p, v) + n = length(v) + mat = matrix(parent_field, 1, n, collect(v)) # collect: workaround for constructor typing + return $T{scalar_type}(mat) + end -# Forward multiplication with oscar matrices. -Base.:*(A::MatElem, v::RayVector) = A*v.p + function $_t(p::Union{scalar_type_or_field, ZZRing}, n::Base.Integer) + parent_field, scalar_type = _determine_parent_and_scalar(p) + mat = zero_matrix(parent_field, 1, n) + return $T{scalar_type}(mat) + end -RayVector{U}(p::Union{Field, ZZRing}, v::AbstractVector) where U<:scalar_types_extended = RayVector{U}(p.(v)) + $_t(x::Union{AbstractVector, Base.Integer}) = $_t(QQ, x) -RayVector{U}(p::Union{Field, ZZRing}, n::Base.Integer) where U<:scalar_types_extended = RayVector{U}(p.(zeros(Int, n))) + # $T{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = $T{U}(undef, length(n)) -RayVector(p::Union{Field, ZZRing}, x::Union{<:AbstractVector, <:Base.Integer}) = RayVector{elem_type(p)}(p, x) + function Base.similar(X::$T, ::Type{S}, dims::Dims{1}) where S <: scalar_types_extended + return $_t(coefficient_field(X), dims...) + end -RayVector{U}(::UndefInitializer, n::Base.Integer) where U<:scalar_types_extended = RayVector{U}(Vector{U}(undef, n)) + Base.BroadcastStyle(::Type{<:$T}) = Broadcast.ArrayStyle{$T}() -RayVector{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = RayVector{U}(undef, length(n)) + _parent_or_coefficient_field(po::$T) = coefficient_field(po) -function Base.similar(X::RayVector, ::Type{S}, dims::Dims{1}) where S <: scalar_types_extended - return RayVector{S}(undef, dims...) -end + function Base.similar(bc::Broadcast.Broadcasted{Broadcast.ArrayStyle{$T}}, ::Type{ElType}) where ElType + U, f = _promote_scalar_field(_parent_or_coefficient_field.(bc.args)...) + return $_t(f, axes(bc)...) + end -Base.BroadcastStyle(::Type{<:RayVector}) = Broadcast.ArrayStyle{RayVector}() + Base.:*(k::scalar_types_extended, po::$T) = k .* po -function Base.similar(bc::Broadcast.Broadcasted{Broadcast.ArrayStyle{RayVector}}, ::Type{ElType}) where ElType - return RayVector{ElType}(axes(bc)...) + Base.:*(A::MatElem, v::$T) = A * v.p + + end end -################################################################################ - -Base.:*(k::scalar_types_extended, po::Union{PointVector, RayVector}) = k .* po - ################################################################################ ######## Halfspaces and Hyperplanes ################################################################################ diff --git a/src/exports.jl b/src/exports.jl index 2ea8cd877866..d05f66f3d168 100644 --- a/src/exports.jl +++ b/src/exports.jl @@ -1196,6 +1196,7 @@ export rational_equivalence_class export rational_solutions export rational_to_continued_fraction_hirzebruch_jung export ray_indices +export ray_vector export rays export rays_modulo_lineality export real_projective_plane diff --git a/test/PolyhedralGeometry/types.jl b/test/PolyhedralGeometry/types.jl index b932d72d273c..726b21f89e1c 100644 --- a/test/PolyhedralGeometry/types.jl +++ b/test/PolyhedralGeometry/types.jl @@ -17,8 +17,7 @@ NF, sr2 = quadratic_field(2) - #TODO: RayVector - @testset "$T" for (T, fun) in ((PointVector, point_vector),) + @testset "$T" for (T, fun) in ((PointVector, point_vector), (RayVector, ray_vector)) @test fun(a) isa T{QQFieldElem} From b600337bafba085a1a35c1a19f53040f82647c3a Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Mon, 17 Jul 2023 21:34:37 +0200 Subject: [PATCH 04/21] remove (not embeddeded) nf_elem from testfiles --- test/PolyhedralGeometry/cone.jl | 2 +- test/PolyhedralGeometry/polyhedral_complex.jl | 4 ++-- test/PolyhedralGeometry/polyhedral_fan.jl | 4 +--- 3 files changed, 4 insertions(+), 6 deletions(-) diff --git a/test/PolyhedralGeometry/cone.jl b/test/PolyhedralGeometry/cone.jl index 24cc9af84c45..2098289e1a6c 100644 --- a/test/PolyhedralGeometry/cone.jl +++ b/test/PolyhedralGeometry/cone.jl @@ -4,7 +4,7 @@ NF, sr2 = quadratic_field(2) Qx, x = QQ["x"] K, (a1, a2) = embedded_number_field([x^2 - 2, x^3 - 5], [(0, 2), (0, 2)]) -for f in (QQ, NF, K) +for f in (QQ, K) T = elem_type(f) diff --git a/test/PolyhedralGeometry/polyhedral_complex.jl b/test/PolyhedralGeometry/polyhedral_complex.jl index e9544c058f92..734852859a14 100644 --- a/test/PolyhedralGeometry/polyhedral_complex.jl +++ b/test/PolyhedralGeometry/polyhedral_complex.jl @@ -1,8 +1,8 @@ NF, sr2 = quadratic_field(2) Qx, x = QQ["x"] -K, (a1, a2) = embedded_number_field([x^2 - 2, x^3 - 5], [(0, 2), (0, 2)]) +K, a = Hecke.embedded_field(NF, real_embeddings(NF)[2]) -for f in (QQ, NF, K) +for f in (QQ, K) T = elem_type(f) @testset "PolyhedralComplex{$T}" begin diff --git a/test/PolyhedralGeometry/polyhedral_fan.jl b/test/PolyhedralGeometry/polyhedral_fan.jl index 4b5e9065c58c..f08d3b5e0b17 100644 --- a/test/PolyhedralGeometry/polyhedral_fan.jl +++ b/test/PolyhedralGeometry/polyhedral_fan.jl @@ -1,9 +1,7 @@ -NF, sr2 = quadratic_field(2) -ENF, sre2 = Hecke.embedded_field(NF, real_embeddings(NF)[2]) Qx, x = QQ["x"] K, (a1, a2) = embedded_number_field([x^2 - 2, x^3 - 5], [(0, 2), (0, 2)]) -for f in (QQ, ENF, K) +for f in (QQ, K) T = elem_type(f) @testset "PolyhedralFan{$T}" begin From 4021a5713896a695d14ba635b1593029703ffef8 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Wed, 19 Jul 2023 23:33:12 +0200 Subject: [PATCH 05/21] improve usage of where T<:scalar_types; adjusted docs --- .../PolyhedralGeometry/Polyhedra/auxiliary.md | 2 +- .../Polyhedra/constructions.md | 6 +- docs/src/PolyhedralGeometry/cones.md | 2 +- .../subdivisions_of_points.md | 6 +- src/PolyhedralGeometry/Cone/properties.jl | 32 ++++--- .../PolyhedralComplex/properties.jl | 96 ++++++++----------- .../PolyhedralFan/properties.jl | 18 ++-- .../Polyhedron/constructors.jl | 2 +- .../Polyhedron/properties.jl | 90 +++++++++-------- .../SubdivisionOfPoints/properties.jl | 8 +- src/PolyhedralGeometry/helpers.jl | 2 + 11 files changed, 127 insertions(+), 137 deletions(-) diff --git a/docs/src/PolyhedralGeometry/Polyhedra/auxiliary.md b/docs/src/PolyhedralGeometry/Polyhedra/auxiliary.md index adebd185c2f3..4567224be2f0 100644 --- a/docs/src/PolyhedralGeometry/Polyhedra/auxiliary.md +++ b/docs/src/PolyhedralGeometry/Polyhedra/auxiliary.md @@ -64,7 +64,7 @@ is_smooth(P::Polyhedron{QQFieldElem}) is_very_ample(P::Polyhedron{QQFieldElem}) lattice_points(P::Polyhedron{QQFieldElem}) lattice_volume(P::Polyhedron{QQFieldElem}) -normalized_volume(P::Polyhedron{T}) where T<:scalar_types +normalized_volume(P::Polyhedron) polarize(P::Polyhedron{T}) where T<:scalar_types project_full(P::Polyhedron{T}) where T<:scalar_types print_constraints(A::AnyVecOrMat, b::AbstractVector; trivial::Bool = false) diff --git a/docs/src/PolyhedralGeometry/Polyhedra/constructions.md b/docs/src/PolyhedralGeometry/Polyhedra/constructions.md index 777d26eef733..e6992c1bd70f 100644 --- a/docs/src/PolyhedralGeometry/Polyhedra/constructions.md +++ b/docs/src/PolyhedralGeometry/Polyhedra/constructions.md @@ -17,8 +17,8 @@ from other objects in OSCAR. ### Intersecting halfspaces: $H$-representation ```@docs -polyhedron(::Type{T}, A::AnyVecOrMat, b::AbstractVector) where T<:scalar_types -polyhedron(::Type{T}, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) where T<:scalar_types +polyhedron(::Oscar.scalar_type_or_field, A::AnyVecOrMat, b::AbstractVector) +polyhedron(::Oscar.scalar_type_or_field, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) ``` The complete $H$-representation can be retrieved using [`facets`](@ref facets) @@ -54,7 +54,7 @@ true ### Computing convex hulls: $V$-representation ```@docs -convex_hull(::Type{T}, ::AnyVecOrMat; non_redundant::Bool=false) where T<:scalar_types +convex_hull(::Oscar.scalar_type_or_field, ::AnyVecOrMat; non_redundant::Bool=false) ``` This is a standard triangle, defined via a (redundant) $V$-representation and diff --git a/docs/src/PolyhedralGeometry/cones.md b/docs/src/PolyhedralGeometry/cones.md index 95ef93ad2c3f..404c7ea5fe2c 100644 --- a/docs/src/PolyhedralGeometry/cones.md +++ b/docs/src/PolyhedralGeometry/cones.md @@ -27,7 +27,7 @@ reason for keeping cones as a distinct type. ## Construction ```@docs -positive_hull(::Type{T}, ::Union{Oscar.MatElem, AbstractMatrix}) where T<:scalar_types +positive_hull(::Oscar.scalar_type_or_field, ::Union{Oscar.MatElem, AbstractMatrix}) cone_from_inequalities cone_from_equations secondary_cone(SOP::SubdivisionOfPoints{T}) where T<:scalar_types diff --git a/docs/src/PolyhedralGeometry/subdivisions_of_points.md b/docs/src/PolyhedralGeometry/subdivisions_of_points.md index 3b2f3b9a0825..4e4e9039af95 100644 --- a/docs/src/PolyhedralGeometry/subdivisions_of_points.md +++ b/docs/src/PolyhedralGeometry/subdivisions_of_points.md @@ -34,8 +34,8 @@ vector, but if it does, it is called `regular`. ```@docs -subdivision_of_points(::Type{T}, Points::Union{Oscar.MatElem,AbstractMatrix}, Incidence::IncidenceMatrix) where T<:scalar_types -subdivision_of_points(::Type{T}, Points::Union{Oscar.MatElem,AbstractMatrix}, Weights::AbstractVector) where T<:scalar_types +subdivision_of_points(::Oscar.scalar_type_or_field, Points::Union{Oscar.MatElem,AbstractMatrix}, Incidence::IncidenceMatrix) +subdivision_of_points(::Oscar.scalar_type_or_field, Points::Union{Oscar.MatElem,AbstractMatrix}, Weights::AbstractVector) ``` From a subdivision of points one can construct the @@ -50,6 +50,6 @@ is_regular(SOP::SubdivisionOfPoints) maximal_cells min_weights n_maximal_cells(SOP::SubdivisionOfPoints) -points(SOP::SubdivisionOfPoints) +points(SOP::SubdivisionOfPoints{T}) where T<:scalar_types secondary_cone ``` diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index 9492a22b64a3..4fa91139b766 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -4,8 +4,8 @@ ############################################################################### ############################################################################### -rays(as::Type{RayVector{T}}, C::Cone) where T<:scalar_types = lineality_dim(C) == 0 ? _rays(as, C) : _empty_subobjectiterator(as, C) -_rays(as::Type{RayVector{T}}, C::Cone) where T<:scalar_types = SubObjectIterator{as}(C, _ray_cone, _nrays(C)) +rays(as::Type{RayVector{T}}, C::Cone{T}) where T<:scalar_types = lineality_dim(C) == 0 ? _rays(as, C) : _empty_subobjectiterator(as, C) +_rays(as::Type{RayVector{T}}, C::Cone{T}) where T<:scalar_types = SubObjectIterator{as}(C, _ray_cone, _nrays(C)) _ray_cone(U::Type{RayVector{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(C), view(pm_object(C).RAYS, i, :))::U @@ -13,8 +13,8 @@ _vector_matrix(::Val{_ray_cone}, C::Cone; homogenized=false) = homogenized ? hom _matrix_for_polymake(::Val{_ray_cone}) = _vector_matrix -rays(::Type{RayVector}, C::Cone{T}) where T<:scalar_types = rays(RayVector{T}, C) -_rays(::Type{RayVector}, C::Cone{T}) where T<:scalar_types = _rays(RayVector{T}, C) +rays(::Type{<:RayVector}, C::Cone{T}) where T<:scalar_types = rays(RayVector{T}, C) +_rays(::Type{<:RayVector}, C::Cone{T}) where T<:scalar_types = _rays(RayVector{T}, C) @doc raw""" rays(C::Cone) @@ -111,13 +111,13 @@ julia> RML.lineality_basis ``` """ rays_modulo_lineality(C::Cone{T}) where T<:scalar_types = rays_modulo_lineality(NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}, C) -function rays_modulo_lineality(as::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, C::Cone) where T<:scalar_types +function rays_modulo_lineality(::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, C::Cone{T}) where T<:scalar_types return ( rays_modulo_lineality = _rays(C), lineality_basis = lineality_space(C) ) end -rays_modulo_lineality(as::Type{RayVector}, C::Cone) = _rays(C) +rays_modulo_lineality(::Type{<:RayVector}, C::Cone) = _rays(C) @doc raw""" @@ -147,8 +147,8 @@ function faces(C::Cone{T}, face_dim::Int) where T<:scalar_types return SubObjectIterator{Cone{T}}(C, _face_cone, size(Polymake.polytope.faces_of_dim(pm_object(C), n), 1), (f_dim = n,)) end -function _face_cone(::Type{Cone{T}}, C::Cone, i::Base.Integer; f_dim::Int = 0) where T<:scalar_types - return Cone{T}(Polymake.polytope.Cone{_scalar_type_to_polymake(T)}(RAYS = pm_object(C).RAYS[collect(Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(C), f_dim)[i])), :], LINEALITY_SPACE = pm_object(C).LINEALITY_SPACE)) +function _face_cone(U::Type{Cone{T}}, C::Cone{T}, i::Base.Integer; f_dim::Int = 0) where T<:scalar_types + return cone(coefficient_field(C), pm_object(C).RAYS[collect(Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(C), f_dim)[i])), :], pm_object(C).LINEALITY_SPACE)::U end function _ray_indices(::Val{_face_cone}, C::Cone; f_dim::Int = 0) @@ -156,8 +156,8 @@ function _ray_indices(::Val{_face_cone}, C::Cone; f_dim::Int = 0) return IncidenceMatrix([collect(f[i]) for i in 1:length(f)]) end -function _face_cone_facet(::Type{Cone{T}}, C::Cone, i::Base.Integer) where T<:scalar_types - return Cone{T}(Polymake.polytope.Cone{_scalar_type_to_polymake(T)}(RAYS = pm_object(C).RAYS[collect(pm_object(C).RAYS_IN_FACETS[i, :]), :], LINEALITY_SPACE = pm_object(C).LINEALITY_SPACE), coefficient_field(C)) +function _face_cone_facet(U::Type{Cone{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types + return cone(coefficient_field(C), pm_object(C).RAYS[collect(pm_object(C).RAYS_IN_FACETS[i, :]), :], pm_object(C).LINEALITY_SPACE)::U end _ray_indices(::Val{_face_cone_facet}, C::Cone) = pm_object(C).RAYS_IN_FACETS @@ -441,9 +441,9 @@ julia> f = facets(Halfspace, c) """ facets(as::Type{<:Union{LinearHalfspace{T}, Cone{T}}}, C::Cone) where T<:scalar_types = SubObjectIterator{as}(C, _facet_cone, nfacets(C)) -_facet_cone(::Type{LinearHalfspace{T}}, C::Cone, i::Base.Integer) where T<:scalar_types = linear_halfspace(coefficient_field(C), -pm_object(C).FACETS[[i], :]) +_facet_cone(U::Type{LinearHalfspace{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = linear_halfspace(coefficient_field(C), -pm_object(C).FACETS[[i], :])::U -_facet_cone(::Type{Cone{T}}, C::Cone, i::Base.Integer) where T<:scalar_types = Cone{T}(Polymake.polytope.facet(pm_object(C), i-1), coefficient_field(C)) +_facet_cone(::Type{Cone{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = Cone{T}(Polymake.polytope.facet(pm_object(C), i-1), coefficient_field(C)) _linear_inequality_matrix(::Val{_facet_cone}, C::Cone) = -pm_object(C).FACETS @@ -455,7 +455,9 @@ _incidencematrix(::Val{_facet_cone}) = _ray_indices facets(C::Cone{T}) where T<:scalar_types = facets(LinearHalfspace{T}, C) -facets(::Type{Halfspace}, C::Cone{T}) where T<:scalar_types = facets(LinearHalfspace{T}, C) +facets(::Type{<:Halfspace}, C::Cone{T}) where T<:scalar_types = facets(LinearHalfspace{T}, C) + +facets(::Type{<:AffineHalfspace}, C::Cone{T}) where T<:scalar_types = facets(AffineHalfspace{T}, C) facets(::Type{Cone}, C::Cone{T}) where T<:scalar_types = facets(Cone{T}, C) @@ -501,7 +503,7 @@ x₃ = 0 """ linear_span(C::Cone{T}) where T<:scalar_types = SubObjectIterator{LinearHyperplane{T}}(C, _linear_span, size(pm_object(C).LINEAR_SPAN, 1)) -_linear_span(::Type{LinearHyperplane{T}}, C::Cone, i::Base.Integer) where T<:scalar_types = linear_hyperplane(coefficient_field(C), view(pm_object(C).LINEAR_SPAN, i, :)) +_linear_span(U::Type{LinearHyperplane{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = linear_hyperplane(coefficient_field(C), view(pm_object(C).LINEAR_SPAN, i, :))::U _linear_equation_matrix(::Val{_linear_span}, C::Cone) = pm_object(C).LINEAR_SPAN @@ -530,7 +532,7 @@ function hilbert_basis(C::Cone{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(C, _hilbert_generator, size(pm_object(C).HILBERT_BASIS_GENERATORS[1], 1)) end -_hilbert_generator(T::Type{PointVector{ZZRingElem}}, C::Cone, i::Base.Integer) = point_vector(ZZ, view(pm_object(C).HILBERT_BASIS_GENERATORS[1], i, :))::T +_hilbert_generator(T::Type{PointVector{ZZRingElem}}, C::Cone{QQFieldElem}, i::Base.Integer) = point_vector(ZZ, view(pm_object(C).HILBERT_BASIS_GENERATORS[1], i, :))::T _generator_matrix(::Val{_hilbert_generator}, C::Cone; homogenized=false) = homogenized ? homogenize(pm_object(C).HILBERT_BASIS_GENERATORS[1], 0) : pm_object(C).HILBERT_BASIS_GENERATORS[1] diff --git a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl index ad32f051816c..839cc70ae4cd 100644 --- a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl @@ -60,7 +60,7 @@ julia> matrix(QQ, vertices(PC)) """ vertices(as::Type{PointVector{T}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} = lineality_dim(PC) == 0 ? _vertices(as, PC) : _empty_subobjectiterator(as, PC) -_vertices(as::Type{PointVector{T}}, PC::PolyhedralComplex) where {T<:scalar_types} = +_vertices(as::Type{PointVector{T}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} = SubObjectIterator{as}(PC, _vertex_complex, length(_vertex_indices(pm_object(PC)))) _vertex_complex(U::Type{PointVector{T}}, PC::PolyhedralComplex{T}, i::Base.Integer) where {T<:scalar_types} = point_vector( @@ -98,11 +98,10 @@ function _vertex_or_ray_complex(::Type{Union{PointVector{T}, RayVector{T}}}, PC: end end -_vertices( - as::Type{Union{RayVector{T},PointVector{T}}}, PC::PolyhedralComplex -) where {T<:scalar_types} = SubObjectIterator{as}( - PC, _vertex_or_ray_complex, length(_all_vertex_indices(pm_object(PC))) -) +_vertices(as::Type{Union{RayVector{T}, PointVector{T}}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} = + SubObjectIterator{as}( + PC, _vertex_or_ray_complex, length(_all_vertex_indices(pm_object(PC))) + ) vertices_and_rays(PC::PolyhedralComplex{T}) where {T<:scalar_types} = _vertices(Union{PointVector{T},RayVector{T}}, PC) @@ -116,7 +115,7 @@ _vertices(::Type{PointVector}, PC::PolyhedralComplex{T}) where {T<:scalar_types} vertices(PC::PolyhedralComplex{T}) where {T<:scalar_types} = vertices(PointVector{T}, PC) _vertices(PC::PolyhedralComplex{T}) where {T<:scalar_types} = _vertices(PointVector{T}, PC) -_rays(as::Type{RayVector{T}}, PC::PolyhedralComplex) where {T<:scalar_types} = +_rays(as::Type{RayVector{T}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} = SubObjectIterator{as}( PC, _ray_polyhedral_complex, length(_ray_indices_polyhedral_complex(pm_object(PC))) ) @@ -161,28 +160,18 @@ julia> RML.lineality_basis """ rays_modulo_lineality(PC::PolyhedralComplex{T}) where {T<:scalar_types} = rays_modulo_lineality( - NamedTuple{ - (:rays_modulo_lineality, :lineality_basis), - Tuple{SubObjectIterator{RayVector{T}},SubObjectIterator{RayVector{T}}}, - }, - PC, - ) -function rays_modulo_lineality( - as::Type{ - NamedTuple{ - (:rays_modulo_lineality, :lineality_basis), - Tuple{SubObjectIterator{RayVector{T}},SubObjectIterator{RayVector{T}}}, - }, - }, - PC::PolyhedralComplex, -) where {T<:scalar_types} - return ( - rays_modulo_lineality=_rays(RayVector{T}, PC), lineality_basis=lineality_space(PC) + NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}, PC ) + +function rays_modulo_lineality(::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} + return ( + rays_modulo_lineality = _rays(RayVector{T}, PC), + lineality_basis = lineality_space(PC) + ) end -rays_modulo_lineality( - as::Type{RayVector{T}}, PC::PolyhedralComplex{T} -) where {T<:scalar_types} = _rays(RayVector{T}, PC) + +rays_modulo_lineality(U::Type{RayVector{T}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} = + _rays(U, PC) @doc raw""" minimal_faces(as, PC::PolyhedralComplex) @@ -219,24 +208,19 @@ julia> MFPC.lineality_basis [0, 0, 1] ``` """ -minimal_faces(PC::PolyhedralComplex{T}) where {T<:scalar_types} = minimal_faces( - NamedTuple{ - (:base_points, :lineality_basis), - Tuple{SubObjectIterator{PointVector{T}},SubObjectIterator{RayVector{T}}}, - }, - PC, -) -function minimal_faces( - as::Type{ - NamedTuple{ - (:base_points, :lineality_basis), - Tuple{SubObjectIterator{PointVector{T}},SubObjectIterator{RayVector{T}}}, - }, - }, - PC::PolyhedralComplex{T}, -) where {T<:scalar_types} - return (base_points=_vertices(PointVector{T}, PC), lineality_basis=lineality_space(PC)) +minimal_faces(PC::PolyhedralComplex{T}) where {T<:scalar_types} = + minimal_faces( + NamedTuple{(:base_points, :lineality_basis), Tuple{SubObjectIterator{PointVector{T}}, SubObjectIterator{RayVector{T}}}}, + PC + ) + +function minimal_faces(::Type{NamedTuple{(:base_points, :lineality_basis), Tuple{SubObjectIterator{PointVector{T}}, SubObjectIterator{RayVector{T}}}}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} + return ( + base_points = _vertices(PointVector{T}, PC), + lineality_basis = lineality_space(PC) + ) end + minimal_faces(as::Type{PointVector{T}}, PC::PolyhedralComplex{T}) where {T<:scalar_types} = _vertices(PointVector{T}, PC) @@ -280,10 +264,11 @@ _vector_matrix(::Val{_ray_polyhedral_complex}, PC::PolyhedralComplex; homogenize _ray_indices_polyhedral_complex(pm_object(PC)), (homogenized ? 1 : 2):end ] -_maximal_polyhedron( - ::Type{Polyhedron{T}}, PC::PolyhedralComplex, i::Base.Integer -) where {T<:scalar_types} = - Polyhedron{T}(Polymake.fan.polytope(pm_object(PC), i - 1), coefficient_field(PC)) +_maximal_polyhedron(::Type{Polyhedron{T}}, PC::PolyhedralComplex{T}, i::Base.Integer) where {T<:scalar_types} = + Polyhedron{T}( + Polymake.fan.polytope(pm_object(PC), i - 1), + coefficient_field(PC) + ) _vertex_indices(::Val{_maximal_polyhedron}, PC::PolyhedralComplex) = pm_object(PC).MAXIMAL_POLYTOPES[:, _vertex_indices(pm_object(PC))] @@ -473,18 +458,15 @@ end function _ith_polyhedron( ::Type{Polyhedron{T}}, - PC::PolyhedralComplex, + PC::PolyhedralComplex{T}, i::Base.Integer; - f_dim::Int=-1, - f_ind::Vector{Int64}=Vector{Int64}(), -) where {T<:scalar_types} + f_dim::Int = -1, + f_ind::Vector{Int64} = Vector{Int64}() + ) where {T<:scalar_types} pface = Polymake.row(Polymake.fan.cones_of_dim(pm_object(PC), f_dim), f_ind[i]) return Polyhedron{T}( - Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(; - VERTICES=pm_object(PC).VERTICES[collect(pface), :], - LINEALITY_SPACE=pm_object(PC).LINEALITY_SPACE, - ), - coefficient_field(PC), + Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(PC).VERTICES[collect(pface), :], LINEALITY_SPACE = pm_object(PC).LINEALITY_SPACE), + coefficient_field(PC) ) end diff --git a/src/PolyhedralGeometry/PolyhedralFan/properties.jl b/src/PolyhedralGeometry/PolyhedralFan/properties.jl index 6db0afa94141..f8455291c210 100644 --- a/src/PolyhedralGeometry/PolyhedralFan/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralFan/properties.jl @@ -52,7 +52,7 @@ _vector_matrix(::Val{_ray_fan}, PF::_FanLikeType; homogenized=false) = homogeniz _matrix_for_polymake(::Val{_ray_fan}) = _vector_matrix -_maximal_cone(::Type{Cone{T}}, PF::_FanLikeType, i::Base.Integer) where T<:scalar_types = Cone{T}(Polymake.fan.cone(pm_object(PF), i - 1), coefficient_field(PF)) +_maximal_cone(::Type{Cone{T}}, PF::_FanLikeType{T}, i::Base.Integer) where T<:scalar_types = Cone{T}(Polymake.fan.cone(pm_object(PF), i - 1), coefficient_field(PF)) @doc raw""" @@ -87,14 +87,20 @@ julia> rays(NF) 0-element SubObjectIterator{RayVector{QQFieldElem}} ``` """ -rays_modulo_lineality(PF::_FanLikeType) = rays_modulo_lineality(NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{_get_scalar_type(PF)}}, SubObjectIterator{RayVector{_get_scalar_type(PF)}}}}, PF) -function rays_modulo_lineality(as::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, F::_FanLikeType) where T<:scalar_types +rays_modulo_lineality(F::_FanLikeType{T}) where {T<:scalar_types} = + rays_modulo_lineality( + NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}, + F + ) + +function rays_modulo_lineality(::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, F::_FanLikeType{T}) where T<:scalar_types return ( rays_modulo_lineality = _rays(F), lineality_basis = lineality_space(F) ) end -rays_modulo_lineality(as::Type{RayVector}, F::_FanLikeType) = _rays(F) + +rays_modulo_lineality(::Type{<:RayVector}, F::_FanLikeType) = _rays(F) @@ -158,8 +164,8 @@ function cones(PF::_FanLikeType, cone_dim::Int) return SubObjectIterator{Cone{_get_scalar_type(PF)}}(PF, _cone_of_dim, size(Polymake.fan.cones_of_dim(pm_object(PF), l), 1), (c_dim = l,)) end -function _cone_of_dim(::Type{Cone{T}}, PF::_FanLikeType, i::Base.Integer; c_dim::Int = 0) where T<:scalar_types - return Cone{T}(Polymake.polytope.Cone{_scalar_type_to_polymake(T)}(RAYS = pm_object(PF).RAYS[collect(Polymake.row(Polymake.fan.cones_of_dim(pm_object(PF), c_dim), i)),:], LINEALITY_SPACE = pm_object(PF).LINEALITY_SPACE), coefficient_field(PF)) +function _cone_of_dim(U::Type{Cone{T}}, PF::_FanLikeType{T}, i::Base.Integer; c_dim::Int = 0) where T<:scalar_types + return cone(coefficient_field(PF), pm_object(PF).RAYS[collect(Polymake.row(Polymake.fan.cones_of_dim(pm_object(PF), c_dim), i)),:], pm_object(PF).LINEALITY_SPACE)::U end _ray_indices(::Val{_cone_of_dim}, PF::_FanLikeType; c_dim::Int = 0) = Polymake.fan.cones_of_dim(pm_object(PF), c_dim) diff --git a/src/PolyhedralGeometry/Polyhedron/constructors.jl b/src/PolyhedralGeometry/Polyhedron/constructors.jl index 0f2bf83b2a90..9b175bebb6c4 100644 --- a/src/PolyhedralGeometry/Polyhedron/constructors.jl +++ b/src/PolyhedralGeometry/Polyhedron/constructors.jl @@ -106,7 +106,7 @@ julia> vertices(P) [0, 0] ``` """ -function polyhedron(f::scalar_type_or_field, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing; parent_field::Union{Nothing, Field} = nothing) +function polyhedron(f::scalar_type_or_field, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) parent_field, scalar_type = _determine_parent_and_scalar(f, I, E) if isnothing(I) || _isempty_halfspace(I) EM = affine_matrix_for_polymake(E) diff --git a/src/PolyhedralGeometry/Polyhedron/properties.jl b/src/PolyhedralGeometry/Polyhedron/properties.jl index 71c949d07b64..81fd0f74d712 100644 --- a/src/PolyhedralGeometry/Polyhedron/properties.jl +++ b/src/PolyhedralGeometry/Polyhedron/properties.jl @@ -35,9 +35,9 @@ function faces(P::Polyhedron{T}, face_dim::Int) where T<:scalar_types return SubObjectIterator{Polyhedron{T}}(P, _face_polyhedron, length(rfaces), (f_dim = n, f_ind = rfaces)) end -function _face_polyhedron(::Type{Polyhedron{T}}, P::Polyhedron, i::Base.Integer; f_dim::Int = -1, f_ind::Vector{Int64} = Vector{Int64}()) where T<:scalar_types - pface = Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(P), f_dim))[f_ind[i]] - return Polyhedron{T}(Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(P).VERTICES[collect(pface),:], LINEALITY_SPACE = pm_object(P).LINEALITY_SPACE), coefficient_field(P)) +function _face_polyhedron(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base.Integer; f_dim::Int = -1, f_ind::Vector{Int64} = Vector{Int64}()) where T<:scalar_types + pface = Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(P), f_dim))[f_ind[i]] + return Polyhedron{T}(Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(P).VERTICES[collect(pface), :], LINEALITY_SPACE = pm_object(P).LINEALITY_SPACE), coefficient_field(P)) end function _vertex_indices(::Val{_face_polyhedron}, P::Polyhedron; f_dim = -1, f_ind::Vector{Int64} = Vector{Int64}()) @@ -54,9 +54,9 @@ end _incidencematrix(::Val{_face_polyhedron}) = _vertex_and_ray_indices -function _face_polyhedron_facet(::Type{Polyhedron{T}}, P::Polyhedron, i::Base.Integer) where T<:scalar_types - pface = pm_object(P).VERTICES_IN_FACETS[_facet_index(pm_object(P), i), :] - return Polyhedron{T}(Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(P).VERTICES[collect(pface),:], LINEALITY_SPACE = pm_object(P).LINEALITY_SPACE), coefficient_field(P)) +function _face_polyhedron_facet(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types + pface = pm_object(P).VERTICES_IN_FACETS[_facet_index(pm_object(P), i), :] + return Polyhedron{T}(Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(P).VERTICES[collect(pface), :], LINEALITY_SPACE = pm_object(P).LINEALITY_SPACE), coefficient_field(P)) end _vertex_indices(::Val{_face_polyhedron_facet}, P::Polyhedron) = vcat(pm_object(P).VERTICES_IN_FACETS[1:(_facet_at_infinity(pm_object(P)) - 1), _vertex_indices(pm_object(P))], pm_object(P).VERTICES_IN_FACETS[(_facet_at_infinity(pm_object(P)) + 1):end, _vertex_indices(pm_object(P))]) @@ -127,13 +127,13 @@ julia> minimal_faces(P) ``` """ minimal_faces(P::Polyhedron{T}) where T<:scalar_types = minimal_faces(NamedTuple{(:base_points, :lineality_basis), Tuple{SubObjectIterator{PointVector{T}}, SubObjectIterator{RayVector{T}}}}, P) -function minimal_faces(as::Type{NamedTuple{(:base_points, :lineality_basis), Tuple{SubObjectIterator{PointVector{T}}, SubObjectIterator{RayVector{T}}}}}, P::Polyhedron{T}) where T<:scalar_types - return ( - base_points = _vertices(PointVector{T}, P), - lineality_basis = lineality_space(P) - ) +function minimal_faces(::Type{NamedTuple{(:base_points, :lineality_basis), Tuple{SubObjectIterator{PointVector{T}}, SubObjectIterator{RayVector{T}}}}}, P::Polyhedron{T}) where T<:scalar_types + return ( + base_points = _vertices(PointVector{T}, P), + lineality_basis = lineality_space(P) + ) end -minimal_faces(as::Type{PointVector{T}}, P::Polyhedron{T}) where T<:scalar_types = _vertices(PointVector{T}, P) +minimal_faces(::Type{<:PointVector}, P::Polyhedron) = _vertices(P) @@ -163,13 +163,13 @@ julia> rmlP.lineality_basis ``` """ rays_modulo_lineality(P::Polyhedron{T}) where T<:scalar_types = rays_modulo_lineality(NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}, P) -function rays_modulo_lineality(as::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, P::Polyhedron) where T<:scalar_types - return ( - rays_modulo_lineality = _rays(P), - lineality_basis = lineality_space(P) - ) +function rays_modulo_lineality(::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, P::Polyhedron{T}) where T<:scalar_types + return ( + rays_modulo_lineality = _rays(P), + lineality_basis = lineality_space(P) + ) end -rays_modulo_lineality(as::Type{RayVector}, P::Polyhedron) = _rays(P) +rays_modulo_lineality(::Type{<:RayVector}, P::Polyhedron) = _rays(P) @doc raw""" @@ -197,8 +197,8 @@ julia> vertices(PointVector, P) [1, 2] ``` """ -vertices(as::Type{PointVector{T}}, P::Polyhedron) where T<:scalar_types = lineality_dim(P) == 0 ? _vertices(as, P) : _empty_subobjectiterator(as, P) -_vertices(as::Type{PointVector{T}}, P::Polyhedron) where T<:scalar_types = SubObjectIterator{as}(P, _vertex_polyhedron, length(_vertex_indices(pm_object(P)))) +vertices(as::Type{PointVector{T}}, P::Polyhedron{T}) where T<:scalar_types = lineality_dim(P) == 0 ? _vertices(as, P) : _empty_subobjectiterator(as, P) +_vertices(as::Type{PointVector{T}}, P::Polyhedron{T}) where T<:scalar_types = SubObjectIterator{as}(P, _vertex_polyhedron, length(_vertex_indices(pm_object(P)))) _vertex_polyhedron(U::Type{PointVector{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types = point_vector(coefficient_field(P), @view pm_object(P).VERTICES[_vertex_indices(pm_object(P))[i], 2:end])::U @@ -206,8 +206,8 @@ _point_matrix(::Val{_vertex_polyhedron}, P::Polyhedron; homogenized=false) = @vi _matrix_for_polymake(::Val{_vertex_polyhedron}) = _point_matrix -vertices(::Type{PointVector}, P::Polyhedron{T}) where T<:scalar_types = vertices(PointVector{T}, P) -_vertices(::Type{PointVector}, P::Polyhedron{T}) where T<:scalar_types = _vertices(PointVector{T}, P) +vertices(::Type{<:PointVector}, P::Polyhedron{T}) where T<:scalar_types = vertices(PointVector{T}, P) +_vertices(::Type{<:PointVector}, P::Polyhedron{T}) where T<:scalar_types = _vertices(PointVector{T}, P) _facet_indices(::Val{_vertex_polyhedron}, P::Polyhedron) = pm_object(P).FACETS_THRU_VERTICES[_vertex_indices(pm_object(P)),_facet_indices(pm_object(P))] @@ -311,8 +311,8 @@ julia> rays(RayVector, PO) [0, 1] ``` """ -rays(as::Type{RayVector{T}}, P::Polyhedron) where T<:scalar_types = lineality_dim(P) == 0 ? _rays(as, P) : _empty_subobjectiterator(as, P) -_rays(as::Type{RayVector{T}}, P::Polyhedron) where T<:scalar_types = SubObjectIterator{as}(P, _ray_polyhedron, length(_ray_indices(pm_object(P)))) +rays(as::Type{RayVector{T}}, P::Polyhedron{T}) where T<:scalar_types = lineality_dim(P) == 0 ? _rays(as, P) : _empty_subobjectiterator(as, P) +_rays(as::Type{RayVector{T}}, P::Polyhedron{T}) where T<:scalar_types = SubObjectIterator{as}(P, _ray_polyhedron, length(_ray_indices(pm_object(P)))) _ray_polyhedron(U::Type{RayVector{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(P), @view pm_object(P).VERTICES[_ray_indices(pm_object(P))[i], 2:end])::U @@ -324,8 +324,8 @@ _matrix_for_polymake(::Val{_ray_polyhedron}) = _vector_matrix _incidencematrix(::Val{_ray_polyhedron}) = _facet_indices -rays(::Type{RayVector}, P::Polyhedron{T}) where T<:scalar_types = rays(RayVector{T}, P) -_rays(::Type{RayVector}, P::Polyhedron{T}) where T<:scalar_types = _rays(RayVector{T}, P) +rays(::Type{<:RayVector}, P::Polyhedron{T}) where T<:scalar_types = rays(RayVector{T}, P) +_rays(::Type{<:RayVector}, P::Polyhedron{T}) where T<:scalar_types = _rays(RayVector{T}, P) @doc raw""" rays(P::Polyhedron) @@ -410,17 +410,17 @@ x₃ ≦ 1 """ facets(as::Type{T}, P::Polyhedron{S}) where {R, S<:scalar_types, T<:Union{AffineHalfspace{S}, Pair{R, S}, Polyhedron{S}}} = SubObjectIterator{as}(P, _facet_polyhedron, nfacets(P)) -function _facet_polyhedron(::Type{AffineHalfspace{S}}, P::Polyhedron, i::Base.Integer) where S<:scalar_types - h = decompose_hdata(view(pm_object(P).FACETS, [_facet_index(pm_object(P), i)], :)) - return affine_halfspace(coefficient_field(P), h[1], h[2][]) +function _facet_polyhedron(U::Type{AffineHalfspace{S}}, P::Polyhedron{S}, i::Base.Integer) where S<:scalar_types + h = decompose_hdata(view(pm_object(P).FACETS, [_facet_index(pm_object(P), i)], :)) + return affine_halfspace(coefficient_field(P), h[1], h[2][])::U end -function _facet_polyhedron(::Type{Pair{R, S}}, P::Polyhedron, i::Base.Integer) where {R, S<:scalar_types} - f = coefficient_field(P) - h = decompose_hdata(view(pm_object(P).FACETS, [_facet_index(pm_object(P), i)], :)) - return Pair{R, S}(f.(view(h[1], :, :)), f(h[2][])) # view_broadcast +function _facet_polyhedron(U::Type{Pair{R, S}}, P::Polyhedron{S}, i::Base.Integer) where {R, S<:scalar_types} + f = coefficient_field(P) + h = decompose_hdata(view(pm_object(P).FACETS, [_facet_index(pm_object(P), i)], :)) + return U(f.(view(h[1], :, :)), f(h[2][])) # view_broadcast end -function _facet_polyhedron(::Type{Polyhedron{T}}, P::Polyhedron, i::Base.Integer) where T<:scalar_types - return Polyhedron{T}(Polymake.polytope.facet(pm_object(P), _facet_index(pm_object(P), i)-1), coefficient_field(P)) +function _facet_polyhedron(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types + return Polyhedron{T}(Polymake.polytope.facet(pm_object(P), _facet_index(pm_object(P), i) - 1), coefficient_field(P)) end _affine_inequality_matrix(::Val{_facet_polyhedron}, P::Polyhedron) = -_remove_facet_at_infinity(pm_object(P)) @@ -435,7 +435,7 @@ _vertex_and_ray_indices(::Val{_facet_polyhedron}, P::Polyhedron) = vcat(pm_objec _incidencematrix(::Val{_facet_polyhedron}) = _vertex_and_ray_indices -facets(::Type{Pair}, P::Polyhedron{T}) where T<:scalar_types = facets(Pair{Matrix{T}, T}, P) +facets(::Type{<:Pair}, P::Polyhedron{T}) where T<:scalar_types = facets(Pair{Matrix{T}, T}, P) facets(::Type{Polyhedron}, P::Polyhedron{T}) where T<:scalar_types = facets(Polyhedron{T}, P) @@ -461,9 +461,7 @@ x₃ ≦ 1 """ facets(P::Polyhedron{T}) where T<:scalar_types = facets(AffineHalfspace{T}, P) -facets(::Type{Halfspace}, P::Polyhedron{T}) where T<:scalar_types = facets(AffineHalfspace{T}, P) - -facets(::Type{AffineHalfspace}, P::Polyhedron{T}) where T<:scalar_types = facets(AffineHalfspace{T}, P) +facets(::Type{<:Halfspace}, P::Polyhedron{T}) where T<:scalar_types = facets(AffineHalfspace{T}, P) function _facet_index(P::Polymake.BigObject, i::Base.Integer) i < _facet_at_infinity(P) && return i @@ -562,7 +560,7 @@ julia> normalized_volume(C) 8 ``` """ -normalized_volume(P::Polyhedron{T}) where T<:scalar_types = coefficient_field(P)(factorial(dim(P))*(pm_object(P)).VOLUME) +normalized_volume(P::Polyhedron) = coefficient_field(P)(factorial(dim(P))*(pm_object(P)).VOLUME) @doc raw""" @@ -615,7 +613,7 @@ function lattice_points(P::Polyhedron{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(P, _lattice_point, size(pm_object(P).LATTICE_POINTS_GENERATORS[1], 1)) end -_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).LATTICE_POINTS_GENERATORS[1][i, 2:end])::T +_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron{QQFieldElem}, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).LATTICE_POINTS_GENERATORS[1][i, 2:end])::T _point_matrix(::Val{_lattice_point}, P::Polyhedron; homogenized=false) = @view pm_object(P).LATTICE_POINTS_GENERATORS[1][:, (homogenized ? 1 : 2):end] @@ -646,7 +644,7 @@ function interior_lattice_points(P::Polyhedron{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(P, _interior_lattice_point, size(pm_object(P).INTERIOR_LATTICE_POINTS, 1)) end -_interior_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).INTERIOR_LATTICE_POINTS[i, 2:end])::T +_interior_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron{QQFieldElem}, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).INTERIOR_LATTICE_POINTS[i, 2:end])::T _point_matrix(::Val{_interior_lattice_point}, P::Polyhedron; homogenized=false) = homogenized ? pm_object(P).INTERIOR_LATTICE_POINTS : @view pm_object(P).INTERIOR_LATTICE_POINTS[:, 2:end] @@ -686,7 +684,7 @@ function boundary_lattice_points(P::Polyhedron{QQFieldElem}) return SubObjectIterator{PointVector{ZZRingElem}}(P, _boundary_lattice_point, size(pm_object(P).BOUNDARY_LATTICE_POINTS, 1)) end -_boundary_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).BOUNDARY_LATTICE_POINTS[i, 2:end])::T +_boundary_lattice_point(T::Type{PointVector{ZZRingElem}}, P::Polyhedron{QQFieldElem}, i::Base.Integer) = point_vector(ZZ, @view pm_object(P).BOUNDARY_LATTICE_POINTS[i, 2:end])::T _point_matrix(::Val{_boundary_lattice_point}, P::Polyhedron; homogenized=false) = homogenized ? pm_object(P).BOUNDARY_LATTICE_POINTS : @view pm_object(P).BOUNDARY_LATTICE_POINTS[:, 2:end] @@ -779,9 +777,9 @@ x₄ = 5 """ affine_hull(P::Polyhedron{T}) where T<:scalar_types = SubObjectIterator{AffineHyperplane{T}}(P, _affine_hull, size(pm_object(P).AFFINE_HULL, 1)) -function _affine_hull(::Type{AffineHyperplane{T}}, P::Polyhedron, i::Base.Integer) where T<:scalar_types - h = decompose_hdata(-view(pm_object(P).AFFINE_HULL, [i], :)) - return affine_hyperplane(coefficient_field(P), h[1], h[2][]) +function _affine_hull(U::Type{AffineHyperplane{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types + h = decompose_hdata(-view(pm_object(P).AFFINE_HULL, [i], :)) + return affine_hyperplane(coefficient_field(P), h[1], h[2][])::U end _affine_equation_matrix(::Val{_affine_hull}, P::Polyhedron) = pm_object(P).AFFINE_HULL diff --git a/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl b/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl index 76d2e8f176e8..21205271f14f 100644 --- a/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl +++ b/src/PolyhedralGeometry/SubdivisionOfPoints/properties.jl @@ -29,11 +29,11 @@ julia> points(MOAE) [1, 1, 2] ``` """ -function points(SOP::SubdivisionOfPoints) - return SubObjectIterator{PointVector{QQFieldElem}}(SOP, _point, size(pm_object(SOP).POINTS, 1)) +function points(SOP::SubdivisionOfPoints{T}) where T<:scalar_types + return SubObjectIterator{PointVector{T}}(SOP, _point, size(pm_object(SOP).POINTS, 1)) end -_point(T::Type{PointVector{QQFieldElem}}, SOP::SubdivisionOfPoints, i::Base.Integer) = point_vector(QQ, pm_object(SOP).POINTS[i, 2:end])::T +_point(U::Type{PointVector{T}}, SOP::SubdivisionOfPoints{T}, i::Base.Integer) where T<:scalar_types = point_vector(coefficient_field(SOP), pm_object(SOP).POINTS[i, 2:end])::U _point_matrix(::Val{_point}, SOP::SubdivisionOfPoints; homogenized=false) = pm_object(SOP).POINTS[:, (homogenized ? 1 : 2):end] @@ -185,7 +185,7 @@ julia> min_weights(SOP) 0 ``` """ -min_weights(SOP::SubdivisionOfPoints{T}) where T<:scalar_types = Vector{Int}(pm_object(SOP).MIN_WEIGHTS) +min_weights(SOP::SubdivisionOfPoints) = Vector{Int}(pm_object(SOP).MIN_WEIGHTS) @doc raw""" diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index c94abdf20fe7..2011352ee3ca 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -521,6 +521,8 @@ end _parent_or_coefficient_field(r::Base.RefValue{<:Union{FieldElem, ZZRingElem}}) = parent(r.x) +_parent_or_coefficient_field(v::AbstractVector{T}) where T<:Union{FieldElem, ZZRingElem} = _determine_parent_and_scalar(T, v)[1] + function _promoted_bigobject(::Type{T}, obj::PolyhedralObject{U}) where {T <: scalar_types, U <: scalar_types} T == U ? pm_object(obj) : Polymake.common.convert_to{_scalar_type_to_polymake(T)}(pm_object(obj)) end From f84c984db03145d740606ff5da7c87a9c27f8e2f Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Fri, 21 Jul 2023 23:57:22 +0200 Subject: [PATCH 06/21] re-factored Halfspace and Hyperplane --- src/PolyhedralGeometry/iterators.jl | 144 +++++++++++----------------- 1 file changed, 56 insertions(+), 88 deletions(-) diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index e672b42b342c..7984406b766d 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -43,8 +43,6 @@ for (T, _t) in ((:PointVector, :point_vector), (:RayVector, :ray_vector)) $_t(x::Union{AbstractVector, Base.Integer}) = $_t(QQ, x) - # $T{U}(n::UR) where {U<:scalar_types_extended, UR<:Base.AbstractUnitRange} = $T{U}(undef, length(n)) - function Base.similar(X::$T, ::Type{S}, dims::Dims{1}) where S <: scalar_types_extended return $_t(coefficient_field(X), dims...) end @@ -69,116 +67,86 @@ end ######## Halfspaces and Hyperplanes ################################################################################ -abstract type Halfspace{T} end - -################################################################################ - -@doc raw""" - Halfspace(a, b) - -One halfspace `H(a,b)` is given by a vector `a` and a value `b` such that -$$H(a,b) = \{ x | ax ≤ b \}.$$ -""" -struct AffineHalfspace{T} <: Halfspace{T} - a::Vector{T} - b::T - parent_field::Field - - AffineHalfspace{T}(p::Field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) where T<:scalar_types = new{T}(p.(vec(a)), p(b), p) -end - -halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_halfspace(a, b) -halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_halfspace(f, a, b) +for (h, comp) in (("alfspace", "≤"), ("yperplane", "=")) -invert(H::AffineHalfspace{T}) where T<:scalar_types = AffineHalfspace{T}(coefficient_field(H), -normal_vector(H), -negbias(H)) + Habs = Symbol("H", h) + Haff = Symbol("AffineH", h) + Hlin = Symbol("LinearH", h) + Fabs = Symbol("h", h) + Faff = Symbol("affine_h", h) + Flin = Symbol("linear_h", h) -function affine_halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) - parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) - return AffineHalfspace{scalar_type}(parent_field, a, b) -end + @eval begin -affine_halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) = affine_halfspace(QQ, a, b) + abstract type $Habs{T} end -################################################################################ + # Affine types + struct $Haff{T} <: $Habs{T} + a::MatElem{T} + b::T -struct LinearHalfspace{T} <: Halfspace{T} - a::Vector{T} - parent_field::Field - - LinearHalfspace{T}(p::Field, a::Union{MatElem, AbstractMatrix, AbstractVector}) where T<:scalar_types = new{T}(p.(vec(a)), p) -end - -halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_halfspace(a) -halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_halfspace(f, a) - -invert(H::LinearHalfspace{T}) where T<:scalar_types = LinearHalfspace{T}(coefficient_field(H), -normal_vector(H)) + $Haff{T}(a::MatElem{T}, b::T) where T<:scalar_types = new{T}(a, b) + end -function linear_halfspace(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) - parent_field, scalar_type = _determine_parent_and_scalar(f, a) - return LinearHalfspace{scalar_type}(parent_field, a) -end + $Fabs(a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = $Faff(a, b) + $Fabs(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = $Faff(f, a, b) -linear_halfspace(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_halfspace(QQ, a) + invert(H::$Haff{T}) where T<:scalar_types = $Haff{T}(-H.a, -negbias(H)) -coefficient_field(h::Halfspace) = h.parent_field + @doc """ + """ * string($Faff) * """([p,] a, b) -################################################################################ +Return the """ * string($Haff) * raw""" `H(a,b)`, which is given by a vector `a` and a value `b` such that +$$H(a,b) = \{ x | ax """ * $comp * raw""" b \}.$$ +`p` specifies the `Field` or `Type` of its coefficient. + """ + function $Faff(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) + parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) + mat = matrix(parent_field, 1, length(a), collect(a)) + return $Haff{scalar_type}(mat, parent_field(b)) + end -abstract type Hyperplane{T} end + $Faff(a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) = $Faff(QQ, a, b) -################################################################################ -@doc raw""" - AffineHyperplane(a, b) + # Linear types -One hyperplane `H(a,b)` is given by a vector `a` and a value `b` such that -$$H(a,b) = \{ x | ax = b \}.$$ -""" -struct AffineHyperplane{T} <: Hyperplane{T} - a::Vector{T} - b::T - parent_field::Field + struct $Hlin{T} <: $Habs{T} + a::MatElem{T} - AffineHyperplane{T}(p::Field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b=0) where T<:scalar_types = new{T}(p.(vec(a)), p(b), p) -end - -hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_hyperplane(a, b) -hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b) = affine_hyperplane(f, a, b) - -function affine_hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) - parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) - return AffineHyperplane{scalar_type}(parent_field, a, b) -end + $Hlin{T}(a::MatElem{T}) where T<:scalar_types = new{T}(a) + end -affine_hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) = affine_hyperplane(QQ, a, b) + $Fabs(a::Union{MatElem, AbstractMatrix, AbstractVector}) = $Flin(a) + $Fabs(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) = $Flin(f, a) -################################################################################ + invert(H::$Hlin{T}) where T<:scalar_types = $Hlin{T}(-H.a) -struct LinearHyperplane{T} <: Hyperplane{T} - a::Vector{T} - parent_field::Field - - LinearHyperplane{T}(p::Field, a::Union{MatElem, AbstractMatrix, AbstractVector}) where T<:scalar_types = new{T}(p.(vec(a)), p) -end + @doc """ + """ * string($Flin) * """([p,] a, b) -hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_hyperplane(a) -hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_hyperplane(f, a) +Return the """ * string($Hlin) * raw""" `H(a,b)`, which is given by a vector `a` such that +$$H(a,b) = \{ x | ax """ * $comp * raw""" 0 \}.$$ +`p` specifies the `Field` or `Type` of its coefficient. + """ + function $Flin(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) + parent_field, scalar_type = _determine_parent_and_scalar(f, a) + mat = matrix(parent_field, 1, length(a), collect(a)) + return $Hlin{scalar_type}(mat) + end -function linear_hyperplane(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) - parent_field, scalar_type = _determine_parent_and_scalar(f, a) - return LinearHyperplane{scalar_type}(parent_field, a) -end + $Flin(a::Union{MatElem, AbstractMatrix, AbstractVector}) = $Flin(QQ, a) -linear_hyperplane(a::Union{MatElem, AbstractMatrix, AbstractVector}) = linear_hyperplane(QQ, a) + coefficient_field(h::$Habs) = base_ring(h.a) -coefficient_field(h::Hyperplane) = h.parent_field + end -################################################################################ +end # Field access -negbias(H::Union{AffineHalfspace{T}, AffineHyperplane{T}}) where T<:scalar_types = H.b -negbias(H::Union{LinearHalfspace{T}, LinearHyperplane{T}}) where T<:scalar_types = coefficient_field(H)(0) -normal_vector(H::Union{Halfspace{T}, Hyperplane{T}}) where T <: scalar_types = Vector{T}(H.a) +negbias(H::Union{AffineHalfspace, AffineHyperplane}) = H.b +negbias(H::Union{LinearHalfspace, LinearHyperplane}) = coefficient_field(H)(0) +normal_vector(H::Union{Halfspace, Hyperplane}) = [H.a[1, i] for i in 1:length(H.a)] _ambient_dim(x::Union{Halfspace, Hyperplane}) = length(x.a) From 53149e7c7b07c57af68b9e7241d89c8a52531e44 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Sat, 22 Jul 2023 01:33:15 +0200 Subject: [PATCH 07/21] adjusted and included docs for point_vector, affine_halfspace etc --- docs/src/PolyhedralGeometry/intro.md | 16 ++++++++++++++++ src/PolyhedralGeometry/iterators.jl | 24 ++++++++++++++++++++---- 2 files changed, 36 insertions(+), 4 deletions(-) diff --git a/docs/src/PolyhedralGeometry/intro.md b/docs/src/PolyhedralGeometry/intro.md index 672a0b804e77..c148398e386e 100644 --- a/docs/src/PolyhedralGeometry/intro.md +++ b/docs/src/PolyhedralGeometry/intro.md @@ -72,6 +72,13 @@ any of the corresponding notions below. ### Vectors +There are two specialized `Vector`-like types, `PointVector` and `RayVector`, which commonly are returned by functions from Polyhedral Geometry. These can also be manually constructed: + +```@docs +point_vector +ray_vector +``` + While `RayVector`s can not be used do describe `PointVector`s (and vice versa), matrices are generally allowed. @@ -98,6 +105,15 @@ Type | A `RayVector` corresponds to... ### Halfspaces and Hyperplanes +Similar to points and rays, there are types `AffineHalfspace`, `LinearHalfspace`, `AffineHyperplane` and `LinearHyperplane`: + +```@docs +affine_halfspace +linear_halfspace +affine_hyperplane +linar_hyperplane +``` + These collections allow to mix up affine halfspaces/hyperplanes and their linear counterparts, but note that an error will be produced when trying to convert an affine description with bias not equal to zero to a linear description. diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 7984406b766d..4a2d7888631e 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -63,6 +63,22 @@ for (T, _t) in ((:PointVector, :point_vector), (:RayVector, :ray_vector)) end end +@doc """ + point_vector(p = QQ, v::AbstractVector) + +Return a `PointVector` resembling a point whose coordinates equal the entries of `v`. +`p` specifies the `Field` or `Type` of its coefficient. +""" +point_vector + +@doc """ + ray_vector(p = QQ, v::AbstractVector) + +Return a `RayVector` resembling a ray from the origin through the point whose coordinates equal the entries of `v`. +`p` specifies the `Field` or `Type` of its coefficient. +""" +ray_vector + ################################################################################ ######## Halfspaces and Hyperplanes ################################################################################ @@ -94,9 +110,9 @@ for (h, comp) in (("alfspace", "≤"), ("yperplane", "=")) invert(H::$Haff{T}) where T<:scalar_types = $Haff{T}(-H.a, -negbias(H)) @doc """ - """ * string($Faff) * """([p,] a, b) + """ * string($Faff) * """(p = QQ, a, b) -Return the """ * string($Haff) * raw""" `H(a,b)`, which is given by a vector `a` and a value `b` such that +Return the `""" * string($Haff) * raw"""` `H(a,b)`, which is given by a vector `a` and a value `b` such that $$H(a,b) = \{ x | ax """ * $comp * raw""" b \}.$$ `p` specifies the `Field` or `Type` of its coefficient. """ @@ -123,9 +139,9 @@ $$H(a,b) = \{ x | ax """ * $comp * raw""" b \}.$$ invert(H::$Hlin{T}) where T<:scalar_types = $Hlin{T}(-H.a) @doc """ - """ * string($Flin) * """([p,] a, b) + """ * string($Flin) * """(p = QQ, a, b) -Return the """ * string($Hlin) * raw""" `H(a,b)`, which is given by a vector `a` such that +Return the `""" * string($Hlin) * raw"""` `H(a,b)`, which is given by a vector `a` such that $$H(a,b) = \{ x | ax """ * $comp * raw""" 0 \}.$$ `p` specifies the `Field` or `Type` of its coefficient. """ From c740f638661daa0c77796973e07dbf556d05e3b1 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Tue, 25 Jul 2023 00:43:34 +0200 Subject: [PATCH 08/21] improved docs for Halfspace and Hyperplane types --- docs/src/PolyhedralGeometry/intro.md | 2 +- src/PolyhedralGeometry/iterators.jl | 21 ++++++++++----------- 2 files changed, 11 insertions(+), 12 deletions(-) diff --git a/docs/src/PolyhedralGeometry/intro.md b/docs/src/PolyhedralGeometry/intro.md index c148398e386e..f185306a969b 100644 --- a/docs/src/PolyhedralGeometry/intro.md +++ b/docs/src/PolyhedralGeometry/intro.md @@ -111,7 +111,7 @@ Similar to points and rays, there are types `AffineHalfspace`, `LinearHalfspace` affine_halfspace linear_halfspace affine_hyperplane -linar_hyperplane +linear_hyperplane ``` These collections allow to mix up affine halfspaces/hyperplanes and their linear diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 4a2d7888631e..182c122174e2 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -13,7 +13,6 @@ for (T, _t) in ((:PointVector, :point_vector), (:RayVector, :ray_vector)) end Base.IndexStyle(::Type{<:$T}) = IndexLinear() - Base.getindex(po::$T, i::Base.Integer) = po.p[1, i] function Base.setindex!(po::$T, val, i::Base.Integer) @@ -67,7 +66,7 @@ end point_vector(p = QQ, v::AbstractVector) Return a `PointVector` resembling a point whose coordinates equal the entries of `v`. -`p` specifies the `Field` or `Type` of its coefficient. +`p` specifies the `Field` or `Type` of its coefficients. """ point_vector @@ -75,7 +74,7 @@ point_vector ray_vector(p = QQ, v::AbstractVector) Return a `RayVector` resembling a ray from the origin through the point whose coordinates equal the entries of `v`. -`p` specifies the `Field` or `Type` of its coefficient. +`p` specifies the `Field` or `Type` of its coefficients. """ ray_vector @@ -110,11 +109,11 @@ for (h, comp) in (("alfspace", "≤"), ("yperplane", "=")) invert(H::$Haff{T}) where T<:scalar_types = $Haff{T}(-H.a, -negbias(H)) @doc """ - """ * string($Faff) * """(p = QQ, a, b) + $($Faff)(p = QQ, a, b) -Return the `""" * string($Haff) * raw"""` `H(a,b)`, which is given by a vector `a` and a value `b` such that -$$H(a,b) = \{ x | ax """ * $comp * raw""" b \}.$$ -`p` specifies the `Field` or `Type` of its coefficient. +Return the `$($Haff)` `H(a,b)`, which is given by a vector `a` and a value `b` such that +\$\$H(a,b) = \\{ x | ax $($comp) b \\}.\$\$ +`p` specifies the `Field` or `Type` of its coefficients. """ function $Faff(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}, b = 0) parent_field, scalar_type = _determine_parent_and_scalar(f, a, b) @@ -139,11 +138,11 @@ $$H(a,b) = \{ x | ax """ * $comp * raw""" b \}.$$ invert(H::$Hlin{T}) where T<:scalar_types = $Hlin{T}(-H.a) @doc """ - """ * string($Flin) * """(p = QQ, a, b) + $($Flin)(p = QQ, a, b) -Return the `""" * string($Hlin) * raw"""` `H(a,b)`, which is given by a vector `a` such that -$$H(a,b) = \{ x | ax """ * $comp * raw""" 0 \}.$$ -`p` specifies the `Field` or `Type` of its coefficient. +Return the `$($Hlin)` `H(a)`, which is given by a vector `a` such that +\$\$H(a,b) = \\{ x | ax $($comp) 0 \\}.\$\$ +`p` specifies the `Field` or `Type` of its coefficients. """ function $Flin(f::scalar_type_or_field, a::Union{MatElem, AbstractMatrix, AbstractVector}) parent_field, scalar_type = _determine_parent_and_scalar(f, a) From ee65ab516f7fd1a6137733cbd7f288f082a62243 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Tue, 1 Aug 2023 01:11:49 +0200 Subject: [PATCH 09/21] wip: adjust handling for sub-objects working with QuadraticExtension --- src/PolyhedralGeometry/Cone/properties.jl | 2 +- .../PolyhedralComplex/properties.jl | 5 ++- .../PolyhedralFan/properties.jl | 7 +++- .../Polyhedron/properties.jl | 10 ++++- src/PolyhedralGeometry/helpers.jl | 4 ++ src/PolyhedralGeometry/type_functions.jl | 14 +++---- test/PolyhedralGeometry/scalar_types.jl | 38 +++++++++++++++++++ 7 files changed, 67 insertions(+), 13 deletions(-) diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index 4fa91139b766..b953fa86283d 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -439,7 +439,7 @@ julia> f = facets(Halfspace, c) -x₂ + x₃ ≦ 0 ``` """ -facets(as::Type{<:Union{LinearHalfspace{T}, Cone{T}}}, C::Cone) where T<:scalar_types = SubObjectIterator{as}(C, _facet_cone, nfacets(C)) +facets(as::Type{<:Union{LinearHalfspace{T}, Cone{T}}}, C::Cone{T}) where T<:scalar_types = SubObjectIterator{as}(C, _facet_cone, nfacets(C)) _facet_cone(U::Type{LinearHalfspace{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types = linear_halfspace(coefficient_field(C), -pm_object(C).FACETS[[i], :])::U diff --git a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl index 839cc70ae4cd..71d0deacc653 100644 --- a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl @@ -464,8 +464,11 @@ function _ith_polyhedron( f_ind::Vector{Int64} = Vector{Int64}() ) where {T<:scalar_types} pface = Polymake.row(Polymake.fan.cones_of_dim(pm_object(PC), f_dim), f_ind[i]) + V = pm_object(PC).VERTICES[collect(pface), :] + L = pm_object(PC).LINEALITY_SPACE + PT = _scalar_type_to_polymake(T, eltype(V)) return Polyhedron{T}( - Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(PC).VERTICES[collect(pface), :], LINEALITY_SPACE = pm_object(PC).LINEALITY_SPACE), + Polymake.polytope.Polytope{PT}(VERTICES = V, LINEALITY_SPACE = L), coefficient_field(PC) ) end diff --git a/src/PolyhedralGeometry/PolyhedralFan/properties.jl b/src/PolyhedralGeometry/PolyhedralFan/properties.jl index f8455291c210..eee3303deb06 100644 --- a/src/PolyhedralGeometry/PolyhedralFan/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralFan/properties.jl @@ -164,8 +164,11 @@ function cones(PF::_FanLikeType, cone_dim::Int) return SubObjectIterator{Cone{_get_scalar_type(PF)}}(PF, _cone_of_dim, size(Polymake.fan.cones_of_dim(pm_object(PF), l), 1), (c_dim = l,)) end -function _cone_of_dim(U::Type{Cone{T}}, PF::_FanLikeType{T}, i::Base.Integer; c_dim::Int = 0) where T<:scalar_types - return cone(coefficient_field(PF), pm_object(PF).RAYS[collect(Polymake.row(Polymake.fan.cones_of_dim(pm_object(PF), c_dim), i)),:], pm_object(PF).LINEALITY_SPACE)::U +function _cone_of_dim(::Type{Cone{T}}, PF::_FanLikeType{T}, i::Base.Integer; c_dim::Int = 0) where T<:scalar_types + R = pm_object(PF).RAYS[collect(Polymake.row(Polymake.fan.cones_of_dim(pm_object(PF), c_dim), i)), :] + L = pm_object(PF).LINEALITY_SPACE + PT = _scalar_type_to_polymake(T, eltype(R)) + return Cone{T}(Polymake.polytope.Cone{PT}(RAYS = R, LINEALITY_SPACE = L), coefficient_field(PF)) end _ray_indices(::Val{_cone_of_dim}, PF::_FanLikeType; c_dim::Int = 0) = Polymake.fan.cones_of_dim(pm_object(PF), c_dim) diff --git a/src/PolyhedralGeometry/Polyhedron/properties.jl b/src/PolyhedralGeometry/Polyhedron/properties.jl index 81fd0f74d712..42083b6cb2de 100644 --- a/src/PolyhedralGeometry/Polyhedron/properties.jl +++ b/src/PolyhedralGeometry/Polyhedron/properties.jl @@ -37,7 +37,10 @@ end function _face_polyhedron(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base.Integer; f_dim::Int = -1, f_ind::Vector{Int64} = Vector{Int64}()) where T<:scalar_types pface = Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(P), f_dim))[f_ind[i]] - return Polyhedron{T}(Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(P).VERTICES[collect(pface), :], LINEALITY_SPACE = pm_object(P).LINEALITY_SPACE), coefficient_field(P)) + V = pm_object(P).VERTICES[collect(pface), :] + L = pm_object(P).LINEALITY_SPACE + PT = _scalar_type_to_polymake(T, eltype(V)) + return Polyhedron{T}(Polymake.polytope.Polytope{PT}(VERTICES = V, LINEALITY_SPACE = L), coefficient_field(P)) end function _vertex_indices(::Val{_face_polyhedron}, P::Polyhedron; f_dim = -1, f_ind::Vector{Int64} = Vector{Int64}()) @@ -56,7 +59,10 @@ _incidencematrix(::Val{_face_polyhedron}) = _vertex_and_ray_indices function _face_polyhedron_facet(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base.Integer) where T<:scalar_types pface = pm_object(P).VERTICES_IN_FACETS[_facet_index(pm_object(P), i), :] - return Polyhedron{T}(Polymake.polytope.Polytope{_scalar_type_to_polymake(T)}(VERTICES = pm_object(P).VERTICES[collect(pface), :], LINEALITY_SPACE = pm_object(P).LINEALITY_SPACE), coefficient_field(P)) + V = pm_object(P).VERTICES[collect(pface), :] + L = pm_object(P).LINEALITY_SPACE + PT = _scalar_type_to_polymake(T, eltype(V)) + return Polyhedron{T}(Polymake.polytope.Polytope{PT}(VERTICES = V, LINEALITY_SPACE = L), coefficient_field(P)) end _vertex_indices(::Val{_face_polyhedron_facet}, P::Polyhedron) = vcat(pm_object(P).VERTICES_IN_FACETS[1:(_facet_at_infinity(pm_object(P)) - 1), _vertex_indices(pm_object(P))], pm_object(P).VERTICES_IN_FACETS[(_facet_at_infinity(pm_object(P)) + 1):end, _vertex_indices(pm_object(P))]) diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index 2011352ee3ca..2851e2e7f42f 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -392,6 +392,10 @@ _scalar_type_to_polymake(::Type{QQFieldElem}) = Polymake.Rational _scalar_type_to_polymake(::Type{<:FieldElem}) = Polymake.OscarNumber _scalar_type_to_polymake(::Type{Float64}) = Float64 +_scalar_type_to_polymake(::Type{EmbeddedElem{nf_elem}}, ::Type{Polymake.OscarNumber}) = Polymake.OscarNumber +_scalar_type_to_polymake(::Type{EmbeddedElem{nf_elem}}, ::Type{<:Polymake.QuadraticExtension}) = Polymake.QuadraticExtension{Polymake.Rational} +_scalar_type_to_polymake(::Type{T}, ::Type{U}) where {T<:scalar_types, U} = _scalar_type_to_polymake(T) + #################################################################### # Parent Fields #################################################################### diff --git a/src/PolyhedralGeometry/type_functions.jl b/src/PolyhedralGeometry/type_functions.jl index ac6af70df478..a2830b40469d 100644 --- a/src/PolyhedralGeometry/type_functions.jl +++ b/src/PolyhedralGeometry/type_functions.jl @@ -28,13 +28,13 @@ end function halfspace_matrix_pair(iter::SubObjectIterator{<:Union{Halfspace{T}, Hyperplane{T}, Polyhedron{T}, Cone{T}, Pair{Matrix{T}, T}}}) where T<:scalar_types - try - f = coefficient_field(iter.Obj) - h = affine_matrix_for_polymake(iter) - return (A = matrix(f, h[:, 2:end]), b = Vector{T}(f.(-h[:, 1]))) - catch e - throw(ArgumentError("Halfspace-Matrix-Pair not defined in this context.")) - end + try + f = coefficient_field(iter.Obj) + h = affine_matrix_for_polymake(iter) + return (A = matrix(f, h[:, 2:end]), b = [f(x) for x in -h[:, 1]]) + catch e + throw(ArgumentError("Halfspace-Matrix-Pair not defined in this context.")) + end end for fun in ( diff --git a/test/PolyhedralGeometry/scalar_types.jl b/test/PolyhedralGeometry/scalar_types.jl index 2134693e2489..e0e52add5138 100644 --- a/test/PolyhedralGeometry/scalar_types.jl +++ b/test/PolyhedralGeometry/scalar_types.jl @@ -186,4 +186,42 @@ end end + @testset "QuadraticExtension-templated sub-objects" begin + + j = johnson_solid(2) + + let f = facets(Polyhedron, j) + for i in 1:nfacets(j) + @test facets(f[i])[] in facets(j) + end + @test halfspace_matrix_pair(f) isa NamedTuple{(:A, :b), Tuple{AbstractAlgebra.Generic.MatSpaceElem{EmbeddedElem{nf_elem}}, Vector{EmbeddedElem{nf_elem}}}} + g = halfspace_matrix_pair(f) + affine_halfspace(coefficient_field(j), g.A[1, :], g.b[1]) in facets(j) + end + for n in (1, 2) # faces which are facets have a different access function + let f = faces(j, n) + for i in 1:Int(f_vector(j)[n + 1]) + @test issubset(vertices(f[i]), vertices(j)) + end + end + end + + k = face_fan(j) + let f = cones(k, 3) + for i in 1:Int(f_vector(k)[3]) + @test issubset(rays(f[i]), rays(k)) + end + l = f[1] + end + + l = cone(Polymake.polytope.Cone{Polymake.QuadraticExtension{Polymake.Rational}}(RAYS = Oscar.pm_object(j).VERTICES)) + + let f = facets(Polyhedron, l) + for i in 1:nfacets(l) + @test facets(f[i])[] in facets(l) + end + end + + end + end From 7ad70af068325b55d775f6be4fd74b66694ed567 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Sun, 6 Aug 2023 18:29:46 +0200 Subject: [PATCH 10/21] extended support for Halfspace and Hyperplane --- src/PolyhedralGeometry/iterators.jl | 26 +++++++++++++++++++------- test/PolyhedralGeometry/types.jl | 9 +++++---- 2 files changed, 24 insertions(+), 11 deletions(-) diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 182c122174e2..3d8ca0ad504f 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -165,14 +165,26 @@ normal_vector(H::Union{Halfspace, Hyperplane}) = [H.a[1, i] for i in 1:length(H. _ambient_dim(x::Union{Halfspace, Hyperplane}) = length(x.a) -# TODO: abstract notion of equality -Base.:(==)(x::AffineHalfspace, y::AffineHalfspace) = x.a == y.a && x.b == y.b - -Base.:(==)(x::LinearHalfspace, y::LinearHalfspace) = x.a == y.a - -Base.:(==)(x::AffineHyperplane, y::AffineHyperplane) = x.a == y.a && x.b == y.b +function Base.:(==)(x::Halfspace, y::Halfspace) + ax = normal_vector(x) + ay = normal_vector(y) + ix = findfirst(a -> !iszero(a), ax) + iy = findfirst(a -> !iszero(a), ay) + ix == iy || return false + r = y.a[iy]//x.a[ix] + r > 0 || return false + return (r .* ax == ay) && (r * negbias(x) == negbias(y)) +end -Base.:(==)(x::LinearHyperplane, y::LinearHyperplane) = x.a == y.a +function Base.:(==)(x::Hyperplane, y::Hyperplane) + ax = normal_vector(x) + ay = normal_vector(y) + ix = findfirst(a -> !iszero(a), ax) + iy = findfirst(a -> !iszero(a), ay) + ix == iy || return false + r = y.a[iy]//x.a[ix] + return (r .* ax == ay) && (r * negbias(x) == negbias(y)) +end Base.hash(x::T, h::UInt) where {T<:Union{AffineHalfspace,AffineHyperplane}} = hash((x.a, x.b), hash(T, h)) diff --git a/test/PolyhedralGeometry/types.jl b/test/PolyhedralGeometry/types.jl index 726b21f89e1c..924075566a6f 100644 --- a/test/PolyhedralGeometry/types.jl +++ b/test/PolyhedralGeometry/types.jl @@ -16,12 +16,13 @@ b = [8, 6, 4] NF, sr2 = quadratic_field(2) + ENF, sre2 = Hecke.embedded_field(NF, real_embeddings(NF)[2]) @testset "$T" for (T, fun) in ((PointVector, point_vector), (RayVector, ray_vector)) @test fun(a) isa T{QQFieldElem} - for f in (ZZ, QQ, NF) + for f in (ZZ, QQ, ENF) U = elem_type(f) @testset "$T{$U}" begin @@ -49,7 +50,7 @@ @test_throws BoundsError A[0] @test_throws BoundsError A[4] - for g in (ZZ, QQ, NF) + for g in (ZZ, QQ, ENF) V = elem_type(g) B = fun(g, b) @@ -79,7 +80,7 @@ @testset "$T" for (T, f) in ((AffineHalfspace, affine_halfspace), (AffineHyperplane, affine_hyperplane)) - for p in [QQ, NF] + for p in [QQ, ENF] U = elem_type(p) @test f(p, a, 0) isa T{U} @test f(p, permutedims(a), 0) isa T{U} @@ -107,7 +108,7 @@ @testset "$T" for (T, f) in ((LinearHalfspace, linear_halfspace), (LinearHyperplane, linear_hyperplane)) - for p in [QQ, NF] + for p in [QQ, ENF] U = elem_type(p) @test f(p, a) isa T{U} @test f(p, permutedims(a)) isa T{U} From 152484043e895c97b86f4fa43913781d8eb7f944 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Tue, 8 Aug 2023 23:17:37 +0200 Subject: [PATCH 11/21] adjust handling for sub-objects working with QuadraticExtension(2) --- src/PolyhedralGeometry/Cone/properties.jl | 14 ++++++++++---- test/PolyhedralGeometry/scalar_types.jl | 8 ++++++++ 2 files changed, 18 insertions(+), 4 deletions(-) diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index b953fa86283d..c34394c0385d 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -147,8 +147,11 @@ function faces(C::Cone{T}, face_dim::Int) where T<:scalar_types return SubObjectIterator{Cone{T}}(C, _face_cone, size(Polymake.polytope.faces_of_dim(pm_object(C), n), 1), (f_dim = n,)) end -function _face_cone(U::Type{Cone{T}}, C::Cone{T}, i::Base.Integer; f_dim::Int = 0) where T<:scalar_types - return cone(coefficient_field(C), pm_object(C).RAYS[collect(Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(C), f_dim)[i])), :], pm_object(C).LINEALITY_SPACE)::U +function _face_cone(::Type{Cone{T}}, C::Cone{T}, i::Base.Integer; f_dim::Int = 0) where T<:scalar_types + R = pm_object(C).RAYS[collect(Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(C), f_dim)[i])), :] + L = pm_object(C).LINEALITY_SPACE + PT = _scalar_type_to_polymake(T, eltype(R)) + return Cone{T}(Polymake.polytope.Cone{PT}(RAYS = R, LINEALITY_SPACE = L), coefficient_field(C)) end function _ray_indices(::Val{_face_cone}, C::Cone; f_dim::Int = 0) @@ -156,8 +159,11 @@ function _ray_indices(::Val{_face_cone}, C::Cone; f_dim::Int = 0) return IncidenceMatrix([collect(f[i]) for i in 1:length(f)]) end -function _face_cone_facet(U::Type{Cone{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types - return cone(coefficient_field(C), pm_object(C).RAYS[collect(pm_object(C).RAYS_IN_FACETS[i, :]), :], pm_object(C).LINEALITY_SPACE)::U +function _face_cone_facet(::Type{Cone{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types + R = pm_object(C).RAYS[collect(pm_object(C).RAYS_IN_FACETS[i, :]), :] + L = pm_object(C).LINEALITY_SPACE + PT = _scalar_type_to_polymake(T, eltype(R)) + return Cone{T}(Polymake.polytope.Cone{PT}(RAYS = R, LINEALITY_SPACE = pm_object(C).LINEALITY_SPACE), coefficient_field(C)) end _ray_indices(::Val{_face_cone_facet}, C::Cone) = pm_object(C).RAYS_IN_FACETS diff --git a/test/PolyhedralGeometry/scalar_types.jl b/test/PolyhedralGeometry/scalar_types.jl index e0e52add5138..3b3b81aa01ed 100644 --- a/test/PolyhedralGeometry/scalar_types.jl +++ b/test/PolyhedralGeometry/scalar_types.jl @@ -221,6 +221,14 @@ @test facets(f[i])[] in facets(l) end end + + for n in (2, 3) # faces which are facets have a different access function + let f = faces(l, n) + for i in 1:Int(f_vector(l)[n]) + @test issubset(rays(f[i]), rays(l)) + end + end + end end From a37ddc92c9878feeff14b918499616a52cf61cf4 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Thu, 17 Aug 2023 11:11:39 +0200 Subject: [PATCH 12/21] re-added accidentally skipped test block; further adjustments for QuadraticExtension/OscarNumber compatibility --- .../Polyhedron/constructors.jl | 45 +++++++++++++++---- .../Polyhedron/properties.jl | 2 +- src/PolyhedralGeometry/iterators.jl | 2 +- test/PolyhedralGeometry/polyhedron.jl | 12 ++--- test/PolyhedralGeometry/scalar_types.jl | 2 +- 5 files changed, 46 insertions(+), 17 deletions(-) diff --git a/src/PolyhedralGeometry/Polyhedron/constructors.jl b/src/PolyhedralGeometry/Polyhedron/constructors.jl index 9b175bebb6c4..6e385526e0a3 100644 --- a/src/PolyhedralGeometry/Polyhedron/constructors.jl +++ b/src/PolyhedralGeometry/Polyhedron/constructors.jl @@ -57,17 +57,17 @@ julia> polyhedron(A,b) Polyhedron in ambient dimension 2 ``` """ -polyhedron(f::scalar_type_or_field, A::AnyVecOrMat, b::AbstractVector) = polyhedron(f, (A, b)) +polyhedron(f::scalar_type_or_field, A::AnyVecOrMat, b::AbstractVector; non_redundant::Bool = false) = polyhedron(f, (A, b); non_redundant = non_redundant) -polyhedron(f::scalar_type_or_field, A::AbstractVector, b::Any) = polyhedron(f, ([A], [b])) +polyhedron(f::scalar_type_or_field, A::AbstractVector, b::Any; non_redundant::Bool = false) = polyhedron(f, ([A], [b]); non_redundant = non_redundant) -polyhedron(f::scalar_type_or_field, A::AbstractVector, b::AbstractVector) = polyhedron(f, ([A], b)) +polyhedron(f::scalar_type_or_field, A::AbstractVector, b::AbstractVector; non_redundant::Bool = false) = polyhedron(f, ([A], b); non_redundant = non_redundant) -polyhedron(f::scalar_type_or_field, A::AbstractVector{<:AbstractVector}, b::Any) = polyhedron(f, (A, [b])) +polyhedron(f::scalar_type_or_field, A::AbstractVector{<:AbstractVector}, b::Any; non_redundant::Bool = false) = polyhedron(f, (A, [b]); non_redundant = non_redundant) -polyhedron(f::scalar_type_or_field, A::AbstractVector{<:AbstractVector}, b::AbstractVector) = polyhedron(f, (A, b)) +polyhedron(f::scalar_type_or_field, A::AbstractVector{<:AbstractVector}, b::AbstractVector; non_redundant::Bool = false) = polyhedron(f, (A, b); non_redundant = non_redundant) -polyhedron(f::scalar_type_or_field, A::AnyVecOrMat, b::Any) = polyhedron(f, A, [b]) +polyhedron(f::scalar_type_or_field, A::AnyVecOrMat, b::Any; non_redundant::Bool = false) = polyhedron(f, A, [b]; non_redundant = non_redundant) @doc raw""" polyhedron(::Union{Type{T}, Field}, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) where T<:scalar_types @@ -106,7 +106,7 @@ julia> vertices(P) [0, 0] ``` """ -function polyhedron(f::scalar_type_or_field, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing) +function polyhedron(f::scalar_type_or_field, I::Union{Nothing, AbstractCollection[AffineHalfspace]}, E::Union{Nothing, AbstractCollection[AffineHyperplane]} = nothing; non_redundant::Bool = false) parent_field, scalar_type = _determine_parent_and_scalar(f, I, E) if isnothing(I) || _isempty_halfspace(I) EM = affine_matrix_for_polymake(E) @@ -116,7 +116,11 @@ function polyhedron(f::scalar_type_or_field, I::Union{Nothing, AbstractCollectio EM = isnothing(E) || _isempty_halfspace(E) ? Polymake.Matrix{_scalar_type_to_polymake(scalar_type)}(undef, 0, size(IM, 2)) : affine_matrix_for_polymake(E) end - return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(INEQUALITIES = remove_zero_rows(IM), EQUATIONS = remove_zero_rows(EM)), parent_field) + if non_redundant + return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(FACETS = remove_zero_rows(IM), AFFINE_HULL = remove_zero_rows(EM)), parent_field) + else + return Polyhedron{scalar_type}(Polymake.polytope.Polytope{_scalar_type_to_polymake(scalar_type)}(INEQUALITIES = remove_zero_rows(IM), EQUATIONS = remove_zero_rows(EM)), parent_field) + end end """ @@ -126,7 +130,7 @@ Get the underlying polymake `Polytope`. """ pm_object(P::Polyhedron) = P.pm_polytope -function ==(P0::Polyhedron, P1::Polyhedron) +function ==(P0::Polyhedron{T}, P1::Polyhedron{T}) where T<:scalar_types # TODO: Remove the following 3 lines, see #758 for pair in Iterators.product([P0, P1], ["RAYS", "FACETS"]) Polymake.give(pm_object(pair[1]),pair[2]) @@ -134,6 +138,29 @@ function ==(P0::Polyhedron, P1::Polyhedron) Polymake.polytope.equal_polyhedra(pm_object(P0), pm_object(P1)) end +function Base.:(==)(P0::Polyhedron{EmbeddedElem{nf_elem}}, P1::Polyhedron{EmbeddedElem{nf_elem}}) + R0 = pm_object(P0).RAYS + R1 = pm_object(P1).RAYS + T0 = eltype(R0) + T1 = eltype(R1) + if T0 != T1 + if T0 <: Polymake.QuadraticExtension + Q0 = convex_hull(coefficient_field(P0), point_matrix(vertices(P0)), vector_matrix(rays(P0)), generator_matrix(lineality_space(P0)); non_redundant = true) + Polymake.give(pm_object(Q0), "FACETS") + Polymake.give(pm_object(P1), "FACETS") + return Polymake.polytope.equal_polyhedra(pm_object(Q0), pm_object(P1)) + else + Q1 = convex_hull(coefficient_field(P1), point_matrix(vertices(P1)), vector_matrix(rays(P1)), generator_matrix(lineality_space(P1)); non_redundant = true) + Polymake.give(pm_object(P0), "FACETS") + Polymake.give(pm_object(Q1), "FACETS") + return Polymake.polytope.equal_polyhedra(pm_object(P0), pm_object(Q1)) + end + end + Polymake.give(pm_object(P0), "FACETS") + Polymake.give(pm_object(P1), "FACETS") + Polymake.polytope.equal_polyhedra(pm_object(P0), pm_object(P1)) +end + ### Construct polyhedron from V-data, as the convex hull of points, rays and lineality. @doc raw""" diff --git a/src/PolyhedralGeometry/Polyhedron/properties.jl b/src/PolyhedralGeometry/Polyhedron/properties.jl index 42083b6cb2de..dbfdd878d9bc 100644 --- a/src/PolyhedralGeometry/Polyhedron/properties.jl +++ b/src/PolyhedralGeometry/Polyhedron/properties.jl @@ -1014,7 +1014,7 @@ julia> [1, -2] in PO false ``` """ -Base.in(v::AbstractVector, P::Polyhedron) = Polymake.polytope.contains(pm_object(P), [1; v])::Bool +Base.in(v::AbstractVector, P::Polyhedron) = Polymake.polytope.contains(pm_object(P), coefficient_field(P).([1; v]))::Bool @doc raw""" diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 3d8ca0ad504f..17ea89289384 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -353,7 +353,7 @@ for f in ("_point_matrix", "_vector_matrix", "_generator_matrix") function $M(::Val{_empty_access}, P::PolyhedralObjectUnion; homogenized=false) scalar_regexp = match(r"[^<]*<(.*)>[^>]*", String(Polymake.type_name(pm_object(P)))) typename = scalar_regexp[1] - T = _scalar_type_to_polymake(scalar_type_to_oscar[typename]) + T = typename == "OscarNumber" ? Polymake.OscarNumber : _scalar_type_to_polymake(scalar_type_to_oscar[typename]) return Polymake.Matrix{T}(undef, 0, Polymake.polytope.ambient_dim(pm_object(P)) + homogenized) end end diff --git a/test/PolyhedralGeometry/polyhedron.jl b/test/PolyhedralGeometry/polyhedron.jl index 594edf872176..c753a2be79cb 100644 --- a/test/PolyhedralGeometry/polyhedron.jl +++ b/test/PolyhedralGeometry/polyhedron.jl @@ -483,11 +483,12 @@ for f in (QQ, ENF) end - if T == nf_elem + if T == EmbeddedElem{nf_elem} @testset "Dodecahedron" begin - R, a = quadratic_field(5) + Q, = quadratic_field(5) + R, a = Hecke.embedded_field(Q, real_embeddings(Q)[2]) V = [[1//2, a//4 + 3//4, 0], [-1//2, a//4 + 3//4, 0], @@ -527,13 +528,14 @@ for f in (QQ, ENF) if S == Pair{Matrix{T}, T} @test facets(S, D) == [Pair(A[i], b[i]) for i in 1:12] elseif S == Polyhedron{T} - @test facets(S, D) == [polyhedron(R, A[i], b[i]) for i in 1:12] + # @test facets(S, D) == [polyhedron(R, A[i], b[i]) for i in 1:12] else @test facets(S, D) == [affine_halfspace(R, A[i], b[i]) for i in 1:12] end @test length(facets(S, D)) == 12 @test affine_inequality_matrix(facets(S, D)) == matrix(R, hcat(-b, vcat(A...))) - @test halfspace_matrix_pair(facets(S, D)).A == vcat(A...) && halfspace_matrix_pair(facets(S, D)).b == b + @test halfspace_matrix_pair(facets(S, D)) == (A = matrix(R, vcat(A...)), b = b) + @test ray_indices(facets(S, D)) == IncidenceMatrix(12, 0) @test vertex_indices(facets(S, D)) == IncidenceMatrix([[1, 3, 5, 9, 10], [1, 2, 3, 4, 6], [1, 2, 5, 7, 8], [11, 12, 16, 18, 20], [5, 8, 10, 15, 17], [2, 4, 7, 11, 12], [3, 6, 9, 13, 14], [4, 6, 11, 13, 16], [13, 14, 16, 19, 20], [7, 8, 12, 15, 18], [9, 10, 14, 17, 19], [15, 17, 18, 19, 20]]) @test vertex_and_ray_indices(facets(S, D)) == IncidenceMatrix([[1, 3, 5, 9, 10], [1, 2, 3, 4, 6], [1, 2, 5, 7, 8], [11, 12, 16, 18, 20], [5, 8, 10, 15, 17], [2, 4, 7, 11, 12], [3, 6, 9, 13, 14], [4, 6, 11, 13, 16], [13, 14, 16, 19, 20], [7, 8, 12, 15, 18], [9, 10, 14, 17, 19], [15, 17, 18, 19, 20]]) @@ -559,7 +561,7 @@ for f in (QQ, ENF) @test isempty(rays(D)) @test lineality_dim(D) == 0 @test isempty(lineality_space(D)) - @test faces(D, 0) == convex_hull.(T, V) + # @test faces(D, 0) == convex_hull.(T, V) @test isempty(affine_hull(D)) @test relative_interior_point(D) == [0, 0, 0] diff --git a/test/PolyhedralGeometry/scalar_types.jl b/test/PolyhedralGeometry/scalar_types.jl index 3b3b81aa01ed..9e933224b699 100644 --- a/test/PolyhedralGeometry/scalar_types.jl +++ b/test/PolyhedralGeometry/scalar_types.jl @@ -196,7 +196,7 @@ end @test halfspace_matrix_pair(f) isa NamedTuple{(:A, :b), Tuple{AbstractAlgebra.Generic.MatSpaceElem{EmbeddedElem{nf_elem}}, Vector{EmbeddedElem{nf_elem}}}} g = halfspace_matrix_pair(f) - affine_halfspace(coefficient_field(j), g.A[1, :], g.b[1]) in facets(j) + @test affine_halfspace(coefficient_field(j), g.A[1, :], g.b[1]) in facets(j) end for n in (1, 2) # faces which are facets have a different access function let f = faces(j, n) From eae7cae44bae6832e136393f84835eeed4fd9cf3 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Mon, 21 Aug 2023 11:16:15 +0200 Subject: [PATCH 13/21] removed unnecessary lines from comparison for polytopes; adjusted human readability for halfspace/hyperplane code generation --- .../Polyhedron/constructors.jl | 28 ------------------- src/PolyhedralGeometry/iterators.jl | 17 +++++------ 2 files changed, 9 insertions(+), 36 deletions(-) diff --git a/src/PolyhedralGeometry/Polyhedron/constructors.jl b/src/PolyhedralGeometry/Polyhedron/constructors.jl index 6e385526e0a3..9911576131ef 100644 --- a/src/PolyhedralGeometry/Polyhedron/constructors.jl +++ b/src/PolyhedralGeometry/Polyhedron/constructors.jl @@ -131,37 +131,9 @@ Get the underlying polymake `Polytope`. pm_object(P::Polyhedron) = P.pm_polytope function ==(P0::Polyhedron{T}, P1::Polyhedron{T}) where T<:scalar_types - # TODO: Remove the following 3 lines, see #758 - for pair in Iterators.product([P0, P1], ["RAYS", "FACETS"]) - Polymake.give(pm_object(pair[1]),pair[2]) - end Polymake.polytope.equal_polyhedra(pm_object(P0), pm_object(P1)) end -function Base.:(==)(P0::Polyhedron{EmbeddedElem{nf_elem}}, P1::Polyhedron{EmbeddedElem{nf_elem}}) - R0 = pm_object(P0).RAYS - R1 = pm_object(P1).RAYS - T0 = eltype(R0) - T1 = eltype(R1) - if T0 != T1 - if T0 <: Polymake.QuadraticExtension - Q0 = convex_hull(coefficient_field(P0), point_matrix(vertices(P0)), vector_matrix(rays(P0)), generator_matrix(lineality_space(P0)); non_redundant = true) - Polymake.give(pm_object(Q0), "FACETS") - Polymake.give(pm_object(P1), "FACETS") - return Polymake.polytope.equal_polyhedra(pm_object(Q0), pm_object(P1)) - else - Q1 = convex_hull(coefficient_field(P1), point_matrix(vertices(P1)), vector_matrix(rays(P1)), generator_matrix(lineality_space(P1)); non_redundant = true) - Polymake.give(pm_object(P0), "FACETS") - Polymake.give(pm_object(Q1), "FACETS") - return Polymake.polytope.equal_polyhedra(pm_object(P0), pm_object(Q1)) - end - end - Polymake.give(pm_object(P0), "FACETS") - Polymake.give(pm_object(P1), "FACETS") - Polymake.polytope.equal_polyhedra(pm_object(P0), pm_object(P1)) -end - - ### Construct polyhedron from V-data, as the convex hull of points, rays and lineality. @doc raw""" convex_hull([::Union{Type{T}, Field} = QQFieldElem,] V [, R [, L]]; non_redundant::Bool = false) diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 17ea89289384..0e797f4879d6 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -82,14 +82,15 @@ ray_vector ######## Halfspaces and Hyperplanes ################################################################################ -for (h, comp) in (("alfspace", "≤"), ("yperplane", "=")) - - Habs = Symbol("H", h) - Haff = Symbol("AffineH", h) - Hlin = Symbol("LinearH", h) - Fabs = Symbol("h", h) - Faff = Symbol("affine_h", h) - Flin = Symbol("linear_h", h) +for (h, comp) in (("halfspace", "≤"), ("hyperplane", "=")) + + H = uppercasefirst(h) + Habs = Symbol(H) + Haff = Symbol("Affine", H) + Hlin = Symbol("Linear", H) + Fabs = Symbol(h) + Faff = Symbol("affine_", h) + Flin = Symbol("linear_", h) @eval begin From f25a7086621f0afe75d25ad629d981022b95218f Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Thu, 24 Aug 2023 13:10:09 +0200 Subject: [PATCH 14/21] convert Polymake.BigObjects that are Polytopes over QuadraticExtension to Polytopes over OscarNumber during construction using polyhedron(::Polymake.BigObject) method --- .../Polyhedron/constructors.jl | 11 ++++++-- src/PolyhedralGeometry/helpers.jl | 25 +++++++++++++++++++ 2 files changed, 34 insertions(+), 2 deletions(-) diff --git a/src/PolyhedralGeometry/Polyhedron/constructors.jl b/src/PolyhedralGeometry/Polyhedron/constructors.jl index 9911576131ef..30e011259b37 100644 --- a/src/PolyhedralGeometry/Polyhedron/constructors.jl +++ b/src/PolyhedralGeometry/Polyhedron/constructors.jl @@ -30,8 +30,15 @@ polyhedron(A) = polyhedron(QQFieldElem, A) Construct a `Polyhedron` corresponding to a `Polymake.BigObject` of type `Polytope`. Scalar type and parent field will be detected automatically. To improve type stability and performance, please use [`Polyhedron{T}(p::Polymake.BigObject, f::Field) where T<:scalar_types`](@ref) instead, where possible. """ function polyhedron(p::Polymake.BigObject) - T, f = _detect_scalar_and_field(Polyhedron, p) - return Polyhedron{T}(p, f) + T, f = _detect_scalar_and_field(Polyhedron, p) + if T == EmbeddedElem{nf_elem} && Hecke.isquadratic_type(number_field(f))[1] + scalar_regexp = match(r"[^<]*<(.*)>[^>]*", String(Polymake.type_name(p))) + typename = scalar_regexp[1] + if typename == "QuadraticExtension" + p = _polyhedron_qe_to_on(p, f) + end + end + return Polyhedron{T}(p, f) end @doc raw""" diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index 2851e2e7f42f..81013fdc12bd 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -540,3 +540,28 @@ end function Polymake._fieldelem_is_rational(e::Hecke.EmbeddedNumFieldElem) return is_rational(e) end + +# convert a Polymake.BigObject's scalar from QuadraticExtension to OscarNumber (Polytope only) + +function _polyhedron_qe_to_on(x::Polymake.BigObject, f::Field) + res = Polymake.polytope.Polytope{Polymake.OscarNumber}() + for pn in Polymake.list_properties(x) + prop = Polymake.give(x, pn) + Polymake.take(res, string(pn), _property_qe_to_on(prop, f)) + end + return res +end + +_property_qe_to_on(x::Polymake.BigObject, f::Field) = Polymake.BigObject(Polymake.bigobject_type(x), x) + +_property_qe_to_on(x::Polymake.PropertyValue, f::Field) = x + +_property_qe_to_on(x::Polymake.QuadraticExtension{Polymake.Rational}, f::Field) = f(x) + +function _property_qe_to_on(x, f::Field) + if hasmethod(length, (typeof(x),)) && eltype(x) <: Polymake.QuadraticExtension{Polymake.Rational} + return f.(x) + else + return x + end +end From f9bc2823bf49a5d30ffad6aad7b674f23eb6b028 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Thu, 24 Aug 2023 13:42:29 +0200 Subject: [PATCH 15/21] wip: adjust functionality and tests for reduced support of QuadraticExtension --- src/PolyhedralGeometry/Cone/properties.jl | 2 +- src/PolyhedralGeometry/Polyhedron/properties.jl | 4 ++-- src/PolyhedralGeometry/helpers.jl | 4 ---- test/PolyhedralGeometry/polyhedron.jl | 12 ++++++++---- 4 files changed, 11 insertions(+), 11 deletions(-) diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index c34394c0385d..ce134dda07a5 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -162,7 +162,7 @@ end function _face_cone_facet(::Type{Cone{T}}, C::Cone{T}, i::Base.Integer) where T<:scalar_types R = pm_object(C).RAYS[collect(pm_object(C).RAYS_IN_FACETS[i, :]), :] L = pm_object(C).LINEALITY_SPACE - PT = _scalar_type_to_polymake(T, eltype(R)) + PT = _scalar_type_to_polymake(T) return Cone{T}(Polymake.polytope.Cone{PT}(RAYS = R, LINEALITY_SPACE = pm_object(C).LINEALITY_SPACE), coefficient_field(C)) end diff --git a/src/PolyhedralGeometry/Polyhedron/properties.jl b/src/PolyhedralGeometry/Polyhedron/properties.jl index dbfdd878d9bc..e2b074d9dc69 100644 --- a/src/PolyhedralGeometry/Polyhedron/properties.jl +++ b/src/PolyhedralGeometry/Polyhedron/properties.jl @@ -39,7 +39,7 @@ function _face_polyhedron(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base.Integ pface = Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(P), f_dim))[f_ind[i]] V = pm_object(P).VERTICES[collect(pface), :] L = pm_object(P).LINEALITY_SPACE - PT = _scalar_type_to_polymake(T, eltype(V)) + PT = _scalar_type_to_polymake(T) return Polyhedron{T}(Polymake.polytope.Polytope{PT}(VERTICES = V, LINEALITY_SPACE = L), coefficient_field(P)) end @@ -61,7 +61,7 @@ function _face_polyhedron_facet(::Type{Polyhedron{T}}, P::Polyhedron{T}, i::Base pface = pm_object(P).VERTICES_IN_FACETS[_facet_index(pm_object(P), i), :] V = pm_object(P).VERTICES[collect(pface), :] L = pm_object(P).LINEALITY_SPACE - PT = _scalar_type_to_polymake(T, eltype(V)) + PT = _scalar_type_to_polymake(T) return Polyhedron{T}(Polymake.polytope.Polytope{PT}(VERTICES = V, LINEALITY_SPACE = L), coefficient_field(P)) end diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index 81013fdc12bd..ace95f853961 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -392,10 +392,6 @@ _scalar_type_to_polymake(::Type{QQFieldElem}) = Polymake.Rational _scalar_type_to_polymake(::Type{<:FieldElem}) = Polymake.OscarNumber _scalar_type_to_polymake(::Type{Float64}) = Float64 -_scalar_type_to_polymake(::Type{EmbeddedElem{nf_elem}}, ::Type{Polymake.OscarNumber}) = Polymake.OscarNumber -_scalar_type_to_polymake(::Type{EmbeddedElem{nf_elem}}, ::Type{<:Polymake.QuadraticExtension}) = Polymake.QuadraticExtension{Polymake.Rational} -_scalar_type_to_polymake(::Type{T}, ::Type{U}) where {T<:scalar_types, U} = _scalar_type_to_polymake(T) - #################################################################### # Parent Fields #################################################################### diff --git a/test/PolyhedralGeometry/polyhedron.jl b/test/PolyhedralGeometry/polyhedron.jl index c753a2be79cb..3a6f85d5278d 100644 --- a/test/PolyhedralGeometry/polyhedron.jl +++ b/test/PolyhedralGeometry/polyhedron.jl @@ -487,8 +487,14 @@ for f in (QQ, ENF) @testset "Dodecahedron" begin - Q, = quadratic_field(5) - R, a = Hecke.embedded_field(Q, real_embeddings(Q)[2]) + D = polyhedron(Polymake.polytope.dodecahedron()) + R = coefficient_field(D) + NF = number_field(R) + let isq = Hecke.isquadratic_type(NF) + @test isq[1] + @test isq[2] == 5 + end + a = R(gens(NF)[]) V = [[1//2, a//4 + 3//4, 0], [-1//2, a//4 + 3//4, 0], @@ -511,9 +517,7 @@ for f in (QQ, ENF) [1//2, -a//4 - 3//4, 0], [-1//2, -a//4 - 3//4, 0]] - D = Polyhedron{T}(Polymake.polytope.dodecahedron(), R) @test D isa Polyhedron{T} - @test polyhedron(Polymake.polytope.dodecahedron()) == D @test nvertices(D) == 20 @test vertices(D) == V From d10e38495722556cabeed1b3335bf79643b9b19b Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Fri, 8 Sep 2023 23:35:43 +0200 Subject: [PATCH 16/21] finished first rebase - tests and doctests successful --- .../PolyhedralComplex/properties.jl | 2 +- .../PolyhedralFan/properties.jl | 18 +++++++------ src/PolyhedralGeometry/helpers.jl | 27 ++++++++++++++++++- src/PolyhedralGeometry/type_functions.jl | 2 +- test/PolyhedralGeometry/polyhedron.jl | 3 --- test/PolyhedralGeometry/scalar_types.jl | 18 +------------ 6 files changed, 39 insertions(+), 31 deletions(-) diff --git a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl index 71d0deacc653..96122284a7c5 100644 --- a/src/PolyhedralGeometry/PolyhedralComplex/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralComplex/properties.jl @@ -466,7 +466,7 @@ function _ith_polyhedron( pface = Polymake.row(Polymake.fan.cones_of_dim(pm_object(PC), f_dim), f_ind[i]) V = pm_object(PC).VERTICES[collect(pface), :] L = pm_object(PC).LINEALITY_SPACE - PT = _scalar_type_to_polymake(T, eltype(V)) + PT = _scalar_type_to_polymake(T) return Polyhedron{T}( Polymake.polytope.Polytope{PT}(VERTICES = V, LINEALITY_SPACE = L), coefficient_field(PC) diff --git a/src/PolyhedralGeometry/PolyhedralFan/properties.jl b/src/PolyhedralGeometry/PolyhedralFan/properties.jl index eee3303deb06..191c1a884a8e 100644 --- a/src/PolyhedralGeometry/PolyhedralFan/properties.jl +++ b/src/PolyhedralGeometry/PolyhedralFan/properties.jl @@ -46,13 +46,15 @@ julia> matrix(QQ, rays(NF)) rays(PF::_FanLikeType) = lineality_dim(PF) == 0 ? _rays(PF) : _empty_subobjectiterator(RayVector{_get_scalar_type(PF)}, PF) _rays(PF::_FanLikeType) = SubObjectIterator{RayVector{_get_scalar_type(PF)}}(PF, _ray_fan, _nrays(PF)) -_ray_fan(U::Type{RayVector{T}}, PF::_FanLikeType{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(PF), view(pm_object(PF).RAYS, i, :))::U +_ray_fan(U::Type{RayVector{T}}, PF::_FanLikeType, i::Base.Integer) where {T<:scalar_types} = + ray_vector(coefficient_field(PF), view(pm_object(PF).RAYS, i, :))::U _vector_matrix(::Val{_ray_fan}, PF::_FanLikeType; homogenized=false) = homogenized ? homogenize(pm_object(PF).RAYS, 0) : pm_object(PF).RAYS _matrix_for_polymake(::Val{_ray_fan}) = _vector_matrix -_maximal_cone(::Type{Cone{T}}, PF::_FanLikeType{T}, i::Base.Integer) where T<:scalar_types = Cone{T}(Polymake.fan.cone(pm_object(PF), i - 1), coefficient_field(PF)) +_maximal_cone(::Type{Cone{T}}, PF::_FanLikeType, i::Base.Integer) where {T<:scalar_types} = + Cone{T}(Polymake.fan.cone(pm_object(PF), i - 1), coefficient_field(PF)) @doc raw""" @@ -87,13 +89,13 @@ julia> rays(NF) 0-element SubObjectIterator{RayVector{QQFieldElem}} ``` """ -rays_modulo_lineality(F::_FanLikeType{T}) where {T<:scalar_types} = +rays_modulo_lineality(F::_FanLikeType) = rays_modulo_lineality( - NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}, + NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{_get_scalar_type(F)}}, SubObjectIterator{RayVector{_get_scalar_type(F)}}}}, F ) -function rays_modulo_lineality(::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, F::_FanLikeType{T}) where T<:scalar_types +function rays_modulo_lineality(::Type{NamedTuple{(:rays_modulo_lineality, :lineality_basis), Tuple{SubObjectIterator{RayVector{T}}, SubObjectIterator{RayVector{T}}}}}, F::_FanLikeType) where {T<:scalar_types} return ( rays_modulo_lineality = _rays(F), lineality_basis = lineality_space(F) @@ -164,10 +166,10 @@ function cones(PF::_FanLikeType, cone_dim::Int) return SubObjectIterator{Cone{_get_scalar_type(PF)}}(PF, _cone_of_dim, size(Polymake.fan.cones_of_dim(pm_object(PF), l), 1), (c_dim = l,)) end -function _cone_of_dim(::Type{Cone{T}}, PF::_FanLikeType{T}, i::Base.Integer; c_dim::Int = 0) where T<:scalar_types +function _cone_of_dim(::Type{Cone{T}}, PF::_FanLikeType, i::Base.Integer; c_dim::Int = 0) where T<:scalar_types R = pm_object(PF).RAYS[collect(Polymake.row(Polymake.fan.cones_of_dim(pm_object(PF), c_dim), i)), :] L = pm_object(PF).LINEALITY_SPACE - PT = _scalar_type_to_polymake(T, eltype(R)) + PT = _scalar_type_to_polymake(T) return Cone{T}(Polymake.polytope.Cone{PT}(RAYS = R, LINEALITY_SPACE = L), coefficient_field(PF)) end @@ -394,7 +396,7 @@ julia> lineality_space(PF) """ lineality_space(PF::_FanLikeType) = SubObjectIterator{RayVector{_get_scalar_type(PF)}}(PF, _lineality_fan, lineality_dim(PF)) -_lineality_fan(U::Type{RayVector{T}}, PF::_FanLikeType{T}, i::Base.Integer) where T<:scalar_types = ray_vector(coefficient_field(PF), view(pm_object(PF).LINEALITY_SPACE, i, :))::U +_lineality_fan(U::Type{RayVector{T}}, PF::_FanLikeType, i::Base.Integer) where {T<:scalar_types} = ray_vector(coefficient_field(PF), view(pm_object(PF).LINEALITY_SPACE, i, :))::U _generator_matrix(::Val{_lineality_fan}, PF::_FanLikeType; homogenized=false) = homogenized ? homogenize(pm_object(PF).LINEALITY_SPACE, 0) : pm_object(PF).LINEALITY_SPACE diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index ace95f853961..f75f440dc9bd 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -498,9 +498,34 @@ function _detect_default_field(::Type{T}, p::Polymake.BigObject) where T<:FieldE throw(ArgumentError("BigObject does not contain information about a parent Field")) end +function _detect_wrapped_type_and_field(p::Polymake.BigObject) + # we only want to check existing properties + f = x -> Polymake.exists(p, string(x)) + propnames = intersect(propertynames(p), [:INPUT_RAYS, :POINTS, :RAYS, :VERTICES, :VECTORS, :INPUT_LINEALITY, :LINEALITY_SPACE, :FACETS, :INEQUALITIES, :EQUATIONS, :LINEAR_SPAN, :AFFINE_HULL]) + i = findfirst(f, propnames) + # find first OscarNumber wrapping a FieldElem + while !isnothing(i) + prop = getproperty(p, propnames[i]) + for el in prop + on = Polymake.unwrap(el) + if on isa FieldElem + f = parent(on) + T = elem_type(f) + return (T, f) + end + end + i = findnext(f, propnames, i + 1) + end + throw(ArgumentError("BigObject does not contain information about a parent Field")) +end + function _detect_scalar_and_field(::Type{U}, p::Polymake.BigObject) where U<:PolyhedralObject - T = detect_scalar_type(U, p) + T = detect_scalar_type(U, p) + if isnothing(T) + return _detect_wrapped_type_and_field(p) + else return (T, _detect_default_field(T, p)) + end end # promotion helpers diff --git a/src/PolyhedralGeometry/type_functions.jl b/src/PolyhedralGeometry/type_functions.jl index a2830b40469d..5ed4480fa391 100644 --- a/src/PolyhedralGeometry/type_functions.jl +++ b/src/PolyhedralGeometry/type_functions.jl @@ -6,7 +6,7 @@ function detect_scalar_type(n::Type{T}, p::Polymake.BigObject) where T<:Union{Polyhedron, Cone, PolyhedralFan, SubdivisionOfPoints, PolyhedralComplex} scalar_regexp = match(r"[^<]*<(.*)>[^>]*", String(Polymake.type_name(p))) typename = scalar_regexp[1] - return scalar_type_to_oscar[typename] + return typename == "OscarNumber" ? nothing : scalar_type_to_oscar[typename] end scalar_type(::Union{Polyhedron{T}, Cone{T}, Hyperplane{T}, Halfspace{T}}) where T<:scalar_types = T diff --git a/test/PolyhedralGeometry/polyhedron.jl b/test/PolyhedralGeometry/polyhedron.jl index 3a6f85d5278d..33b63f15f2bd 100644 --- a/test/PolyhedralGeometry/polyhedron.jl +++ b/test/PolyhedralGeometry/polyhedron.jl @@ -568,9 +568,6 @@ for f in (QQ, ENF) # @test faces(D, 0) == convex_hull.(T, V) @test isempty(affine_hull(D)) @test relative_interior_point(D) == [0, 0, 0] - - @test platonic_solid("dodecahedron") == D - end end diff --git a/test/PolyhedralGeometry/scalar_types.jl b/test/PolyhedralGeometry/scalar_types.jl index 9e933224b699..3ca61579c473 100644 --- a/test/PolyhedralGeometry/scalar_types.jl +++ b/test/PolyhedralGeometry/scalar_types.jl @@ -70,7 +70,7 @@ @test isq[2] == 2 end let j = johnson_solid(1) - jj = polyhedron(Polymake.polytope.Polytope{Polymake.QuadraticExtension{Polymake.Rational}}(POINTS=Oscar.pm_object(j).VERTICES)) + jj = polyhedron(Polymake.polytope.Polytope{Polymake.OscarNumber}(POINTS=Oscar.pm_object(j).VERTICES)) @test number_field(coefficient_field(j)) == number_field(coefficient_field(jj)) end end @@ -213,22 +213,6 @@ end l = f[1] end - - l = cone(Polymake.polytope.Cone{Polymake.QuadraticExtension{Polymake.Rational}}(RAYS = Oscar.pm_object(j).VERTICES)) - - let f = facets(Polyhedron, l) - for i in 1:nfacets(l) - @test facets(f[i])[] in facets(l) - end - end - - for n in (2, 3) # faces which are facets have a different access function - let f = faces(l, n) - for i in 1:Int(f_vector(l)[n]) - @test issubset(rays(f[i]), rays(l)) - end - end - end end From 1ac23c4b9b628135658e7fa8a7a4e935b94ff85d Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Mon, 11 Sep 2023 23:23:17 +0200 Subject: [PATCH 17/21] further adjusted to recent changes in master --- src/AlgebraicGeometry/ToricVarieties/Proj/constructors.jl | 4 ++-- src/PolyhedralGeometry/iterators.jl | 2 +- test/PolyhedralGeometry/scalar_types.jl | 4 +++- 3 files changed, 6 insertions(+), 4 deletions(-) diff --git a/src/AlgebraicGeometry/ToricVarieties/Proj/constructors.jl b/src/AlgebraicGeometry/ToricVarieties/Proj/constructors.jl index 90f7fe73a8f1..f7ac0894c643 100644 --- a/src/AlgebraicGeometry/ToricVarieties/Proj/constructors.jl +++ b/src/AlgebraicGeometry/ToricVarieties/Proj/constructors.jl @@ -74,14 +74,14 @@ function _proj_and_total_space(is_proj::Bool, E::Vector{T}) where T <: Union{Tor end end - total_rays_gens = vcat([modified_ray_gens[ray] for ray in rays(v)], [vcat(RayVector(zeros(Int64, dim(v))), l[i]) for i in eachindex(E)]) + total_rays_gens = vcat([modified_ray_gens[ray] for ray in rays(v)], [vcat(ray_vector(zeros(Int64, dim(v))), l[i]) for i in eachindex(E)]) return normal_toric_variety(polyhedral_fan(total_rays_gens, IncidenceMatrix(new_maximal_cones))) end function _m_sigma(sigma::Cone{QQFieldElem}, pol_sigma::Cone{QQFieldElem}, D::Union{ToricDivisor, ToricLineBundle})::RayVector{QQFieldElem} - ans = RayVector(zeros(QQFieldElem, dim(sigma))) + ans = ray_vector(zeros(QQFieldElem, dim(sigma))) coeff = coefficients(isa(D, ToricDivisor) ? D : toric_divisor(D)) dual_ray = QQFieldElem[] diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index 0e797f4879d6..af03742d2ebc 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -57,7 +57,7 @@ for (T, _t) in ((:PointVector, :point_vector), (:RayVector, :ray_vector)) Base.:*(k::scalar_types_extended, po::$T) = k .* po - Base.:*(A::MatElem, v::$T) = A * v.p + Base.:*(A::MatElem, v::$T) = A * transpose(v.p) end end diff --git a/test/PolyhedralGeometry/scalar_types.jl b/test/PolyhedralGeometry/scalar_types.jl index 3ca61579c473..bca6a0799252 100644 --- a/test/PolyhedralGeometry/scalar_types.jl +++ b/test/PolyhedralGeometry/scalar_types.jl @@ -191,8 +191,10 @@ j = johnson_solid(2) let f = facets(Polyhedron, j) + fj = normal_vector.(facets(j)) + fj = [fj; -fj] for i in 1:nfacets(j) - @test facets(f[i])[] in facets(j) + @test normal_vector(affine_hull(f[i])[]) in fj end @test halfspace_matrix_pair(f) isa NamedTuple{(:A, :b), Tuple{AbstractAlgebra.Generic.MatSpaceElem{EmbeddedElem{nf_elem}}, Vector{EmbeddedElem{nf_elem}}}} g = halfspace_matrix_pair(f) From 0a60487fa5eb5291333db886993a7e51f12fee71 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Thu, 14 Sep 2023 13:51:44 +0200 Subject: [PATCH 18/21] fixed error from borcherds_method; small fixes --- src/AlgebraicGeometry/Surfaces/K3Auto.jl | 2 +- src/PolyhedralGeometry/Cone/properties.jl | 4 +--- 2 files changed, 2 insertions(+), 4 deletions(-) diff --git a/src/AlgebraicGeometry/Surfaces/K3Auto.jl b/src/AlgebraicGeometry/Surfaces/K3Auto.jl index e2806317765a..c1a00dde74fa 100644 --- a/src/AlgebraicGeometry/Surfaces/K3Auto.jl +++ b/src/AlgebraicGeometry/Surfaces/K3Auto.jl @@ -1594,7 +1594,7 @@ function span_in_S(L, S, weyl) if N==0 M = zero_matrix(QQ, 0, rank(S)) else - M = reduce(vcat, (matrix(QQ, 1, rank(S), v.a) for v in spanC)) + M = linear_equation_matrix(spanC) end k, K = kernel(M) gensN = transpose(K)[1:k,:] diff --git a/src/PolyhedralGeometry/Cone/properties.jl b/src/PolyhedralGeometry/Cone/properties.jl index ce134dda07a5..9fcb39394b45 100644 --- a/src/PolyhedralGeometry/Cone/properties.jl +++ b/src/PolyhedralGeometry/Cone/properties.jl @@ -150,7 +150,7 @@ end function _face_cone(::Type{Cone{T}}, C::Cone{T}, i::Base.Integer; f_dim::Int = 0) where T<:scalar_types R = pm_object(C).RAYS[collect(Polymake.to_one_based_indexing(Polymake.polytope.faces_of_dim(pm_object(C), f_dim)[i])), :] L = pm_object(C).LINEALITY_SPACE - PT = _scalar_type_to_polymake(T, eltype(R)) + PT = _scalar_type_to_polymake(T) return Cone{T}(Polymake.polytope.Cone{PT}(RAYS = R, LINEALITY_SPACE = L), coefficient_field(C)) end @@ -463,8 +463,6 @@ facets(C::Cone{T}) where T<:scalar_types = facets(LinearHalfspace{T}, C) facets(::Type{<:Halfspace}, C::Cone{T}) where T<:scalar_types = facets(LinearHalfspace{T}, C) -facets(::Type{<:AffineHalfspace}, C::Cone{T}) where T<:scalar_types = facets(AffineHalfspace{T}, C) - facets(::Type{Cone}, C::Cone{T}) where T<:scalar_types = facets(Cone{T}, C) @doc raw""" From 31bb3d4a8f6fd67790b4f50d7edc3d9ff24b979a Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Sat, 16 Sep 2023 11:19:12 +0200 Subject: [PATCH 19/21] applied suggested changes from PR on github --- .../Polyhedron/constructors.jl | 8 +--- src/PolyhedralGeometry/helpers.jl | 41 ++++++++----------- src/PolyhedralGeometry/iterators.jl | 3 +- src/PolyhedralGeometry/type_functions.jl | 3 +- test/PolyhedralGeometry/polyhedron.jl | 4 +- test/PolyhedralGeometry/types.jl | 2 - 6 files changed, 24 insertions(+), 37 deletions(-) diff --git a/src/PolyhedralGeometry/Polyhedron/constructors.jl b/src/PolyhedralGeometry/Polyhedron/constructors.jl index 30e011259b37..d5347eeb6780 100644 --- a/src/PolyhedralGeometry/Polyhedron/constructors.jl +++ b/src/PolyhedralGeometry/Polyhedron/constructors.jl @@ -31,12 +31,8 @@ Construct a `Polyhedron` corresponding to a `Polymake.BigObject` of type `Polyto """ function polyhedron(p::Polymake.BigObject) T, f = _detect_scalar_and_field(Polyhedron, p) - if T == EmbeddedElem{nf_elem} && Hecke.isquadratic_type(number_field(f))[1] - scalar_regexp = match(r"[^<]*<(.*)>[^>]*", String(Polymake.type_name(p))) - typename = scalar_regexp[1] - if typename == "QuadraticExtension" - p = _polyhedron_qe_to_on(p, f) - end + if T == EmbeddedElem{nf_elem} && Hecke.isquadratic_type(number_field(f))[1] && Polymake.bigobject_eltype(p) == "QuadraticExtension" + p = _polyhedron_qe_to_on(p, f) end return Polyhedron{T}(p, f) end diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index f75f440dc9bd..ab406a4a859d 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -381,8 +381,9 @@ _get_scalar_type(::NormalToricVarietyType) = QQFieldElem const scalar_types = Union{FieldElem, Float64} const scalar_type_to_oscar = Dict{String, Type}([("Rational", QQFieldElem), - ("QuadraticExtension", Hecke.EmbeddedNumFieldElem{nf_elem}), - ("Float", Float64)]) + ("QuadraticExtension", Hecke.EmbeddedNumFieldElem{nf_elem}), + ("QuadraticExtension", Hecke.EmbeddedNumFieldElem{nf_elem}), + ("Float", Float64)]) const scalar_types_extended = Union{scalar_types, ZZRingElem} @@ -480,32 +481,27 @@ _detect_default_field(::Type{QQFieldElem}, p::Polymake.BigObject) = QQ _detect_default_field(::Type{Float64}, p::Polymake.BigObject) = AbstractAlgebra.Floats{Float64}() function _detect_default_field(::Type{T}, p::Polymake.BigObject) where T<:FieldElem - # we only want to check existing properties - f = x -> Polymake.exists(p, string(x)) - propnames = intersect(propertynames(p), [:INPUT_RAYS, :POINTS, :RAYS, :VERTICES, :VECTORS, :INPUT_LINEALITY, :LINEALITY_SPACE, :FACETS, :INEQUALITIES, :EQUATIONS, :LINEAR_SPAN, :AFFINE_HULL]) - i = findfirst(f, propnames) - # find first OscarNumber wrapping a FieldElem - while !isnothing(i) - prop = getproperty(p, propnames[i]) - for el in prop - on = Polymake.unwrap(el) - if on isa T - return parent(on) - end - end - i = findnext(f, propnames, i + 1) + # we only want to check existing properties + propnames = intersect(Polymake.list_properties(p), ["INPUT_RAYS", "POINTS", "RAYS", "VERTICES", "VECTORS", "INPUT_LINEALITY", "LINEALITY_SPACE", "FACETS", "INEQUALITIES", "EQUATIONS", "LINEAR_SPAN", "AFFINE_HULL"]) + # find first OscarNumber wrapping a FieldElem + for pn in propnames + prop = getproperty(p, pn) + for el in prop + on = Polymake.unwrap(el) + if on isa T + return parent(on) + end end - throw(ArgumentError("BigObject does not contain information about a parent Field")) + end + throw(ArgumentError("BigObject does not contain information about a parent Field")) end function _detect_wrapped_type_and_field(p::Polymake.BigObject) # we only want to check existing properties - f = x -> Polymake.exists(p, string(x)) - propnames = intersect(propertynames(p), [:INPUT_RAYS, :POINTS, :RAYS, :VERTICES, :VECTORS, :INPUT_LINEALITY, :LINEALITY_SPACE, :FACETS, :INEQUALITIES, :EQUATIONS, :LINEAR_SPAN, :AFFINE_HULL]) - i = findfirst(f, propnames) + propnames = intersect(Polymake.list_properties(p), ["INPUT_RAYS", "POINTS", "RAYS", "VERTICES", "VECTORS", "INPUT_LINEALITY", "LINEALITY_SPACE", "FACETS", "INEQUALITIES", "EQUATIONS", "LINEAR_SPAN", "AFFINE_HULL"]) # find first OscarNumber wrapping a FieldElem - while !isnothing(i) - prop = getproperty(p, propnames[i]) + for pn in propnames + prop = getproperty(p, pn) for el in prop on = Polymake.unwrap(el) if on isa FieldElem @@ -514,7 +510,6 @@ function _detect_wrapped_type_and_field(p::Polymake.BigObject) return (T, f) end end - i = findnext(f, propnames, i + 1) end throw(ArgumentError("BigObject does not contain information about a parent Field")) end diff --git a/src/PolyhedralGeometry/iterators.jl b/src/PolyhedralGeometry/iterators.jl index af03742d2ebc..fe0f1154b60c 100644 --- a/src/PolyhedralGeometry/iterators.jl +++ b/src/PolyhedralGeometry/iterators.jl @@ -352,8 +352,7 @@ for f in ("_point_matrix", "_vector_matrix", "_generator_matrix") M = Symbol(f) @eval begin function $M(::Val{_empty_access}, P::PolyhedralObjectUnion; homogenized=false) - scalar_regexp = match(r"[^<]*<(.*)>[^>]*", String(Polymake.type_name(pm_object(P)))) - typename = scalar_regexp[1] + typename = Polymake.bigobject_eltype(pm_object(P)) T = typename == "OscarNumber" ? Polymake.OscarNumber : _scalar_type_to_polymake(scalar_type_to_oscar[typename]) return Polymake.Matrix{T}(undef, 0, Polymake.polytope.ambient_dim(pm_object(P)) + homogenized) end diff --git a/src/PolyhedralGeometry/type_functions.jl b/src/PolyhedralGeometry/type_functions.jl index 5ed4480fa391..73bfb92b6375 100644 --- a/src/PolyhedralGeometry/type_functions.jl +++ b/src/PolyhedralGeometry/type_functions.jl @@ -4,8 +4,7 @@ function detect_scalar_type(n::Type{T}, p::Polymake.BigObject) where T<:Union{Polyhedron, Cone, PolyhedralFan, SubdivisionOfPoints, PolyhedralComplex} - scalar_regexp = match(r"[^<]*<(.*)>[^>]*", String(Polymake.type_name(p))) - typename = scalar_regexp[1] + typename = Polymake.bigobject_eltype(p) return typename == "OscarNumber" ? nothing : scalar_type_to_oscar[typename] end diff --git a/test/PolyhedralGeometry/polyhedron.jl b/test/PolyhedralGeometry/polyhedron.jl index 33b63f15f2bd..05d189d3e60a 100644 --- a/test/PolyhedralGeometry/polyhedron.jl +++ b/test/PolyhedralGeometry/polyhedron.jl @@ -532,7 +532,7 @@ for f in (QQ, ENF) if S == Pair{Matrix{T}, T} @test facets(S, D) == [Pair(A[i], b[i]) for i in 1:12] elseif S == Polyhedron{T} - # @test facets(S, D) == [polyhedron(R, A[i], b[i]) for i in 1:12] + @test nvertices.(facets(S, D)) == repeat([5], 12) else @test facets(S, D) == [affine_halfspace(R, A[i], b[i]) for i in 1:12] end @@ -565,7 +565,7 @@ for f in (QQ, ENF) @test isempty(rays(D)) @test lineality_dim(D) == 0 @test isempty(lineality_space(D)) - # @test faces(D, 0) == convex_hull.(T, V) + @test faces(D, 0) == convex_hull.(T, V) @test isempty(affine_hull(D)) @test relative_interior_point(D) == [0, 0, 0] end diff --git a/test/PolyhedralGeometry/types.jl b/test/PolyhedralGeometry/types.jl index 924075566a6f..3256888f9dfb 100644 --- a/test/PolyhedralGeometry/types.jl +++ b/test/PolyhedralGeometry/types.jl @@ -65,8 +65,6 @@ @test *(g(3), A) == 3 * a - # @test [A; B] == [a; b] - end end From 24d8d716c704615b87884ed0a79cbdf442ab7a48 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Sun, 17 Sep 2023 17:54:37 +0200 Subject: [PATCH 20/21] workaround for accessing property using String --- src/PolyhedralGeometry/helpers.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/PolyhedralGeometry/helpers.jl b/src/PolyhedralGeometry/helpers.jl index ab406a4a859d..9b6ec0955080 100644 --- a/src/PolyhedralGeometry/helpers.jl +++ b/src/PolyhedralGeometry/helpers.jl @@ -485,7 +485,7 @@ function _detect_default_field(::Type{T}, p::Polymake.BigObject) where T<:FieldE propnames = intersect(Polymake.list_properties(p), ["INPUT_RAYS", "POINTS", "RAYS", "VERTICES", "VECTORS", "INPUT_LINEALITY", "LINEALITY_SPACE", "FACETS", "INEQUALITIES", "EQUATIONS", "LINEAR_SPAN", "AFFINE_HULL"]) # find first OscarNumber wrapping a FieldElem for pn in propnames - prop = getproperty(p, pn) + prop = getproperty(p, convert(String, pn)) for el in prop on = Polymake.unwrap(el) if on isa T @@ -501,7 +501,7 @@ function _detect_wrapped_type_and_field(p::Polymake.BigObject) propnames = intersect(Polymake.list_properties(p), ["INPUT_RAYS", "POINTS", "RAYS", "VERTICES", "VECTORS", "INPUT_LINEALITY", "LINEALITY_SPACE", "FACETS", "INEQUALITIES", "EQUATIONS", "LINEAR_SPAN", "AFFINE_HULL"]) # find first OscarNumber wrapping a FieldElem for pn in propnames - prop = getproperty(p, pn) + prop = getproperty(p, convert(String, pn)) for el in prop on = Polymake.unwrap(el) if on isa FieldElem From ebf08ab9314e1c6bc39563b215c059bfca812761 Mon Sep 17 00:00:00 2001 From: Alexej Jordan Date: Mon, 18 Sep 2023 21:45:02 +0200 Subject: [PATCH 21/21] adjusted doctest print to new PointVector multiplication return type --- src/PolyhedralGeometry/solving_integrally.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/PolyhedralGeometry/solving_integrally.jl b/src/PolyhedralGeometry/solving_integrally.jl index 7a497471b753..85c0ec4d2f56 100644 --- a/src/PolyhedralGeometry/solving_integrally.jl +++ b/src/PolyhedralGeometry/solving_integrally.jl @@ -57,7 +57,7 @@ SubObjectIterator{PointVector{ZZRingElem}} julia> for x in it print(A*x," ") end -ZZRingElem[7] ZZRingElem[7] ZZRingElem[7] ZZRingElem[7] ZZRingElem[7] ZZRingElem[7] ZZRingElem[7] ZZRingElem[7] +[7] [7] [7] [7] [7] [7] [7] [7] ``` """ solve_mixed(as::Type{T}, A::ZZMatrix, b::ZZMatrix, C::ZZMatrix, d::ZZMatrix) where {T} = solve_mixed(T, A, b, C, d) @@ -106,7 +106,7 @@ SubObjectIterator{PointVector{ZZRingElem}} julia> for x in it print(A*x," ") end -ZZRingElem[3] ZZRingElem[3] ZZRingElem[3] ZZRingElem[3] +[3] [3] [3] [3] ``` """ solve_mixed(as::Type{T}, A::ZZMatrix, b::ZZMatrix, C::ZZMatrix) where {T} = solve_mixed(T, A, b, C, zero_matrix(FlintZZ, nrows(C), 1))