From f58232064087c6376fbc359d0a0ffbe3802f7374 Mon Sep 17 00:00:00 2001 From: Lorenzo Van Munoz <66997677+lxvm@users.noreply.github.com> Date: Thu, 1 Feb 2024 08:15:37 -0500 Subject: [PATCH] [BREAKING] Move functionality to extension packages (#328) * move jump to extension * support for regular dependency in julia 1.6 * reorganize JuMP ext * move RecipesBase and GeometryBasics to ext * uncomment verbose statement * add error message for Mesh * add news file describing breaking change * update docs * update examples * simplify doc * fix verbose print of model * remove exported JuMP function * import functions missing in GeometryBasicsExt * add missing compat entry --- NEWS.md | 14 +++++ Project.toml | 16 +++++- docs/src/plot.md | 21 ++++---- ...a projection of the 4D permutahedron.ipynb | 19 +------ examples/Polyhedral Function.ipynb | 10 +--- .../PolyhedraGeometryBasicsExt.jl | 10 +++- ext/PolyhedraJuMPExt.jl | 53 +++++++++++++++++++ .../PolyhedraRecipesBaseExt.jl | 9 +++- src/Polyhedra.jl | 27 +++++++--- src/center.jl | 1 - src/jump.jl | 11 ---- src/linearity.jl | 3 +- src/opt.jl | 32 +---------- src/projection_opt.jl | 5 -- test/decompose.jl | 3 +- 15 files changed, 137 insertions(+), 97 deletions(-) create mode 100644 NEWS.md rename src/decompose.jl => ext/PolyhedraGeometryBasicsExt.jl (97%) create mode 100644 ext/PolyhedraJuMPExt.jl rename src/recipe.jl => ext/PolyhedraRecipesBaseExt.jl (86%) delete mode 100644 src/jump.jl diff --git a/NEWS.md b/NEWS.md new file mode 100644 index 00000000..e4b544ee --- /dev/null +++ b/NEWS.md @@ -0,0 +1,14 @@ +# Polyhedra.jl v0.8 Release Notes + +## Load time improvements +- Dependencies on JuMP.jl, RecipesBase.jl, and GeometryBasics.jl were moved to + weak dependencies on Julia versions supporting package extensions, i.e. v1.9 + and above. On v1.10 this reduces installation time by 15% and load time by + 18% (see [#328]). + +## Breaking changes +- `JuMP.optimizer_with_attributes` is no longer exported. Call it from JuMP.jl instead. +- The following change is only breaking on Julia v1.9 and above: + `Polyhedra.Mesh` is now implemented in a package extension requiring + GeometryBasics.jl. It is sufficient to load your plotting package, i.e. + Makie.jl or MeshCat.jl, **before** calling `Polyhedra.Mesh` \ No newline at end of file diff --git a/Project.toml b/Project.toml index 99b32e04..c3ffbbb1 100644 --- a/Project.toml +++ b/Project.toml @@ -6,18 +6,32 @@ version = "0.7.7" [deps] GenericLinearAlgebra = "14197337-ba66-59df-a3e3-ca00e7dcff7a" GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326" -JuMP = "4076af6c-e467-56ae-b986-b466b2749572" LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +JuMP = "4076af6c-e467-56ae-b986-b466b2749572" +MathOptInterface = "b8f27783-ece8-5eb3-8dc8-9495eed66fee" MutableArithmetics = "d8a4904e-b15c-11e9-3269-09a3773c0cb0" RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" +[weakdeps] +GeometryBasics = "5c1252a2-5f33-56bf-86c9-59e7332b4326" +JuMP = "4076af6c-e467-56ae-b986-b466b2749572" +RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01" + +[extensions] +PolyhedraGeometryBasicsExt = "GeometryBasics" +PolyhedraJuMPExt = "JuMP" +PolyhedraRecipesBaseExt = "RecipesBase" + [compat] GenericLinearAlgebra = "0.2, 0.3" GeometryBasics = "0.2, 0.3, 0.4" JuMP = "0.23, 1" +LinearAlgebra = "1.6" +MathOptInterface = "1" MutableArithmetics = "1" RecipesBase = "0.7, 0.8, 1.0" +SparseArrays = "1.6" StaticArrays = "0.12, 1.0" julia = "1.6" diff --git a/docs/src/plot.md b/docs/src/plot.md index a3903794..4770dbe1 100644 --- a/docs/src/plot.md +++ b/docs/src/plot.md @@ -63,26 +63,25 @@ Polyhedron DefaultPolyhedron{Rational{BigInt}, Polyhedra.Intersection{Rational{B Ray(Rational{BigInt}[0, 0, 1]) ``` -Then, we need to create a mesh from the polyhedron: -```jldoctest plots3 -julia> m = Polyhedra.Mesh(p) -Polyhedra.Mesh{3, Rational{BigInt}, DefaultPolyhedron{Rational{BigInt}, Polyhedra.Intersection{Rational{BigInt}, Vector{Rational{BigInt}}, Int64}, Polyhedra.Hull{Rational{BigInt}, Vector{Rational{BigInt}}, Int64}}}(convexhull([0, 0, 0]) + convexhull(Ray(Rational{BigInt}[1, 0, 0]), Ray(Rational{BigInt}[0, 1, 0]), Ray(Rational{BigInt}[0, 0, 1])), nothing, nothing, nothing) -``` - -```@docs -Polyhedra.Mesh -``` - -The polyhedron can be plotted with [MeshCat](https://github.com/rdeits/MeshCat.jl) as follows +The polyhedron can then be plotted with [MeshCat](https://github.com/rdeits/MeshCat.jl) as follows ```julia julia> using MeshCat +julia> m = Polyhedra.Mesh(p) + julia> vis = Visualizer() julia> setobject!(vis, m) julia> open(vis) ``` + +Note that the `Mesh` object should be created **after** loading the plotting package: + +```@docs +Polyhedra.Mesh +``` + To plot it in a notebook, replace `open(vis)` with `IJuliaCell(vis)`. To plot it with [Makie](https://github.com/JuliaPlots/Makie.jl) instead, you can use for instance `mesh` or `wireframe`. diff --git a/examples/3D Plotting a projection of the 4D permutahedron.ipynb b/examples/3D Plotting a projection of the 4D permutahedron.ipynb index be2d91e3..f8c18ea0 100644 --- a/examples/3D Plotting a projection of the 4D permutahedron.ipynb +++ b/examples/3D Plotting a projection of the 4D permutahedron.ipynb @@ -77,23 +77,7 @@ "cell_type": "markdown", "metadata": {}, "source": [ - "To get a plottable object, we transform the polyhedron into a mesh as follows." - ] - }, - { - "cell_type": "code", - "execution_count": 5, - "metadata": {}, - "outputs": [], - "source": [ - "m = Polyhedra.Mesh(p3);" - ] - }, - { - "cell_type": "markdown", - "metadata": {}, - "source": [ - "We can now plot this mesh with [MeshCat](https://github.com/rdeits/MeshCat.jl) or [Makie](https://github.com/JuliaPlots/Makie.jl) as follows:" + "We can now plot this polyhedron with [MeshCat](https://github.com/rdeits/MeshCat.jl) or [Makie](https://github.com/JuliaPlots/Makie.jl) as follows:" ] }, { @@ -129,6 +113,7 @@ ], "source": [ "using MeshCat\n", + "m = Polyhedra.Mesh(p3)\n", "vis = Visualizer()\n", "setobject!(vis, m)\n", "IJuliaCell(vis)" diff --git a/examples/Polyhedral Function.ipynb b/examples/Polyhedral Function.ipynb index 49358a5a..4f951718 100644 --- a/examples/Polyhedral Function.ipynb +++ b/examples/Polyhedral Function.ipynb @@ -840,15 +840,6 @@ "The top of the shape will have no face as it is unbounded in this direction." ] }, - { - "cell_type": "code", - "execution_count": null, - "metadata": {}, - "outputs": [], - "source": [ - "m = Polyhedra.Mesh(p)" - ] - }, { "cell_type": "code", "execution_count": null, @@ -856,6 +847,7 @@ "outputs": [], "source": [ "using MeshCat\n", + "m = Polyhedra.Mesh(p)\n", "vis = Visualizer()\n", "setobject!(vis, m)\n", "IJuliaCell(vis)" diff --git a/src/decompose.jl b/ext/PolyhedraGeometryBasicsExt.jl similarity index 97% rename from src/decompose.jl rename to ext/PolyhedraGeometryBasicsExt.jl index 949f72a3..7c4f0319 100644 --- a/src/decompose.jl +++ b/ext/PolyhedraGeometryBasicsExt.jl @@ -1,4 +1,10 @@ +module PolyhedraGeometryBasicsExt + +using LinearAlgebra import GeometryBasics +using Polyhedra +using Polyhedra: FullDim, isapproxzero, _planar_hull, counterclockwise, rotate +using StaticArrays """ struct Mesh{N, T, PT <: Polyhedron{T}} <: GeometryBasics.GeometryPrimitive{N, T} @@ -29,7 +35,7 @@ function Mesh(polyhedron::Polyhedron, N::Int) # use polyhedron built from StaticArrays vector to avoid that. return Mesh{N}(polyhedron) end -function Mesh(polyhedron::Polyhedron) +function Polyhedra.Mesh(polyhedron::Polyhedron) return Mesh(polyhedron, FullDim(polyhedron)) end @@ -223,3 +229,5 @@ GeometryBasics.coordinates(poly::Mesh) = (fulldecompose!(poly); poly.coordinates GeometryBasics.faces(poly::Mesh) = (fulldecompose!(poly); poly.faces) GeometryBasics.texturecoordinates(poly::Mesh) = nothing GeometryBasics.normals(poly::Mesh) = (fulldecompose!(poly); poly.normals) + +end diff --git a/ext/PolyhedraJuMPExt.jl b/ext/PolyhedraJuMPExt.jl new file mode 100644 index 00000000..dfd92467 --- /dev/null +++ b/ext/PolyhedraJuMPExt.jl @@ -0,0 +1,53 @@ +module PolyhedraJuMPExt + +import JuMP +import Polyhedra: hrep, LPHRep, polyhedron, _optimize! +using Polyhedra: Rep, Projection, _moi_set, fulldim, dimension_names, PolyhedraToLPBridge, ProjectionBridge + +""" + hrep(model::JuMP.Model) + +Builds an H-representation from the feasibility set of the JuMP model `model`. +Note that if non-linear constraint are present in the model, they are ignored. +""" +hrep(model::JuMP.Model) = LPHRep(model) +LPHRep(model::JuMP.Model) = LPHRep(JuMP.backend(model)) +polyhedron(model::JuMP.Model, args...) = polyhedron(hrep(model), args...) +_optimize!(model::JuMP.Model) = JuMP.optimize!(model) + +struct VariableInSet{V <: JuMP.ScalarVariable, S <: Union{Rep, Projection}} <: JuMP.AbstractVariable + variables::Vector{V} + set::S +end +function JuMP.build_variable(error_fun::Function, variables::Vector{<:JuMP.ScalarVariable}, set::Union{Rep, Projection}) + if length(variables) != fulldim(set) + error("Number of variables ($(length(variables))) does not match the full dimension of the polyhedron ($(fulldim(set))).") + end + return VariableInSet(variables, set) +end +function JuMP.add_variable(model::JuMP.AbstractModel, v::VariableInSet, names) + dim_names = dimension_names(v.set) + if dim_names !== nothing + names = copy(names) + for i in eachindex(names) + if isempty(names[i]) && !isempty(dim_names[i]) + names[i] = dim_names[i] + end + end + end + JuMP.add_bridge(model, PolyhedraToLPBridge) + JuMP.add_bridge(model, ProjectionBridge) + return JuMP.add_variable(model, JuMP.VariablesConstrainedOnCreation(v.variables, _moi_set(v.set)), names) +end +function JuMP.build_constraint(error_fun::Function, func, set::Rep) + return JuMP.BridgeableConstraint( + JuMP.build_constraint(error_fun, func, _moi_set(set)), + PolyhedraToLPBridge) +end +function JuMP.build_constraint(error_fun::Function, func, set::Projection) + return JuMP.BridgeableConstraint(JuMP.BridgeableConstraint( + JuMP.build_constraint(error_fun, func, _moi_set(set)), + ProjectionBridge), PolyhedraToLPBridge) +end + +end diff --git a/src/recipe.jl b/ext/PolyhedraRecipesBaseExt.jl similarity index 86% rename from src/recipe.jl rename to ext/PolyhedraRecipesBaseExt.jl index 2919a043..dd341923 100644 --- a/src/recipe.jl +++ b/ext/PolyhedraRecipesBaseExt.jl @@ -1,4 +1,9 @@ +module PolyhedraRecipesBaseExt + +using LinearAlgebra using RecipesBase +using Polyhedra +using Polyhedra: basis, _semi_hull function planar_contour(p::Polyhedron) if fulldim(p) != 2 @@ -17,7 +22,7 @@ function planar_contour(p::Polyhedron) error("Plotting empty polyhedron is not supported.") end sort!(ps, by = x -> x[1]) - counterclockwise(p1, p2) = dot(cross([p1; 0], [p2; 0]), [0, 0, 1]) + counterclockwise = (p1, p2) -> dot(cross([p1; 0], [p2; 0]), [0, 0, 1]) sweep_norm = basis(eltype(ps), fulldim(p), 1) top = _semi_hull(ps, 1, counterclockwise, sweep_norm) bot = _semi_hull(ps, -1, counterclockwise, sweep_norm) @@ -39,3 +44,5 @@ end legend --> false planar_contour(p) end + +end diff --git a/src/Polyhedra.jl b/src/Polyhedra.jl index 5f1c0ff7..bec47008 100644 --- a/src/Polyhedra.jl +++ b/src/Polyhedra.jl @@ -12,14 +12,14 @@ import GenericLinearAlgebra import MutableArithmetics const MA = MutableArithmetics +import MathOptInterface as MOI +import MathOptInterface.Utilities as MOIU + export Polyhedron abstract type Library end abstract type Polyhedron{T} end -using JuMP -export optimizer_with_attributes - coefficient_type(::Union{AbstractVector{T}, Type{<:AbstractVector{T}}}) where T = T similar_type(::Type{<:Vector}, ::Int, ::Type{T}) where T = Vector{T} similar_type(::Type{SparseVector{S, IT}}, ::Int, ::Type{T}) where {S, IT, T} = SparseVector{T, IT} @@ -81,7 +81,6 @@ include("extended.jl") include("vecrep.jl") include("mixedrep.jl") include("lphrep.jl") -include("jump.jl") include("matrep.jl") include("liftedrep.jl") include("doubledescription.jl") # FIXME move it after projection.jl once it stops depending on LiftedRep @@ -100,7 +99,23 @@ include("projection_opt.jl") # Visualization include("show.jl") -include("recipe.jl") -include("decompose.jl") + +""" + Mesh(p::Polyhedron) + +Returns wrapper of a polyhedron suitable for plotting with MeshCat.jl and Makie.jl. + +!!! note "Extension in Julia 1.9 and above" + Although we require `using GeometryBasics` to use this function in Julia 1.9 and above, + in most use cases this extension dependency is loaded by the plotting package and no + further action is required. +""" +Mesh(p) = p isa Polyhedron ? error("this method requires using GeometryBasics") : throw(MethodError(Mesh, p)) + +if !isdefined(Base, :get_extension) + include("../ext/PolyhedraJuMPExt.jl") + include("../ext/PolyhedraRecipesBaseExt.jl") + include("../ext/PolyhedraGeometryBasicsExt.jl") +end end # module diff --git a/src/center.jl b/src/center.jl index 7a23254e..3ac6b476 100644 --- a/src/center.jl +++ b/src/center.jl @@ -1,5 +1,4 @@ export maximum_radius_with_center, hchebyshevcenter, vchebyshevcenter, chebyshevcenter -using JuMP """ maximum_radius_with_center(h::HRep, center) diff --git a/src/jump.jl b/src/jump.jl deleted file mode 100644 index c93538a8..00000000 --- a/src/jump.jl +++ /dev/null @@ -1,11 +0,0 @@ -using JuMP - -""" - hrep(model::JuMP.Model) - -Builds an H-representation from the feasibility set of the JuMP model `model`. -Note that if non-linear constraint are present in the model, they are ignored. -""" -hrep(model::JuMP.Model) = LPHRep(model) -LPHRep(model::JuMP.Model) = LPHRep(backend(model)) -polyhedron(model::JuMP.Model, args...) = polyhedron(hrep(model), args...) diff --git a/src/linearity.jl b/src/linearity.jl index 607a27fc..3fcab565 100644 --- a/src/linearity.jl +++ b/src/linearity.jl @@ -151,7 +151,6 @@ function detect_new_linearities(rep::Representation, solver; verbose=0, kws...) nonlins = _nonlinearity(rep) isempty(nonlins) && return Int[] model, T = layered_optimizer(solver) - verbose >= 2 && (_model = JuMP.direct_model(model)) is_lin = falses(length(nonlins)) active = trues(length(nonlins)) # We pass `true` as we break homogeneity of the problem with `sum(λ) == 1`. @@ -173,7 +172,7 @@ function detect_new_linearities(rep::Representation, solver; verbose=0, kws...) end for i in 1:fulldim(rep) # safer than `while ...` (count(is_lin) == length(is_lin) || iszero(count(active))) && break - verbose >= 2 && println(_model) + verbose >= 2 && println(model) # FIXME stopping when we have enough hyperplanes to prove that it's empty # should also resolve the issue with presolve. is_feasible(model, "detecting new linearity (you may need to activate presolve for some solvers, e.g. by replacing `GLPK.Optimizer` with `optimizer_with_attributes(GLPK.Optimizer, \"presolve\" => GLPK.GLP_ON)`).") || break diff --git a/src/opt.jl b/src/opt.jl index 2e78a427..623029a3 100644 --- a/src/opt.jl +++ b/src/opt.jl @@ -9,36 +9,6 @@ MOI.dimension(set::PolyhedraOptSet) = fulldim(set.rep) # will be modified so we don't copy it as it would be expensive. Base.copy(set::PolyhedraOptSet) = set -struct VariableInSet{V <: JuMP.ScalarVariable, S <: Union{Rep, Projection}} <: JuMP.AbstractVariable - variables::Vector{V} - set::S -end -function JuMP.build_variable(error_fun::Function, variables::Vector{<:JuMP.ScalarVariable}, set::Union{Rep, Projection}) - if length(variables) != fulldim(set) - _error("Number of variables ($(length(variables))) does not match the full dimension of the polyhedron ($(fulldim(set))).") - end - return VariableInSet(variables, set) -end -function JuMP.add_variable(model::JuMP.AbstractModel, v::VariableInSet, names) - dim_names = dimension_names(v.set) - if dim_names !== nothing - names = copy(names) - for i in eachindex(names) - if isempty(names[i]) && !isempty(dim_names[i]) - names[i] = dim_names[i] - end - end - end - JuMP.add_bridge(model, PolyhedraToLPBridge) - JuMP.add_bridge(model, ProjectionBridge) - return JuMP.add_variable(model, JuMP.VariablesConstrainedOnCreation(v.variables, _moi_set(v.set)), names) -end -_moi_set(set::Rep) = PolyhedraOptSet(set) -function JuMP.build_constraint(error_fun::Function, func, set::Rep) - return JuMP.BridgeableConstraint( - JuMP.build_constraint(error_fun, func, _moi_set(set)), - PolyhedraToLPBridge) -end abstract type AbstractPolyhedraOptimizer{T} <: MOI.AbstractOptimizer end @@ -153,7 +123,6 @@ function layered_optimizer(solver) return optimizer, T end -_optimize!(model::JuMP.Model) = JuMP.optimize!(model) _optimize!(model::MOI.ModelLike) = MOI.optimize!(model) function _unknown_status(model, status, message) @@ -195,3 +164,4 @@ function Base.isempty(p::Rep, solver=Polyhedra.linear_objective_solver(p)) x, cx = MOI.add_constrained_variables(model, PolyhedraOptSet(p)) return !is_feasible(model, "trying to determine whether the polyhedron is empty.") end +_moi_set(set::Rep) = PolyhedraOptSet(set) diff --git a/src/projection_opt.jl b/src/projection_opt.jl index 240a12b5..d1f819a8 100644 --- a/src/projection_opt.jl +++ b/src/projection_opt.jl @@ -58,8 +58,3 @@ function MOI.delete(model::MOI.ModelLike, b::ProjectionBridge) end end _moi_set(set::Projection) = ProjectionOptSet(set) -function JuMP.build_constraint(error_fun::Function, func, set::Projection) - return JuMP.BridgeableConstraint(JuMP.BridgeableConstraint( - JuMP.build_constraint(error_fun, func, _moi_set(set)), - ProjectionBridge), PolyhedraToLPBridge) -end diff --git a/test/decompose.jl b/test/decompose.jl index a5899d1c..ce8de6f6 100644 --- a/test/decompose.jl +++ b/test/decompose.jl @@ -23,7 +23,8 @@ end nfaces(d::Dict{<:Any, Face}) = length(d) nfaces(d::Dict{<:Any, <:Vector}) = sum(map(length, values(d))) -function test_decompose(p::Polyhedra.Mesh{N}, d::Dict) where N +function test_decompose(p, d::Dict) + N = ndims(p) P = GeometryBasics.Point{N, Float64} NR = GeometryBasics.Point{3, Float64} points = GeometryBasics.decompose(P, p)