Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

#193 #217

Merged
merged 6 commits into from
Feb 6, 2018
Merged

#193 #217

Show file tree
Hide file tree
Changes from 4 commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
1 change: 1 addition & 0 deletions src/LazySets.jl
Original file line number Diff line number Diff line change
Expand Up @@ -63,6 +63,7 @@ include("HPolygon.jl")
include("HPolygonOpt.jl")
include("HPolytope.jl")
include("VPolygon.jl")
include("VPolytope.jl")
include("Zonotope.jl")
include("Ellipsoid.jl")
include("Hyperplane.jl")
Expand Down
224 changes: 224 additions & 0 deletions src/VPolytope.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,224 @@
using MathProgBase, GLPKMathProgInterface

export VPolytope,
vertices_list

"""
VPolytope{N<:Real} <: AbstractPolytope{N}

Type that represents a convex polytope in V-representation.

### Fields

- `vertices` -- the list of vertices
"""
struct VPolytope{N<:Real} <: AbstractPolytope{N}
vertices::Vector{Vector{N}}
end
# constructor for a VPolytope with empty vertices list
VPolytope{N}() where {N<:Real} = VPolytope{N}(Vector{N}(0))

# constructor for a VPolytope with no constraints of type Float64
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

"constraints" → "vertices"

VPolytope() = VPolytope{Float64}()

# constructor from a polygon in V-representation
VPolytope(P::VPolygon) = VPolytope(P.vertices_list)


# --- LazySet interface functions ---


"""
dim(P::VPolytope)::Int

Return the dimension of a polytope in V-representation.

### Input

- `P` -- polytope in V-representation

### Output

The ambient dimension of the polytope in V-representation.
If it is empty, the result is ``-1``.

### Examples

```jldoctest
julia> v = VPolytope()
julia> dim(v) > 0
false
julia> v = VPolytope([ones(3)])
LazySets.VPolytope{Float64}(Array{Float64,1}[[1.0, 1.0, 1.0]])
julia> dim(v) == 3
true
```
"""
function dim(P::VPolytope)::Int
return length(P.vertices) == 0 ? -1 : length(P.vertices[1])
end

"""
σ(d::AbstractVector{<:Real}, P::VPolytope; algorithm="hrep")::Vector{<:Real}

Return the support vector of a polyhedron (in V-representation) in a given
direction.

### Input

- `d` -- direction
- `P` -- polyhedron in V-representation
- `algorithm` -- (optional, default: `'hrep'`) method to compute the support vector

### Output

The support vector in the given direction.
"""
function σ(d::AbstractVector{<:Real}, P::VPolytope; algorithm="hrep")::Vector{<:Real}
if algorithm == "hrep"
@assert isdefined(:Polyhedra) "you need to load Polyhedra for this algorithm"
return σ(d, tohrep(P))
else
error("the algorithm $(hrep) is not known")
end
end

# --- VPolytope functions ---

"""
vertices_list(P::VPolytope{N})::Vector{Vector{N}} where {N<:Real}

Return the list of vertices of a polytope in V-representation.

### Input

- `P` -- polytope in vertex representation

### Output

List of vertices.
"""
function vertices_list(P::VPolytope{N})::Vector{Vector{N}} where {N<:Real}
return P.vertices
end


@require Polyhedra begin

using CDDLib # default backend
import Polyhedra:polyhedron, SimpleHRepresentation, SimpleVRepresentation,
HRep, VRep,
removehredundancy!, removevredundancy!,
hreps, vreps,
intersect,
convexhull,
hcartesianproduct

export intersect, convex_hull, cartesian_product, vertices_list

# VPolytope from a VRep
function VPolytope(P::VRep{N, T}, backend=CDDLib.CDDLibrary()) where {N, T}
vertices = Vector{Vector{T}}()
for vi in vreps(P)
push!(vertices, vi)
end
return VPolytope(vertices)
end

"""
polyhedron(P::VPolytope{N}, [backend]=CDDLib.CDDLibrary()) where {N}

Return an `VRep` polyhedron from `Polyhedra.jl` given a polytope in V-representation.

### Input

- `P` -- polytope
- `backend` -- (optional, default: `CDDLib.CDDLibrary()`) the polyhedral
computations backend, see [Polyhedra's documentation](https://juliapolyhedra.github.io/Polyhedra.jl/latest/installation.html#Getting-Libraries-1)
for further information

### Output

A `VRep` polyhedron.
"""
function polyhedron(P::VPolytope{N}, backend=CDDLib.CDDLibrary()) where {N}
return polyhedron(SimpleVRepresentation(hcat(vertices_list(P)...)'), backend)
end

"""
intersect(P1::VPolytope{N}, P2::VPolytope{N};
[backend]=CDDLib.CDDLibrary(),
[prunefunc]=removehredundancy!)::VPolytope{N} where {N<:Real}

Compute the intersection of two polytopes in V-representation.

### Input

- `P1` -- polytope
- `P2` -- another polytope
- `backend` -- (optional, default: `CDDLib.CDDLibrary()`) the polyhedral
computations backend, see [Polyhedra's documentation](https://juliapolyhedra.github.io/Polyhedra.jl/latest/installation.html#Getting-Libraries-1)
for further information
- `prunefunc` -- (optional, default: `removehredundancy!`) function to post-process
the output of `intersect`

### Output

The `VPolytope` obtained by the intersection of `P1` and `P2`.
"""
function intersect(P1::VPolytope{N}, P2::VPolytope{N};
backend=CDDLib.CDDLibrary(),
prunefunc=removehredundancy!)::VPolytope{N} where {N<:Real}

P1 = polyhedron(P1, backend)
P2 = polyhedron(P2, backend)
Pint = intersect(P1, P2)
prunefunc(Pint)
return VPolytope(Pint)
end

"""
convex_hull(P1::VPolytope, P2::VPolytope; [backend]=CDDLib.CDDLibrary())

Compute the convex hull of the set union of two polytopes in V-representation.

### Input

- `P1` -- polytope
- `P2` -- another polytope
- `backend` -- (optional, default: `CDDLib.CDDLibrary()`) the polyhedral
computations backend, see [Polyhedra's documentation](https://juliapolyhedra.github.io/Polyhedra.jl/latest/installation.html#Getting-Libraries-1)
for further information

### Output

The `VPolytope` obtained by the concrete convex hull of `P1` and `P2`.
"""
function convex_hull(P1::VPolytope, P2::VPolytope; backend=CDDLib.CDDLibrary())
Pch = convexhull(polyhedron(P1, backend), polyhedron(P2, backend))
return VPolytope(Pch)
end

"""
cartesian_product(P1::VPolytope, P2::VPolytope; [backend]=CDDLib.CDDLibrary())

Compute the Cartesian product of two polytopes in V-representaion.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

"representation"


### Input

- `P1` -- polytope
- `P2` -- another polytope
- `backend` -- (optional, default: `CDDLib.CDDLibrary()`) the polyhedral
computations backend, see [Polyhedra's documentation](https://juliapolyhedra.github.io/Polyhedra.jl/latest/installation.html#Getting-Libraries-1)
for further information

### Output

The `VPolytope` obtained by the concrete cartesian product of `P1` and `P2`.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

"Cartesian"

"""
function cartesian_product(P1::VPolytope, P2::VPolytope; backend=CDDLib.CDDLibrary())
Pcp = hcartesianproduct(polyhedron(P1, backend), polyhedron(P2, backend))
return VPolytope(Pcp)
end

end