Drake
Drake C++ Documentation

Detailed Description

Abstract base class for defining a convex set.

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

Public Member Functions

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...
 
std::optional< std::pair< std::vector< double >, Eigen::MatrixXd > > Projection (const Eigen::Ref< const Eigen::MatrixXd > &points) const
 Computes in the L₂ norm the distance and the nearest point in this convex set to every column of points. More...
 
bool has_exact_volume () const
 Returns true if the exact volume can be computed for this convex set instance. More...
 

Protected Member Functions

 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...
 
virtual std::unique_ptr< ConvexSetDoClone () const =0
 Non-virtual interface implementation for Clone(). More...
 
virtual std::optional< bool > DoIsBoundedShortcut () const
 Non-virtual interface implementation for DoIsBoundedShortcut(). More...
 
virtual std::vector< std::optional< double > > DoProjectionShortcut (const Eigen::Ref< const Eigen::MatrixXd > &points, EigenPtr< Eigen::MatrixXd > projected_points) const
 Non-virtual interface implementation for DoProjectionShortcut(). More...
 
virtual bool DoIsEmpty () const
 Non-virtual interface implementation for IsEmpty(). More...
 
virtual std::optional< Eigen::VectorXd > DoMaybeGetPoint () const
 Non-virtual interface implementation for MaybeGetPoint(). More...
 
virtual std::optional< Eigen::VectorXd > DoMaybeGetFeasiblePoint () const
 Non-virtual interface implementation for MaybeGetFeasiblePoint(). More...
 
virtual bool DoPointInSet (const Eigen::Ref< const Eigen::VectorXd > &x, double tol) const
 Non-virtual interface implementation for PointInSet(). More...
 
virtual std::optional< bool > DoPointInSetShortcut (const Eigen::Ref< const Eigen::VectorXd > &x, double tol) const
 A non-virtual interface implementation for PointInSet() that should be used when the PointInSet() can be computed more efficiently than solving a convex program. More...
 
virtual std::pair< VectorX< symbolic::Variable >, std::vector< solvers::Binding< solvers::Constraint > > > DoAddPointInSetConstraints (solvers::MathematicalProgram *prog, const Eigen::Ref< const solvers::VectorXDecisionVariable > &vars) const =0
 Non-virtual interface implementation for AddPointInSetConstraints(). More...
 
virtual std::vector< solvers::Binding< solvers::Constraint > > DoAddPointInNonnegativeScalingConstraints (solvers::MathematicalProgram *prog, const Eigen::Ref< const solvers::VectorXDecisionVariable > &x, const symbolic::Variable &t) const =0
 Non-virtual interface implementation for AddPointInNonnegativeScalingConstraints(). More...
 
virtual std::vector< solvers::Binding< solvers::Constraint > > DoAddPointInNonnegativeScalingConstraints (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 =0
 Non-virtual interface implementation for AddPointInNonnegativeScalingConstraints(). More...
 
virtual std::pair< std::unique_ptr< Shape >, math::RigidTransformd > DoToShapeWithPose () const =0
 Non-virtual interface implementation for ToShapeWithPose(). More...
 
virtual double DoCalcVolume () const
 Non-virtual interface implementation for CalcVolume(). 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...
 
Implements CopyConstructible, CopyAssignable, MoveConstructible, MoveAssignable
 ConvexSet (const ConvexSet &)=default
 
ConvexSetoperator= (const ConvexSet &)=default
 
 ConvexSet (ConvexSet &&)=default
 
ConvexSetoperator= (ConvexSet &&)=default
 

Constructor & Destructor Documentation

◆ ~ConvexSet()

virtual ~ConvexSet ( )
virtual

◆ ConvexSet() [1/3]

ConvexSet ( const ConvexSet )
protecteddefault

◆ ConvexSet() [2/3]

ConvexSet ( ConvexSet &&  )
protecteddefault

◆ ConvexSet() [3/3]

ConvexSet ( int  ambient_dimension,
bool  has_exact_volume 
)
explicitprotected

For use by derived classes to construct a ConvexSet.

Parameters
has_exact_volumeDerived classes should pass true if they've implemented DoCalcVolume() to return a value (at least sometimes).

Member Function Documentation

◆ AddPointInNonnegativeScalingConstraints() [1/2]

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.

When S is bounded, this method adds the convex constraints to imply

x ∈ t S,
t ≥ 0,

where x is a point in ℜⁿ (with n the ambient_dimension) and t is a scalar.

When S is unbounded, then the behavior is almost identical, except when t=0. In this case, the constraints imply t ≥ 0, x ∈ t S ⊕ rec(S), where rec(S) is the recession cone of S (the asymptotic directions in which S is not bounded) and ⊕ is the Minkowski sum. For t > 0, this is equivalent to x ∈ t S, but for t = 0, we have only x ∈ rec(S).

Exceptions
std::exceptionif ambient_dimension() == 0

◆ AddPointInNonnegativeScalingConstraints() [2/2]

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.

When S is bounded, this method adds the convex constraints to imply

A * x + b ∈ (c' * t + d) S,
c' * t + d ≥ 0,

where A is an n-by-m matrix (with n the ambient_dimension), b is a vector of size n, c is a vector of size p, x is a point in ℜᵐ, and t is a point in ℜᵖ.

When S is unbounded, then the behavior is almost identical, except when c' * t+d=0. In this case, the constraints imply

A * x + b ∈ (c' * t + d) S ⊕ rec(S),
c' * t + d ≥ 0,

where rec(S) is the recession cone of S (the asymptotic directions in which S is not bounded) and ⊕ is the Minkowski sum. For c' * t + d > 0, this is equivalent to A * x + b ∈ (c' * t + d) S, but for c' * t + d = 0, we have only A * x + b ∈ rec(S).

Exceptions
std::exceptionif ambient_dimension() == 0

◆ AddPointInSetConstraints()

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.

Returns
(new_vars, new_constraints) Some of the derived class will add new decision variables to enforce this constraint, we return all the newly added decision variables as new_vars. The meaning of these new decision variables differs in each subclass. If no new variables are added, then we return an empty Eigen vector. Also we return all the newly added constraints to prog through this function.
Exceptions
std::exceptionif ambient_dimension() == 0

◆ ambient_dimension()

int ambient_dimension ( ) const

Returns the dimension of the vector space in which the elements of this set are evaluated.

Contrast this with the affine dimension: the dimension of the smallest affine subset of the ambient space that contains our set. For example, if we define a set using A*x = b, where A has linearly independent rows, then the ambient dimension is the dimension of x, but the affine dimension of the set is ambient_dimension() - rank(A).

◆ CalcVolume()

double CalcVolume ( ) const

Computes the exact volume for the convex set.

Note
Not every convex set can report an exact volume. In that case, use CalcVolumeViaSampling() instead.
Exceptions
std::exceptionif has_exact_volume() returns false.
ifambient_dimension() == 0.

◆ CalcVolumeViaSampling()

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.

Note
this method is intended to be used for low to moderate dimensions (d<15). For larger dimensions, a telescopic product approach has yet to be implemented. See, e.g., https://proceedings.mlr.press/v151/chevallier22a/chevallier22a.pdf
Parameters
generatora random number generator.
desired_rel_accuracythe desired relative accuracy of the volume estimate in the sense that the estimated volume is likely to be within the interval defined by (1±2*desired_rel_accuracy)*true_volume with probability of at least* 0.95 according to the Law of Large Numbers. https://people.math.umass.edu/~lr7q/ps_files/teaching/math456/Chapter6.pdf The computation will terminate when the relative error is less than rel_accuracy or when the maximum number of samples is reached.
max_num_samplesthe maximum number of samples to use.
Precondition
desired_rel_accuracy is in the range [0,1].
Returns
a pair the estimated volume of the set and an upper bound for the relative accuracy
Exceptions
ifambient_dimension() == 0.
ifthe minimum axis-aligned bounding box of the set cannot be computed.

◆ Clone()

std::unique_ptr<ConvexSet> Clone ( ) const

Creates a unique deep copy of this set.

◆ DoAddPointInNonnegativeScalingConstraints() [1/2]

virtual std::vector<solvers::Binding<solvers::Constraint> > DoAddPointInNonnegativeScalingConstraints ( solvers::MathematicalProgram prog,
const Eigen::Ref< const solvers::VectorXDecisionVariable > &  x,
const symbolic::Variable t 
) const
protectedpure virtual

Non-virtual interface implementation for AddPointInNonnegativeScalingConstraints().

Precondition
x.size() == ambient_dimension()
ambient_dimension() > 0

◆ DoAddPointInNonnegativeScalingConstraints() [2/2]

virtual std::vector<solvers::Binding<solvers::Constraint> > DoAddPointInNonnegativeScalingConstraints ( 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
protectedpure virtual

Non-virtual interface implementation for AddPointInNonnegativeScalingConstraints().

Subclasses must override to add the constraints needed to keep the point A * x + b in the non-negative scaling of the set. Note that subclasses do not need to add the constraint c * t + d ≥ 0 as it is already added.

Precondition
ambient_dimension() > 0
A.rows() == ambient_dimension()
A.rows() == b.rows()
A.cols() == x.size()
c.rows() == t.size()

◆ DoAddPointInSetConstraints()

virtual std::pair<VectorX<symbolic::Variable>, std::vector<solvers::Binding<solvers::Constraint> > > DoAddPointInSetConstraints ( solvers::MathematicalProgram prog,
const Eigen::Ref< const solvers::VectorXDecisionVariable > &  vars 
) const
protectedpure virtual

Non-virtual interface implementation for AddPointInSetConstraints().

Precondition
vars.size() == ambient_dimension()
ambient_dimension() > 0

◆ DoCalcVolume()

virtual double DoCalcVolume ( ) const
protectedvirtual

Non-virtual interface implementation for CalcVolume().

This will only be called if has_exact_volume() returns true and ambient_dimension() > 0

◆ DoClone()

virtual std::unique_ptr<ConvexSet> DoClone ( ) const
protectedpure virtual

Non-virtual interface implementation for Clone().

◆ DoIsBoundedShortcut()

virtual std::optional<bool> DoIsBoundedShortcut ( ) const
protectedvirtual

Non-virtual interface implementation for DoIsBoundedShortcut().

Trivially returns std::nullopt. This allows a derived class to implement its own boundedness checks, to potentially avoid the more expensive base class checks.

Precondition
ambient_dimension() >= 0

◆ DoIsEmpty()

virtual bool DoIsEmpty ( ) const
protectedvirtual

Non-virtual interface implementation for IsEmpty().

The default implementation solves a feasibility optimization problem, but derived classes can override with a custom (more efficient) implementation. Zero-dimensional sets are considered to be nonempty by default. Sets which can be zero-dimensional and empty must handle this behavior in their derived implementation of DoIsEmpty.

◆ DoMaybeGetFeasiblePoint()

virtual std::optional<Eigen::VectorXd> DoMaybeGetFeasiblePoint ( ) const
protectedvirtual

Non-virtual interface implementation for MaybeGetFeasiblePoint().

The default implementation solves a feasibility optimization problem, but derived classes can override with a custom (more efficient) implementation.

◆ DoMaybeGetPoint()

virtual std::optional<Eigen::VectorXd> DoMaybeGetPoint ( ) const
protectedvirtual

Non-virtual interface implementation for MaybeGetPoint().

The default implementation returns nullopt. Sets that can model a single point should override with a custom implementation.

Precondition
ambient_dimension() >= 0.

◆ DoPointInSet()

virtual bool DoPointInSet ( const Eigen::Ref< const Eigen::VectorXd > &  x,
double  tol 
) const
protectedvirtual

Non-virtual interface implementation for PointInSet().

Precondition
x.size() == ambient_dimension()
ambient_dimension() >= 0

◆ DoPointInSetShortcut()

virtual std::optional<bool> DoPointInSetShortcut ( const Eigen::Ref< const Eigen::VectorXd > &  x,
double  tol 
) const
protectedvirtual

A non-virtual interface implementation for PointInSet() that should be used when the PointInSet() can be computed more efficiently than solving a convex program.

Returns
Returns true if and only if x is known to be in the set. Returns false if and only if x is known to not be in the set. Returns std::nullopt if a shortcut implementation is not provided (i.e. the method has not elected to decide whether the point x is in the set).

For example, membership in a VPolytope cannot be verified without solving a linear program and so no shortcut implementation should be provided. On the other hand, membership in an HPolyhedron can be checked by checking the inequality Ax ≤ b and so a shortcut is possible.

◆ DoProjectionShortcut()

virtual std::vector<std::optional<double> > DoProjectionShortcut ( const Eigen::Ref< const Eigen::MatrixXd > &  points,
EigenPtr< Eigen::MatrixXd >  projected_points 
) const
protectedvirtual

Non-virtual interface implementation for DoProjectionShortcut().

This allows a derived class to implement a method which computes the projection of some, but not necessarily all, of the points more efficiently than the generic implementation.

The default implementation checks whether each column of points is in the set using DoPointInSetShortcut. Points in the set are given a distance of 0 and are projected to themselves.

Parameters
[in]pointsare the points which we wish to project to the convex set.
[in/out]projected_points are the projection of points onto the convex set.
Returns
A vector distances which is the same size as points.cols().These are the distances from points to the convex set. If distances[i] has a value, then projected_points->col(i) is the projection of points.col(i) onto the set. If distances[i] is nullopt, then the projection of points.col(i) has not yet been computed, and so the value at projected_points->col(i) is meaningless.
Precondition
ambient_dimension() >= 0
distances.size() == points.cols()

◆ DoToShapeWithPose()

virtual std::pair<std::unique_ptr<Shape>, math::RigidTransformd> DoToShapeWithPose ( ) const
protectedpure virtual

Non-virtual interface implementation for ToShapeWithPose().

Precondition
ambient_dimension() == 3

◆ HandleZeroAmbientDimensionConstraints()

std::optional<symbolic::Variable> HandleZeroAmbientDimensionConstraints ( solvers::MathematicalProgram prog,
const ConvexSet set,
std::vector< solvers::Binding< solvers::Constraint >> *  constraints 
) const
protected

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.

If the set is empty, a trivially infeasible constraint must be added. We also warn the user when this happens, since they probably didn't intend it to occur. If the set is nonempty, then it's the unique zero-dimensional vector space {0}, and no additional variables or constraints are needed. If a new variable is created, return it, to optionally be stored (as in AddPointInSetConstraints), or not be stored (as in DoAddPointInNonnegativeScalingConstraints).

◆ has_exact_volume()

bool has_exact_volume ( ) const

Returns true if the exact volume can be computed for this convex set instance.

Note
This value reasons about to the generic case of the convex set class rather than the specific instance of the convex set. For example, the exact volume of a box is trivival to compute, but if the box is created as a HPolyhedron, then the exact volume cannot be computed.

◆ IntersectsWith()

bool IntersectsWith ( const ConvexSet other) const

Returns true iff the intersection between this and other is non-empty.

Exceptions
std::exceptionif the ambient dimension of other is not the same as that of this.

◆ IsBounded()

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.

Note: for some derived classes, this check is trivial, but for others it can require solving an (typically small) optimization problem. Check the derived class documentation for any notes.

◆ IsEmpty()

bool IsEmpty ( ) const

Returns true iff the set is empty.

Note: for some derived classes, this check is trivial, but for others, it can require solving a (typically small) optimization problem. Check the derived class documentation for any notes. Zero-dimensional sets must be handled specially. There are two possible sets in a zero-dimensional space – the empty set, and the whole set (which is simply the "zero vector space", {0}.) For more details, see: https://en.wikipedia.org/wiki/Examples_of_vector_spaces#Trivial_or_zero_vector_space Zero-dimensional sets are considered to be nonempty by default. Sets which can be zero-dimensional and empty must handle this behavior in their derived implementation of DoIsEmpty. An example of such a subclass is VPolytope.

◆ MaybeGetFeasiblePoint()

std::optional<Eigen::VectorXd> MaybeGetFeasiblePoint ( ) const

Returns a feasible point within this convex set if it is nonempty, and nullopt otherwise.

◆ MaybeGetPoint()

std::optional<Eigen::VectorXd> MaybeGetPoint ( ) const

If this set trivially contains exactly one point, returns the value of that point.

Otherwise, returns nullopt. By "trivially", we mean that representation of the set structurally maps to a single point; if checking for point-ness would require solving an optimization program, returns nullopt. In other words, this is a relatively cheap function to call.

◆ operator=() [1/2]

ConvexSet& operator= ( ConvexSet &&  )
protecteddefault

◆ operator=() [2/2]

ConvexSet& operator= ( const ConvexSet )
protecteddefault

◆ PointInSet()

bool PointInSet ( const Eigen::Ref< const Eigen::VectorXd > &  x,
double  tol = 0 
) const

Returns true iff the point x is contained in the set.

If the ambient dimension is zero, then if the set is nonempty, the point is trivially in the set, and if the set is empty, the point is trivially not in the set.

◆ Projection()

std::optional<std::pair<std::vector<double>, Eigen::MatrixXd> > Projection ( const Eigen::Ref< const Eigen::MatrixXd > &  points) const

Computes in the L₂ norm the distance and the nearest point in this convex set to every column of points.

If this set is empty, we return nullopt.

Precondition
points.rows() == ambient_dimension().
Exceptions
ifthe internal convex optimization solver fails.

◆ Serialize()

void Serialize ( Archive *  a)
protected

Implements non-virtual base class serialization.

◆ ToShapeWithPose()

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.

Exceptions
std::exceptionif ambient_dimension() != 3 or if the functionality for a particular set has not been implemented yet.

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