Drake
Drake C++ Documentation

Detailed Description

A polytope described using the vertex representation.

The set is defined as the convex hull of the vertices. The vertices are not guaranteed to be in any particular order, nor to be minimal (some vertices could be strictly in the interior of the set).

Note: Unlike the half-space representation, this definition means the set is always bounded (hence the name polytope, instead of polyhedron).

A VPolytope is empty if and only if it is composed of zero vertices, i.e., if vertices_.cols() == 0. This includes the zero-dimensional case. If vertices_.rows() == 0 but vertices_.cols() > 0, we treat this as having one or more copies of 0 in the zero-dimensional vector space {0}. If vertices_.rows() and vertices_.cols() are zero, we treat this as no points in {0}, which is empty.

#include <drake/geometry/optimization/vpolytope.h>

Public Member Functions

 VPolytope ()
 Constructs a set with no vertices in the zero-dimensional space, which is empty (by convention). More...
 
 VPolytope (const Eigen::Ref< const Eigen::MatrixXd > &vertices)
 Constructs the polytope from a d-by-n matrix, where d is the ambient dimension, and n is the number of vertices. More...
 
 VPolytope (const HPolyhedron &H, const double tol=1e-9)
 Constructs the polytope from a bounded polyhedron (using Qhull). More...
 
 VPolytope (const QueryObject< double > &query_object, GeometryId geometry_id, std::optional< FrameId > reference_frame=std::nullopt)
 Constructs the polytope from a SceneGraph geometry. More...
 
 ~VPolytope () final
 
VPolytope GetMinimalRepresentation () const
 Creates a new VPolytope whose vertices are guaranteed to be minimal, i.e., if we remove any point from its vertices, then the convex hull of the remaining vertices is a strict subset of the polytope. More...
 
const Eigen::MatrixXd & vertices () const
 Returns the vertices in a d-by-n matrix, where d is the ambient dimension, and n is the number of vertices. More...
 
void WriteObj (const std::filesystem::path &filename) const
 Uses qhull to compute the Delaunay triangulation and then writes the vertices and faces to filename in the Wavefront Obj format. More...
 
bool PointInSet (const Eigen::Ref< const Eigen::VectorXd > &x, double tol=0) const
 Returns true if the point is within tol of the set under the L∞-norm. More...
 
double CalcVolume () const
 Computes the volume of this V-Polytope. More...
 
Implements CopyConstructible, CopyAssignable, MoveConstructible, MoveAssignable
 VPolytope (const VPolytope &)=default
 
VPolytopeoperator= (const VPolytope &)=default
 
 VPolytope (VPolytope &&)=default
 
VPolytopeoperator= (VPolytope &&)=default
 
- Public Member Functions inherited from ConvexSet
virtual ~ConvexSet ()
 
std::unique_ptr< ConvexSetClone () const
 Creates a unique deep copy of this set. More...
 
int ambient_dimension () const
 Returns the dimension of the vector space in which the elements of this set are evaluated. More...
 
bool IntersectsWith (const ConvexSet &other) const
 Returns true iff the intersection between this and other is non-empty. More...
 
bool IsBounded () const
 Returns true iff the set is bounded, e.g., there exists an element-wise finite lower and upper bound for the set. More...
 
bool IsEmpty () const
 Returns true iff the set is empty. More...
 
std::optional< Eigen::VectorXd > MaybeGetPoint () const
 If this set trivially contains exactly one point, returns the value of that point. More...
 
std::optional< Eigen::VectorXd > MaybeGetFeasiblePoint () const
 Returns a feasible point within this convex set if it is nonempty, and nullopt otherwise. More...
 
bool PointInSet (const Eigen::Ref< const Eigen::VectorXd > &x, double tol=0) const
 Returns true iff the point x is contained in the set. More...
 
std::pair< VectorX< symbolic::Variable >, std::vector< solvers::Binding< solvers::Constraint > > > AddPointInSetConstraints (solvers::MathematicalProgram *prog, const Eigen::Ref< const solvers::VectorXDecisionVariable > &vars) const
 Adds a constraint to an existing MathematicalProgram enforcing that the point defined by vars is inside the set. More...
 
std::vector< solvers::Binding< solvers::Constraint > > AddPointInNonnegativeScalingConstraints (solvers::MathematicalProgram *prog, const Eigen::Ref< const solvers::VectorXDecisionVariable > &x, const symbolic::Variable &t) const
 Let S be this convex set. More...
 
std::vector< solvers::Binding< solvers::Constraint > > AddPointInNonnegativeScalingConstraints (solvers::MathematicalProgram *prog, const Eigen::Ref< const Eigen::MatrixXd > &A, const Eigen::Ref< const Eigen::VectorXd > &b, const Eigen::Ref< const Eigen::VectorXd > &c, double d, const Eigen::Ref< const solvers::VectorXDecisionVariable > &x, const Eigen::Ref< const solvers::VectorXDecisionVariable > &t) const
 Let S be this convex set. More...
 
std::pair< std::unique_ptr< Shape >, math::RigidTransformd > ToShapeWithPose () const
 Constructs a Shape and a pose of the set in the world frame for use in the SceneGraph geometry ecosystem. More...
 
double CalcVolume () const
 Computes the exact volume for the convex set. More...
 
SampledVolume CalcVolumeViaSampling (RandomGenerator *generator, const double desired_rel_accuracy=1e-2, const int max_num_samples=1e4) const
 Calculates an estimate of the volume of the convex set using sampling and performing Monte Carlo integration. More...
 
bool has_exact_volume () const
 Returns true if the exact volume can be computed for this convex set instance. More...
 

Static Public Member Functions

static VPolytope MakeBox (const Eigen::Ref< const Eigen::VectorXd > &lb, const Eigen::Ref< const Eigen::VectorXd > &ub)
 Constructs a polyhedron as an axis-aligned box from the lower and upper corners. More...
 
static VPolytope MakeUnitBox (int dim)
 Constructs the L∞-norm unit box in dim dimensions, {x | |x|∞ <= 1 }. More...
 

Additional Inherited Members

- Protected Member Functions inherited from ConvexSet
 ConvexSet (int ambient_dimension, bool has_exact_volume)
 For use by derived classes to construct a ConvexSet. More...
 
template<typename Archive >
void Serialize (Archive *a)
 Implements non-virtual base class serialization. More...
 
std::optional< symbolic::VariableHandleZeroAmbientDimensionConstraints (solvers::MathematicalProgram *prog, const ConvexSet &set, std::vector< solvers::Binding< solvers::Constraint >> *constraints) const
 Instances of subclasses such as CartesianProduct and MinkowskiSum can have constituent sets with zero ambient dimension, which much be handled in a special manner when calling methods such as DoAddPointInSetConstraints. More...
 
 ConvexSet (const ConvexSet &)=default
 
ConvexSetoperator= (const ConvexSet &)=default
 
 ConvexSet (ConvexSet &&)=default
 
ConvexSetoperator= (ConvexSet &&)=default
 

Constructor & Destructor Documentation

◆ VPolytope() [1/6]

VPolytope ( const VPolytope )
default

◆ VPolytope() [2/6]

VPolytope ( VPolytope &&  )
default

◆ VPolytope() [3/6]

VPolytope ( )

Constructs a set with no vertices in the zero-dimensional space, which is empty (by convention).

◆ VPolytope() [4/6]

VPolytope ( const Eigen::Ref< const Eigen::MatrixXd > &  vertices)
explicit

Constructs the polytope from a d-by-n matrix, where d is the ambient dimension, and n is the number of vertices.

The vertices do not have to be ordered, nor minimal (they can contain points inside their convex hull).

◆ VPolytope() [5/6]

VPolytope ( const HPolyhedron H,
const double  tol = 1e-9 
)
explicit

Constructs the polytope from a bounded polyhedron (using Qhull).

If the HPolyhedron is not full-dimensional, we perform computations in a coordinate system of its affine hull. tol specifies the numerical tolerance used in the computation of the affine hull. See the documentation of AffineSubspace for more details. A loose tolerance is necessary for the built-in solvers, but a tighter tolerance can be used with commercial solvers (e.g. Gurobi and Mosek).

Exceptions
std::runtime_errorif H is unbounded or if Qhull terminates with an error.

◆ VPolytope() [6/6]

VPolytope ( const QueryObject< double > &  query_object,
GeometryId  geometry_id,
std::optional< FrameId reference_frame = std::nullopt 
)

Constructs the polytope from a SceneGraph geometry.

◆ ~VPolytope()

~VPolytope ( )
final

Member Function Documentation

◆ CalcVolume()

double CalcVolume

Computes the volume of this V-Polytope.

Note
this function calls qhull to compute the volume.

◆ GetMinimalRepresentation()

VPolytope GetMinimalRepresentation ( ) const

Creates a new VPolytope whose vertices are guaranteed to be minimal, i.e., if we remove any point from its vertices, then the convex hull of the remaining vertices is a strict subset of the polytope.

In the 2D case the vertices of the new VPolytope are ordered counter-clockwise from the negative X axis. For all other cases an order is not guaranteed.

◆ MakeBox()

static VPolytope MakeBox ( const Eigen::Ref< const Eigen::VectorXd > &  lb,
const Eigen::Ref< const Eigen::VectorXd > &  ub 
)
static

Constructs a polyhedron as an axis-aligned box from the lower and upper corners.

◆ MakeUnitBox()

static VPolytope MakeUnitBox ( int  dim)
static

Constructs the L∞-norm unit box in dim dimensions, {x | |x|∞ <= 1 }.

This is an axis-aligned box, centered at the origin, with edge length 2.

◆ operator=() [1/2]

VPolytope& operator= ( VPolytope &&  )
default

◆ operator=() [2/2]

VPolytope& operator= ( const VPolytope )
default

◆ PointInSet()

bool PointInSet

Returns true if the point is within tol of the set under the L∞-norm.

Note: This requires the solution of a linear program; the achievable tolerance may be dependent on your specific solver and solver parameters.

See also
ConvexSet::set_solver

◆ vertices()

const Eigen::MatrixXd& vertices ( ) const

Returns the vertices in a d-by-n matrix, where d is the ambient dimension, and n is the number of vertices.

◆ WriteObj()

void WriteObj ( const std::filesystem::path &  filename) const

Uses qhull to compute the Delaunay triangulation and then writes the vertices and faces to filename in the Wavefront Obj format.

Note that the extension .obj is not automatically added to the filename.

Precondition
ambient_dimension() == 3.

The documentation for this class was generated from the following file: