mpt4py package

Contents

mpt4py package#

Subpackages#

Submodules#

mpt4py.base module#

class mpt4py.base.HData(A: ndarray[tuple[Any, ...], dtype[float64]] | None = None, b: ndarray[tuple[Any, ...], dtype[float64]] | None = None, Ae: ndarray[tuple[Any, ...], dtype[float64]] | None = None, be: ndarray[tuple[Any, ...], dtype[float64]] | None = None, is_minimal: bool = False, preprocess: bool = True)[source]#

Bases: object

The halfspace representation (H-representation) of a polyhedron, defined as:

\[\{ x\in\mathbb{R}^n \mid A x \leq b, A_e x = b_e \}.\]

Note

If preprocess is set as True, the __init__ method performs pre-processing on the provided arguments. This might cause the A, b, Ae, and be attributes to be different from what the user provided in arguments.

A: ndarray[tuple[Any, ...], dtype[float64]]#
Ae: ndarray[tuple[Any, ...], dtype[float64]]#
b: ndarray[tuple[Any, ...], dtype[float64]]#
be: ndarray[tuple[Any, ...], dtype[float64]]#
property dehomogenize: HData#
property dim: int#
property homogenize: HData#
property is_cone#

Returns True if the Polyhedron is a cone

is_minimal: bool = False#
preprocess()[source]#

Preprocess the H-representation by removing redundant constraints and normalizing the constraints.

  • Deal with zero rows in \(A\):

    If at least one row \(a_i\) in \(A\) is almost zero and its corresponding \(b_i\) is strictly negative, the provided H-presentation will be treated as an empty set, represented by \(\{ x\in\mathbb{R}^n \mid 0_{1\times n}^\top x \leq -1\}\).

    If at least one row \(a_i\) in \(A\) is almost zero and its corresponding \(b_i\) is non-negative, this row will be treated as redundant and removed, unless it is the only row in \(A\). In this case, the provided H-presentation will be treated as the full space, represented by \(\{ x\in\mathbb{R}^n \mid 0_{1\times n}^\top x \leq 1\}\).

  • Deal with zero rows in \(A_e\):

    If at least one row \(a_{e,i}\) in \(A_e\) is almost zero and its corresponding \(b_{e,i}\) is non-zero, the provided H-presentation will be treated as an empty set, represented by \(\{ x\in\mathbb{R}^n \mid 0_{1\times n}^\top x \leq -1\}\).

    If at least one row \(a_{e,i}\) in \(A_e\) is almost zero and its corresponding \(b_{e,i}\) is zero, this row will be treated as redundant and removed.

  • Deal with rows that are multiples of other rows in \(A\) and \(A_e\).

    For \(A\), if there exist \(i<j\), s.t. \(a_i = t a_j, ~b_i = t b_j\) where \(t \neq 0\), then \(a_j\) is treated as redundant and will be removed. Similar for \(A_e\).

Note

This method modifies the H-representation in place.

class mpt4py.base.VData(V: ndarray[tuple[Any, ...], dtype[float64]] | None = None, R: ndarray[tuple[Any, ...], dtype[float64]] | None = None, is_minimal: bool = False, preprocess: bool = True)[source]#

Bases: object

The vertex representation (V-representation) of a polyhedron, defined as:

\[\{ x\in\mathbb{R}^n \mid x = V^\top \lambda + R^\top \mu, \lambda,\mu \geq 0, 1^\top \lambda = 1 \}.\]

Note

If preprocess is set as True, the __init__ method performs pre-processing on the provided arguments. This might cause the V and R attributes to be different from what the user provided in arguments.

R: ndarray[tuple[Any, ...], dtype[float64]]#
V: ndarray[tuple[Any, ...], dtype[float64]]#
property dehomogenize: VData#

Be careful when calling this function. It only makes sense if the V-representation is homogenous.

property dim: int#
property homogenize: VData#
property is_cone#

Returns True if the Polyhedron is a cone

is_minimal: bool = False#
preprocess()[source]#

Preprocess the V-representation by removing obviously redundant vertices and rays.

  • Deal with repeated rows in \(V\):

    If there are repeated rows in \(V\), i.e., \(v_i = v_j\) where \(i<j\), then \(v_j\) will be treated as a redundant vertex and removed.

  • Deal with zero rows in \(R\):

    If there are zero rows in \(R\), these zeros rows will be removed.

  • Deal with rows that are positive multiples of others in \(R\):

    If there exist \(i<j\), s.t. \(r_i = t r_j\) where \(t > 0\), \(r_j\) will be treated as a redundant ray and removed.

mpt4py.exceptions module#

exception mpt4py.exceptions.DimMismatchError(dim_expected, dim_received)[source]#

Bases: Exception

Raised for dimension mismatches.

exception mpt4py.exceptions.GeometryError[source]#

Bases: Exception

exception mpt4py.exceptions.InfeasibleError[source]#

Bases: Exception

Raised when a set is empty.

exception mpt4py.exceptions.UnboundedError[source]#

Bases: Exception

Raised when a set is unbounded.

mpt4py.tolerances module#

mpt4py.tolerances.get_tolerance(name: str)[source]#

Retrieve a tolerance value.

mpt4py.tolerances.print_tolerances()[source]#

Print all current tolerance values.

mpt4py.tolerances.set_tolerance(name: str, value: float)[source]#

Modify a tolerance at runtime.

mpt4py.util module#

class mpt4py.util.Basis(indices: List[int], n: int | None = None)[source]#

Bases: object

Define the basis of a Linear Complementarity Problem (LCP)

Any set \(B\) which is a subset of \(\{1, ..., 2n\}\) such that \(|B| = n\) and \(\text{rank} A_{*,B} = n\) is called a basis.

property indices#
is_almost_complementary() → bool[source]#

Determine if the current basis is almost complementary

is_complementary() → bool[source]#

Determine if the current basis B is a complementary basis, i.e., for all i in B, i+n and i-n are not in B

is_feasible(M: ndarray[tuple[Any, ...], dtype[float64]], q: ndarray[tuple[Any, ...], dtype[float64]]) → bool[source]#

Determine if the basis is feasible regrading a LCP (M, q)

property n#
mpt4py.util.is_lexico_positive(matrix: ndarray[tuple[Any, ...], dtype[float64]], zero_tol: float = tolerances['zero']) → bool[source]#

Tell if a matrix or vector is lexico positive.

  • A vector is lexico positive if its first non-zero elements is nonnegative (following the style in MPT3).

  • A matrix is lexico positive if all its rows are lexico positive.

Parameters:

  • matrix – a ndarray of size [m,n]. The array can be a 1d array representing a vector.

  • zero_tol – tolerance for a number to be recognized as non-zero.

Returns:

a bool variable implying the matrix is lexico positive or not.

mpt4py.util.is_p_matrix(A: ndarray[tuple[Any, ...], dtype[float64]]) → bool[source]#

Check if a matrix is a P-matrix.

A matrix is called P-matrix if its all principal minors are strictly positive. Refer to https://www.tuhh.de/ti3/paper/rump/Rump03b.pdf.

mpt4py.util.is_positive_definite(A: ndarray[tuple[Any, ...], dtype[float64]]) → bool[source]#

Check if a real matrix is positive definite (P.D.) through Cholesky decomposition.

Since \(x^\top Ax = (x^\top Ax)^\top = x^\top A^\top x, x^\top (A+A^\top)x = 2x^\top Ax\). Therefore, \(A\) is P.D. if and only if \(A^\top+A\) is P.D.

mpt4py.util.is_positive_semidefinite(A: ndarray[tuple[Any, ...], dtype[float64]], hermitian: bool = False) → bool[source]#

Check if a real matrix is positive semi-definite (P.S.D.) through Cholesky decomposition.

Since \(x^\top Ax = (x^\top Ax)^\top = x^\top A^\top x, x^\top (A+A^\top)x = 2x^\top Ax\). Therefore, \(A\) is P.S.D. if and only if \(A^\top+A\) is P.S.D.

mpt4py.util.is_symmetric(A: ndarray[tuple[Any, ...], dtype[float64]], zero_tol: float = tolerances['zero']) → bool[source]#

Check if a real matrix is symmetric or not.

mpt4py.util.lexico_argmax(matrix: ndarray[tuple[Any, ...], dtype[float64]], zero_tol: float = tolerances['zero']) → List[int][source]#

Find the (row) index of the lexico maximum of a matrix.

The lexico maximum of a set of vectors \(\{a_1, a_2, ..., a_m\}\) is the vector \(a_j\) satisfying that \(a_j - a_i\) is lexico positive for all \(i \in \{1,2,...,m\}/\{j\}\).

Parameters:

  • matrix – ndarray.

  • zero_tol – the tolerance to determine if two real number are identical or not.

Returns:

An integer vector consisting of index of lexico maximum row (the index could be non-unique).

mpt4py.util.lexico_argmin(matrix: ndarray[tuple[Any, ...], dtype[float64]], zero_tol: float = tolerances['zero']) → List[int][source]#

Find the (row) index of the lexico minimum of a matrix.

The lexico minimum of a set of vectors \(\{a_1, a_2, ..., a_m\}\) is the vector \(a_j\) satisfying that \(a_i - a_j\) is lexico positive for all \(i \in \{1,2,...,m\}/\{j\}\).

Parameters:

  • matrix – ndarray of size [m,n].

  • zero_tol – the tolerance to determine if two real number are identical or not.

Returns:

A list of indices corresponding to the rows that are

Return type:

List[int]

lexicographical minima (the indices may be non-unique).

mpt4py.util.range_space(A: ndarray[tuple[Any, ...], dtype[float64]])[source]#

Compute the range-space of a matrix.

Module contents#

The Multi-Parametric Toolbox for Python (mpt4py) Provides computational geometry, multi-parametric programming, and MPC synthesis tools.

class mpt4py.AffineSet(A: ndarray[tuple[Any, ...], dtype[float64]], b: ndarray[tuple[Any, ...], dtype[float64]])[source]#

Bases: ConvexSet

Represents an affine set:

\[\{ x \in \mathbb{R}^n \mid Ax=b \}\]
property A: ndarray[tuple[Any, ...], dtype[float64]]#
property b: ndarray[tuple[Any, ...], dtype[float64]]#
property basis: Tuple[ndarray[tuple[Any, ...], dtype[float64]], ndarray[tuple[Any, ...], dtype[float64]]]#

Returns a basis and solution for the affine set such that

\[Ax = b \iff x = B y + x_0\]
property coefficients: Tuple[ndarray[tuple[Any, ...], dtype[float64]], ndarray[tuple[Any, ...], dtype[float64]]]#
property dim: int#

Returns the dimension of the convex set.

static from_basis(basis: ndarray[tuple[Any, ...], dtype[float64]], x0: ndarray[tuple[Any, ...], dtype[float64]]) → AffineSet[source]#

Construct an affine set from a basis and a point in the set.

Parameters: - basis (array-like): Basis for the affine set. - x0 (array-like): Point in the affine set.

Returns: - AffineSet: Affine set defined by Ax == b.

homogenize()[source]#

Homogenize the affine set

minimal_representation() → AffineSet[source]#

Returns a minimal representation of the affine set

property x0: ndarray[tuple[Any, ...], dtype[float64]]#

Returns a point in the affine set

class mpt4py.ConvexSet[source]#

Bases: ABC

add_function(func: FunctionBase | Callable[[ndarray[tuple[Any, ...], dtype[float64]]], float], func_name: str) → None[source]#

Adds a function associated to the convex set.

abstract property dim: int#

Returns the dimension of the convex set.

distance(other: ndarray[tuple[Any, ...], dtype[float64]] | ConvexSet) → floating[source]#

Compute the distance between this convex set \(\mathcal{C}\) and the given point \(z\) or convex set \(\mathcal{S}\).

By providing real vector \(z\), the distance between \(z\) and \(\mathcal{C}\) is computed by solving the optimization problem

\[\begin{split}\begin{align*} \min\limits_y &\quad \|z-y\|_2^2 \\ \text{s.t.} &\quad y \in \mathcal{C}. \end{align*}\end{split}\]

By providing a convex set \(\mathcal{S}\), the distance between \(\mathcal{C}\) and \(\mathcal{S}\) is computed by solving the optimization problem

\[\begin{split}\begin{align*} \min\limits_{x,y} &\quad \|x-y\|_2^2 \\ \text{s.t.} &\quad x \in \mathcal{C}, \\ &\quad y \in \mathcal{S}. \end{align*}\end{split}\]
extreme(x: ndarray[tuple[Any, ...], dtype[float64]]) → ndarray[tuple[Any, ...], dtype[float64]][source]#

Computes an extreme point \(y\) of the set \(\mathcal{C}\) in the direction \(x\):

\[\begin{split}\begin{align*} \mathop{\arg\min}\limits_y &\quad x^\top y \\ \text{s.t.} &\quad y \in \mathcal{C} \\ \end{align*}\end{split}\]
get_function(func_name: str) → FunctionBase[source]#

Returns a function associated to the convex set.

has_function(func_name: str | None = None) → bool[source]#

Returns whether the convex set has a function specified by the name.

list_functions()[source]#
project(x: ndarray[tuple[Any, ...], dtype[float64]]) → ndarray[tuple[Any, ...], dtype[float64]][source]#

Projects the given point \(y\) onto the set \(\mathcal{C}\):

\[\begin{split}\begin{aligned} \min_x &\quad \|x - y\|_2^2 \\ \text{s.t.} &\quad x \in \mathcal{C} \end{aligned}\end{split}\]
remove_all_functions() → None[source]#

Removes all functions associated to the convex set.

remove_function(func_name: str) → None[source]#

Removes a function associated to the convex set.

seperate(x: ndarray[tuple[Any, ...], dtype[float64]]) → Hyperplane[source]#

Compute the maximal seperating hyperplane between the point x and the Set.

shoot(r: ndarray[tuple[Any, ...], dtype[float64]], x0: ndarray[tuple[Any, ...], dtype[float64]]) → Tuple[float, ndarray[tuple[Any, ...], dtype[float64]]][source]#

Computes the solution to the optimization problem:

\[\begin{split}\begin{align*} \max\limits_\alpha &\quad \alpha \\ \text{s.t.} &\quad x_0 + \alpha r \in \mathcal{C} \end{align*}\end{split}\]

The problem can be explained as to find the point that lies along the line given by the ray \(r\) which starts from the point \(x_0\) such that it still lies inside the set \(\mathcal{C}\). Typically, the solution of this set is the point that lies on the boundary of the set, or can be Inf if the set is unbounded.

support(x: ndarray[tuple[Any, ...], dtype[float64]]) → float[source]#

Compute the support of the set \(\mathcal{C}\) in the direction \(x\):

\[\begin{split}\begin{align*} \max\limits_y &\quad x^\top y \\ \text{s.t.} &\quad y \in \mathcal{C} \\ \end{align*}\end{split}\]
class mpt4py.ConvexSetUnion(*convex_sets: ConvexSet)[source]#

Bases: object

Represents a general union of convex sets.

add(set: ConvexSet)[source]#

Adds a convex set to the union.

contains(x: ndarray[tuple[Any, ...], dtype[float64]], fastbreak: bool = False) → Tuple[bool, List[int]][source]#

Test if the point x is in the union of convex sets.

Parameters:

  • x (Vector) – A point with the same dimension as the convex sets.

  • fastbreak (bool) – Do a quick stop in the consecutive search when x is contained in the first set it founds. Defaults to False.

Returns:

whether x is in the union of convex sets or not, and the indices of the sets that contain x.

Return type:

Tuple[bool, List[int]]

property dim: int#

Returns the dimension of the convex sets in the union.

property num_sets: int#

Returns the number of convex sets in the union.

remove(index: int)[source]#

Removes the convex set at the specified index.

exception mpt4py.DimMismatchError(dim_expected, dim_received)[source]#

Bases: Exception

Raised for dimension mismatches.

class mpt4py.Ellipsoid(P: ndarray[tuple[Any, ...], dtype[float64]], center: ndarray[tuple[Any, ...], dtype[float64]] | None = None, radius: float = 1.)[source]#

Bases: ConvexSet

Represents an ellipsoid:

\[\{ x \in \mathbb{R}^n \mid (x - c)^\top P (x - c) \leq r^2 \}\]

where \(P\in\mathbb{R}^{n\times n}\) is positive semi-definite.

property dim: int#

Returns the dimension of the convex set.

plot(fig: PlotterProtocol, **kwargs)[source]#

Plot the ellipsoid.

Parameters:

fig (PlotterProtocol) – _description_

exception mpt4py.GeometryError[source]#

Bases: Exception

class mpt4py.Hexagon(center: ndarray[tuple[Any, ...], dtype[float64]] = np.zeros(2), radius: float = 1.0, rotation: float = 0.0)[source]#

Bases: Polyhedron

A class to represent a regular hexagon in 2D.

class mpt4py.Hyperplane(a: ndarray[tuple[Any, ...], dtype[float64]], b: floating | None = None, x0: ndarray[tuple[Any, ...], dtype[float64]] | None = None)[source]#

Bases: AffineSet

Represents the hyperplane: \(\langle a, x \rangle = b\) or \(\langle a, x-x_0 \rangle = 0\).

property normal#
property offset#
exception mpt4py.InfeasibleError[source]#

Bases: Exception

Raised when a set is empty.

class mpt4py.PolyUnion(*polyhedra: Polyhedron, convex: bool | None = None, overlaps: bool | None = None, connected: bool | None = None, bounded: bool | None = None, fulldim: bool | None = None)[source]#

Bases: ConvexSetUnion

Represent a union of polyhedra (represented in the same dimension).

add(polyhedron: Polyhedron)[source]#

Add a polyhedron to the union.

convex_hull() → Polyhedron[source]#

Compute the convex hull for union of polyhedra.

The convex hull of the union of polyhedra is defined as the minimal convex set that contains all polyhedra.

Returns:

the convex hull of the union of polyhedra.

Return type:

Polyhedron

fplot(ax: Axes | Plotter, func_name: str | None = None, **kwargs)[source]#

Plot the functions associated with the polyunion.

property is_bounded: bool#

Determine if the union is built from bounded polyhedra.

property is_connected: bool#

Determine if the union of polyhedra form a connected union.

property is_convex: bool#

Determine if the union of polyhedra is convex.

Warning

This method is very computationally demanding and is suitable for unions with small number of polyhedra.

property is_full_dim: bool#

Determine if the union is built from full-dimensional polyhedra.

property is_overlapping: bool#

Determine if the union of polyhedra is overlapping.

Note

This function considers following two cases to detect overlaps:

1. If two full-dimensional polyhedra overlap, then the intersection of these polyhedra must be full-dimensional.

2. If low-dimensional and full-dimensional polyhedra overlap, then the intersection of these polyhedra must not be empty.

Warning

This method is computationally demanding and is suitable for unions with small number of polyhedra.

merge()[source]#

Simplify the union of polyhedra by merging the neighboring polyhedra if their union is convex. The algorithm cycles through the regions and checks if any two regions form a convex union. If so, the algorithm combines them in one region, and continues checking the remaining regions. To improve the solution, multiple merging loops can be enabled in options.

plot(ax: Axes | Plotter, **kwargs)[source]#

Plot the Polyhedron.

Note

Require V-representation. Will be computed if necessary.

remove(index: int)[source]#

Remove the polyhedron at the specified index.

class mpt4py.Polyhedron(H: HData | None = None, V: VData | None = None)[source]#

Bases: ConvexSet

Represents a polyhedron - a set defined by the intersection of a finite number of inequalities.

property A: ndarray[tuple[Any, ...], dtype[float64]]#

Return the \(A\) matrix in H-representation. Will be computed if needed.

property Ae: ndarray[tuple[Any, ...], dtype[float64]]#

Return the \(A_e\) matrix in H-representation. Will be computed if needed.

property H: ndarray[tuple[Any, ...], dtype[float64]]#

Return the contatenated \([A \mid b]\) matrix in H-representation. Will be computed if needed.

property He: ndarray[tuple[Any, ...], dtype[float64]]#

Return the contatenated \([A_e \mid b_e]\) matrix in H-representation. Will be computed if needed.

property Hrep: HData#

Returns the H-representation of the polyhedron. Computes it if needed.

property R: ndarray[tuple[Any, ...], dtype[float64]]#

Return the \(R\) matrix in V-representation. Will be computed if needed.

property V: ndarray[tuple[Any, ...], dtype[float64]]#

Return the \(V\) matrix in V-representation. Will be computed if needed.

property Vrep: VData#

Returns the V-representation of the polyhedron. Computes it if needed.

property affine_dimension: int#

Compute the affine dimension of the polyhedron, i.e., the number of degrees of freedom you have when moving within the affine set.

Math

The implementation is based on the method proposed in Section 7.3 of Fukuda’s Polyhedral Computation.

  • For a V-polyhedron \(P = \{ x \in \mathbb{R}^n \mid x = V^\top \lambda + R^\top\mu, 1^\top \lambda = 1, \mu \geq 0, \lambda \geq 0\}\),

\[\begin{split}\text{affdim}(P) = \text{rank} \begin{bmatrix} V^\top & R^\top \\ 1_{1\times n} & 0_{1\times n} \end{bmatrix} -1.\end{split}\]

  • For a H-polyhedron \(P = \{ x \in \mathbb{R}^n \mid A x \leq b, A_e x = b_e \}\):

\[\text{affdim}(P) = n - \text{rank}\tilde{A}_e.\]

where \(\tilde{A}_e\) is the matrix containing all explicit and implicit equalities.

Note

If a polyhedron is empty, its affine dimension is defined as -1.

affine_hull(inplace: bool = False) → AffineSet[source]#

Compute the affine hull of the polyhedron, which is the smallest affine set that contains it.

This function detects if the polyhedron has an implicit affine set, e.g., if the inequalities contain \(\alpha^\top x \leq \beta\) and \(\alpha^\top x \geq \beta\).

If inplace == True, then the HData of this object is modified.

affine_map(T: ndarray[tuple[Any, ...], dtype[float64]], t: ndarray[tuple[Any, ...], dtype[float64]] | None = None) → Polyhedron[source]#

Computes an affine map \(P \rightarrow Q\) of polyhedron \(P\) to polyhedron \(Q\) based on the transformation matrix \(T\). The polyhedron \(Q\) is given by

\[Q = \{ y \in \mathbb{R}^n \mid y=Tx+t, x \in P \subseteq \mathbb{R}^d \}\]

The matrix \(T\) must be real with \(n\) rows and \(d\) columns.

  • If \(n < d\) then this operation is referred to as projection.

  • If \(n = d\) then this operation is referred to as rotation/skew.

  • If \(n > d\) then this operation is referred to as lifting.

Note

This operation is also overloaded via __matmul__ and __add__. So for a Polyhedron object P, P.affine_map(T, t) is equivalent to T @ P + t.

See also

property b: ndarray[tuple[Any, ...], dtype[float64]]#

Return the \(b\) vector in H-representation. Will be computed if needed.

property be: ndarray[tuple[Any, ...], dtype[float64]]#

Return the \(b_e\) vector in H-representation. Will be computed if needed.

static bounding_box(obj: ConvexSet) → Polyhedron[source]#

Compute the minimum bounding box of the polyhedron.

If the polyhedron is unbounded, then the box will be unbounded in those directions.

Note that if the polyhedron has an affine hull, then the bounding box will take the form

\[B = \{ x\in\mathbb{R}^n \mid x_{lb} \leq x \leq x_{ub}, A_e x = b_e \}.\]
cartesian_product(other: Polyhedron) → Polyhedron[source]#

Compute the Cartesian product with another polyhedron. The Cartesian product of two polyhedra \(P\subseteq\mathbb{R}^m\) and \(Q\subseteq\mathbb{R}^n\) is defined as:

\[P \times Q := \{ (x,y) \in\mathbb{R}^{m+n} \mid x\in P, y \in Q \}.\]
cheby_center(raise_on_empty_set: bool | None = True) → Tuple[ndarray[tuple[Any, ...], dtype[float64]], float][source]#

Compute the center and radius of the polyhedron’s Chebyshev ball.

Given a polyhedron with H-representation, the center \(c\) and radius \(r\) of the Chebyshev ball is computed by solving the following Linear Programming:

\[\begin{split}\begin{align*} \max_{c,r} &\quad r \\ \text{s.t.} &\quad a_i c + \|a_i\|_2 r \leq b_i, \quad \forall i\in \{1,...,m\} \\ &\quad A_e c = b_e \end{align*}\end{split}\]

where \(a_i\) is the \(i\) th row of \(A\) and \(m\) is the number of rows of \(A\).

Note

  • If the polyhedron is empty, the Chebyshev center is set to \(-\infty\) and the radius is set to \(-\infty\).

  • If the polyhedron is unbounded, the Chebyshev center is set to \(+\infty\) and the radius is set to \(+\infty\).

Warning

  • If the polyhedron is lower-dimensional and it only has a V-representation, then its H-representation is non-unique. Different H-representations give different Chebyshev balls.

property codimension: int#

Returns the codimension of the polyhedron, i.e., the number of constraints/dimensions by which the affine set is “lower” than the ambient space.

For example, suppose we have a non-empty polyhedron

\[P = \{ x \in \mathbb{R}^n \mid A x \leq b, A_e x = b_e \},\]

where the rank of \(A_e\) (number of rows) is \(m\). The codimension of \(P\) is \(m\).

contains(inner: ndarray[tuple[Any, ...], dtype[float64]] | ConvexSet) → bool[source]#

Determine if the given point or convex set is a subset of this polyhedron or not.

Note

  • The H-representation of the self polyhedron is required to test if one set is contained in another. Will be computed if needed.

  • If the inner set is unbounded and not a polyhedron, then this will likely raise an error.

Math

For two homogenized polyhedra \(P = \{x \mid Ax\leq 0, A_e x=0 \}\) and \(Q=\{x\mid Bx\leq 0, B_e x=0 \}\), to test if \(Q\subseteq P\):

  • If \(Q\) has V-rep, just test that all vertices and rays of \(Q\) are contained in \(P\).

  • If \(Q\) only has H-rep, use Farkas’ lemma: \(Q\subseteq P\) if and only if the following LP is feasible:

\[\begin{split}\begin{align*} \min_{\lambda, \mu, \mu_e} &\quad 0 \\ \text{s.t.} &\quad B = \lambda A + \mu A_e \\ &\quad B_e = \mu_e A_e \end{align*}\end{split}\]
property dim: int#

Return the polyhedron’s dimension of representation.

See also

classmethod empty_set(dim: int) → Polyhedron[source]#

Create an empty polyhedron of given dimension which has a H-rep \(P = \{x \mid \mathbf{0}^{\top} x \leq -1\}\) and a V-rep (no vertices or rays).

fplot(ax: Axes | Plotter, func_name: str | None = None, **kwargs)[source]#

Plot a single function over a polyhedron.

Note

Require H-representation. Will be computed if necessary.

classmethod from_Hrep(*args, **kwargs) → Polyhedron[source]#

Create a polyhedron from its H-representation.

This method initializes a polyhedron object using the H-representation data provided as arguments, which is encapsulated in an HData object.

Parameters:

  • *args – Positional arguments to initialize the H-representation.

  • **kwargs – Keyword arguments to initialize the H-representation.

Returns:

An instance of the Polyhedron class initialized with the given H-representation data.

See also

  • from_Vrep() - Create a polyhedron instance from its V-representation.

classmethod from_Vrep(*args, **kwargs) → Polyhedron[source]#

Create a polyhedron from its V-representation.

This method initializes a polyhedron object using the V-representation data provided as arguments, which is encapsulated in an VData object.

Parameters:

  • *args – Positional arguments to initialize the V-representation.

  • **kwargs – Keyword arguments to initialize the V-representation.

Returns:

An instance of the Polyhedron class initialized with the given V-representation data.

See also

  • from_Hrep() - Create a polyhedron instance from its H-representation.

classmethod from_bounds(lb: ndarray[tuple[Any, ...], dtype[float64]], ub: ndarray[tuple[Any, ...], dtype[float64]]) → Polyhedron[source]#

Create a polyhedron from upper and lower bounds.

If the upper bounds contain numpy.inf or the lower bounds contain -numpy.inf, the corresponding dimensions are treated as unbounded.

The resulting polyhedron has both H-representation and V-representation. The H-rep is \(P = \{x \mid x_{lb} \leq x \leq x_{ub}\}\) and the minimal V-rep is computed accordingly.

get_facet(index: int | list[int] | None = None) → Polyhedron | list[Polyhedron][source]#

Extract the facets of the polyhedron specified by the inequality indices. Each facet is a polyhedron formed by rewriting the indexed ineqality as an equality constraint.

For a polyhedron \(P\) in H-representation:

\[P = \{ x\in\mathbb{R}^n \mid a_i^\top x \leq b_i, i \in \{1, ..., m\}, A_e x = b_e \}\]

given an index \(k\), the facet polyhedron \(Q\) is built as:

\[Q = \{ x\in\mathbb{R}^n \mid a_i^\top x \leq b_i, i \in \{1, ..., m\} \backslash \{k\}, a_k^\top x = b_k, A_e x = b_e \}\]
Parameters:

index (int) – an integer or list of integers. If not provided, a list of all facets will be returned.

Note

The polyhedron must have minimal (irredundant) H-representation, otherwise an error will be thrown. Consider calling minHrep() beforehand.

property hasHrep: bool#

Tell if the polyhedron has a H-representation or not.

property hasVrep: bool#

Tell if the polyhedron has a V-representation or not.

hausdorff_distance(other: Polyhedron) → float[source]#

Compute the Hausdorff distance between two polyhedra.

For two polyhedra \(P\) and \(Q\), the Hausdorff distance is defined as:

\[d_H(P, Q) = \max\left\{ \max\limits_{x\in P} \min\limits_{y\in Q} \|x-y\|, \max\limits_{y\in Q} \min\limits_{x\in P} \|x-y\| \right\}.\]

Note

This method is computationally expensive, as it requires to compute the distances of all vertices to another set for both polyhedra.

homogenize() → Polyhedron[source]#

Homogenize the Polyhedron into a cone.

  • For an H-representation \(P = \{ x\in\mathbb{R}^n \mid A x \leq b, A_e x = b_e \}\), its homogenization is:

\[\bar{P} = \{ y\in\mathbb{R}^{n+1} \mid [A ~-b] \leq 0, [A_e ~-b_e] x = 0 \}.\]

  • For a V-representation: \(P = \{ x\in\mathbb{R}^n \mid x = V^\top v + R^\top r, r \geq 0, v\geq 0, 1^\top v = 1 \}\), its homogenization is:

\[\begin{split}P_h = \{ y\in\mathbb{R}^{n+1} \mid y=\begin{bmatrix} V^\top & R^\top \\ 1 & 0\end{bmatrix} l, l \geq 0 \}.\end{split}\]
property incidence_map: dict#

Compute the incidence map of the polyhedron.

The incidence map describes containment of vertices in facets.

Note

This function computes both irredundant H-rep and V-rep and can be time consuming.

Warning

This function is not implemented yet. Calling it will raise NotImplementedError.

interior_point(facet_index: int | None = None) → Tuple[ndarray[tuple[Any, ...], dtype[float64]], bool][source]#

Compute an interior point of the polyhedron.

Parameters:

facet_index (int) – the index of the facet to compute the interior point for. If None, compute the interior point of the whole polyhedron.

Returns:

An interior point of the polyhedron and a boolean indicating whether the point is strictly inside the polyhedron.

intersect(other: Polyhedron) → Polyhedron[source]#

Compute the intersection of this polyhedron with another Polyhedron.

inv_affine_map(T: ndarray[tuple[Any, ...], dtype[float64]], t: ndarray[tuple[Any, ...], dtype[float64]] | None = None) → Polyhedron[source]#

Computes an inverse affine map \(Q\) of polyhedron \(P\) based on the transformation matrix \(T\) and vector \(t\). The polyhedron \(Q\) is given by:

\[Q = \{ x \in \mathbb{R}^n \mid T x + t \in P \subseteq \mathbb{R}^n \}\]

The matrix \(T\) must be a square real matrix. The vector \(t\), if omitted, defaults to a zero vector of corresponding dimension.

See also

is_adjacent(other: Polyhedron) → bool[source]#

Test if a polyhedron shares a facet with another polyhedron.

Return true if the polyhedron P has a facet-to-facet property with the polyhedron Q. Both polyhedra must be in H-representation. If they are not, the irredundant H-representations will be computed.

The function tests if polyhedra P and Q are adjacent by checking if their intersection is of dimension \(d-1\) and if the intersection is a facet of both polyhedra.

See also

property is_bounded#

Returns True if the polyhedron is bounded and False otherwise.

property is_empty: bool#

Returns True if the polyhedron is empty and False otherwise.

property is_full_dim: bool#

Determine if the polyhedron is full dimensional or not.

A polyhedron \(P\in\mathbb{R}^n\) is full dimensional if and only if the dimension of its affine dimension is the same as \(n\) (its dimension of representation), or equivalently if there exists a non-empty ball of dimension \(n\) that is contained within it.

is_neighbor(other: Polyhedron) → Tuple[bool, int | None, int | None][source]#

Test if a polyhedron touches another polyhedron along a given facet.

Return True if the polyhedron \(P\) shares with another polyhedron \(Q\) a part of a common facet. Both polyhedra must be in H-representation. If they are not, the irredundant H-representations will be computed.

The function tests if the intersection of two polyhedra \(P\) and S in dimension d is nonempty and is of dimension \(d-1\). If this holds, then the two polyhedra are touching themself along one facet.

Note

If one of the polyhedra is empty, the function will return False.

See also

lift(num_dims: int) → Polyhedron[source]#

Lift the polyhedron \(P\in\mathbb{R}^n\) to a higher dimension in \(\mathbb{R}^{n+m}\):

\[\{ (x, y) \mid x \in P, y \in \mathbb{R}^m \}.\]
minHrep(inplace: bool = True, use_cdd: bool = False) → HData[source]#

Compute the minimal H-representation of the polyhedron.

Parameters:

  • inplace – If True, the current object is modified.

  • use_cdd – If True, the pycddlib is used for the redundancy elimination. Otherwise the built-in implementation is used.

Returns:

The minimal H-representation of the polyhedron.

Note

  • If an H-representation is already known, then this function performs redundancy elimination. Otherwise, facet enumeration will be performed first.

  • If the polyhedron is empty, no redundancy elimination will be performed.

Math

The redundancy elimination of H-representation performs the following steps.

First, attempt to identify some irredundant inequalities using ray shooting heuristics if the polyhedron is full-dimensional. Draw a number of rays from a strict interior point along random directions. For each ray, the inequality it intersects first must be irredundant.

Second, identify redundant inequalities using Linear Programming. For each inequality \(a_j x \leq b_j\), solve the following LP:

\[\begin{split}\begin{align*} \max_x &\quad a_j x \\ \text{s.t.} &\quad a_i x \leq b_i, \quad \forall i\neq j \\ &\quad a_j x \leq b_j + 1, \\ &\quad A_e x = b_e. \\ \end{align*}\end{split}\]

This LP is guaranteed to be bounded. The inequality \(a_j x \leq b_j\) is irredundant if the LP’s optimal objective value is strictly greater than \(b_j\) and redundant otherwise.

minVrep(inplace: bool = True, use_cdd: bool = False) → VData[source]#

Returns the minimal V-representation of the polyhedron.

If inplace is True, then the current object is modified.

Math

The redundancy elimination of V-representation performs the following steps:

1. Eliminate redundant rays. An extreme ray \(r_j\) is redundant if it can be written as the non-negative linear combination of other rays \(r_i\), i.e., if the following LP is feasible:

\[\begin{split}\begin{align*} \min_\mu &\quad 0 \\ \text{s.t.} &\quad r_j = \sum\limits_{i\neq j} \mu_i r_i \\ &\quad \mu_i \geq 0 \end{align*}\end{split}\]

2. Eliminate redundant vertices. The method we use is an extension of Section 2.19 in Fukuda’s PolyFAQ. A vertex \(v_j\) is redundant if the following LP is infeasible:

\[\begin{split}\begin{align*} \min\limits_{z, z_0} &\quad 0 \\ \text{s.t.} &\quad z^\top v_j > z_0 \\ &\quad z^\top v_i \leq z_0 \\ &\quad z^\top r_k < 0 \end{align*}\end{split}\]
minkowski_sum(other: ndarray[tuple[Any, ...], dtype[float64]] | Polyhedron) → Polyhedron[source]#

Compute the Minkowski sum of the polyhedron with a point or another polyhedron.

  • The Minkowski sum of a polyhedron \(P\subseteq\mathbb{R}^n\) and a point \(z\in\mathbb{R}^n\) is defined as:

\[P \oplus z = \{ x+z \mid x\in P \}.\]

  • The Minkowski sum of two polyhedra \(P\subseteq\mathbb{R}^n\) and \(Q\subseteq\mathbb{R}^n\) is defined as:

\[P \oplus Q = \{ x+y \mid x\in P, y\in Q \}.\]
plot(ax: Axes | Plotter, **kwargs)[source]#

Plot the Polyhedron.

Note

Require V-representation. Will be computed if necessary.

pontryagin_difference(other: ndarray[tuple[Any, ...], dtype[float64]] | Polyhedron) → Polyhedron[source]#

Compute the Pontryagin difference between a polyhedron with a point or another polyhedron.

  • The Pontryagin difference of a polyhedron \(P\subseteq\mathbb{R}^n\) with a point \(z\in\mathbb{R}^n\) is defined as:

\[P \ominus z = \{ x-z \mid x\in P \} = P \oplus (-z).\]

  • The Pontryagin difference of a polyhedron \(P\subseteq\mathbb{R}^n\) with another polyhedron \(Q\subseteq\mathbb{R}^n\) is defined as:

\[P \ominus Q = \{ x\in P \mid \forall y\in Q, x+y \in P \}.\]
projection(dims: Sequence[int], method: str | None = None) → Polyhedron[source]#

Compute the orthogonal projection of the polyhedron on given dimensions.

sample(n_samples: int, strictly_interior: bool = False, method: Literal['hit_and_run', 'reject'] = 'reject', *, rng: Generator | None = None, batch_size: int = 1024, max_attempts: int | None = None, whitening: bool = True, verbose: bool = False) → ndarray[source]#

Uniformuly sample random points inside the polyhedron. Requires the polyhedron to be bounded and non-empty.

classmethod singleton(point: ndarray[tuple[Any, ...], dtype[float64]]) → Polyhedron[source]#

Create a polyhedron that contains only a single point. The H-rep is \(P = \{x \mid x = v\}\) where \(v\) is the given point. The V-rep contains the point as the only vertex and no rays.

property volume: float#

Compute the volume of this polyhedron.

Note

The volume of a polyhedron is defined as 0 if it is lower-dimensional or empty and \(\infty\) if it is full-dimensional and unbounded.

class mpt4py.Simplex(k: int)[source]#

Bases: Polyhedron

Represent a standard \(k\)-simplex.

property is_bounded#

Returns True if the polyhedron is bounded and False otherwise.

exception mpt4py.UnboundedError[source]#

Bases: Exception

Raised when a set is unbounded.

mpt4py.get_tolerance(name: str)[source]#

Retrieve a tolerance value.

mpt4py.new_figure(backend: str | None = 'pyvista')[source]#

Create a new figure using the specified backend.

mpt4py.print_tolerances()[source]#

Print all current tolerance values.

mpt4py.set_tolerance(name: str, value: float)[source]#

Modify a tolerance at runtime.

Contents