Cone Class
This class handles all computations relating to rational polyhedral cones, such cone duality and extremal ray computations. It is mainly used for the study of Kähler and Mori cones.
This class is primarily tailored to pointed (i.e. strongly convex) cones. There are a few computations, such as finding extremal rays, that may produce some unexpected results when working with non-pointed cones. Additionally, cones that are not pointed, and whose dual is also not pointed, are not supported since they are uncommon and difficult to deal with.
Constructor
cytools.cone.Cone
Description:
Constructs a Cone
object. This is handled by the hidden
__init__
function.
Arguments:
rays
(array_like, optional): A list of rays that generates the cone. If it is not specified then the hyperplane normals must be specified.hyperplanes
(array_like, optional): A list of inward-pointing hyperplane normals that define the cone. If it is not specified then the generating rays must be specified.check
(bool, optional, default=True): Whether to check the input. Recommended if constructing a cone directly.
Exactly one of rays
or hyperplanes
must be specified. Otherwise an
exception is raised.
Example
We construct a cone in two different ways. First from a list of rays then from a list of hyperplane normals. We verify that the two inputs result in the same cone.
from cytools import Cone
c1 = Cone([[0,1],[1,1]]) # Create a cone using rays. It can also be done with Cone(rays=[[0,1],[1,1]])
c2 = Cone(hyperplanes=[[1,0],[-1,1]]) # Create a cone using hyperplane normals.
c1 == c2 # We verify that the two cones are the same.
# True
Functions
ambient_dimension
Description: Returns the dimension of the ambient lattice.
Arguments: None.
Returns: (int) The dimension of the ambient lattice.
Aliases:
ambient_dim
.
Example
We construct a cone and find the dimension of the ambient lattice.
c = Cone([[0,1,0],[1,1,0]])
c.ambient_dimension()
# 3
clear_cache
Description: Clears the cached results of any previous computation.
Arguments: None.
Returns: Nothing.
Example
We construct a cone, compute its extremal rays, clear the cache and then compute them again.
c = Cone([[1,0],[1,1],[0,1]])
c.extremal_rays()
# array([[0, 1],
# [1, 0]])
c.clear_cache() # Clears the cached result
c.extremal_rays() # The extremal rays recomputed
# array([[0, 1],
# [1, 0]])
contains
Description: Checks if a point is in the (strict) interior.
Arguments:
other
: The object to check containment of. Can be a 1D array, which is treated as a point. Can be a 2D array, which is treated as a list of points. Can be a Cone.eps
: Check H@pt >= eps.
Returns: Whether pt is in the (strict) interior.
dimension
Description: Returns the dimension of the cone.
Arguments: None.
Returns: (int) The dimension of the cone.
Aliases:
dim
.
Example
We construct a cone and find its dimension.
c = Cone([[0,1,0],[1,1,0]])
c.dimension()
# 2
dual_cone
Description: Returns the dual cone.
Arguments: None.
Returns: (Cone) The dual cone.
Aliases:
dual
.
Example
We construct a cone and find its dual cone.
c = Cone([[0,1],[1,1]])
c.dual_cone()
# A rational polyhedral cone in RR^2 defined by 2 hyperplanes normals
c.dual_cone().rays()
# array([[ 1, 0],
# [-1, 1]])
extremal_hyperplanes
Description: Returns the extremal hyperplanes of the cone.
Arguments:
tol
(float, optional, default=1e-4): Specifies the tolerance for deciding whether a ray is extremal or not.- verbose (bool, optional, default=False): When set to True it show the progress while finding the extremal rays.
Returns: (numpy.ndarray) The list of extremal hyperplanes of the cone.
extremal_rays
Description: Returns the extremal rays of the cone.
By default, this function will use as many CPU threads as there are
available. To fix the number of threads, you can set the n_threads
variable in the config
submodule.
Arguments:
tol
(float, optional, default=1e-4): Specifies the tolerance for deciding whether a ray is extremal or not.- verbose (bool, optional, default=False): When set to True it show the progress while finding the extremal rays.
Returns: (numpy.ndarray) The list of extremal rays of the cone.
Example
We construct a cone and find its extremal_rays.
c = Cone([[0,1],[1,1],[1,0]])
c.extremal_rays()
# array([[0, 1],
# [1, 0]])
find_grading_vector
Description: Finds a grading vector for the cone, i.e. a vector such that any non-zero point in the cone has a positive dot product . Thus, the grading vector must be strictly interior to the dual cone, so it is only defined for pointed cones. This function returns an integer grading vector.
Arguments:
backend
(str, optional, default=None): String that specifies the optimizer to use. The options are the same as for thefind_interior_point
function.
Returns: (numpy.ndarray) A grading vector. If it could not be found then None is returned.
Example
We construct a cone and find a grading vector.
c = Cone([[3,2],[5,3]])
c.find_grading_vector()
# array([-1, 2])
find_interior_point
Description: Finds a point in the strict interior of the cone. If no point is found then None is returned.
Arguments:
c
(float, optional, default=1): A real positive number specifying the stretching of the cone (i.e. the minimum distance to the defining hyperplanes).integral
(bool, optional, default=False): A flag that specifies whether the point should have integral coordinates.backend
(str, optional, default=None): String that specifies the optimizer to use. Options are "glop", "scip", "cpsat", "mosek", "osqp", and "cvxopt". If it is not specified then "glop" is used by default. For it uses "mosek" if it is activated.show_hints
: Whether to show hints about odd backend behavior.
Returns: (numpy.ndarray) A point in the strict interior of the cone. If no point is found then None is returned.
Example
We construct a cone and find some interior points.
c = Cone([[3,2],[5,3]])
c.find_interior_point()
# array([4. , 2.5])
c.find_interior_point(integral=True)
# array([8, 5])
find_lattice_points
Description: Finds lattice points in the cone. The points are found in the region bounded by the cone, and by a cutoff surface given by the grading vector. Note that this requires the cone to be pointed. The minimum number of points to find can be specified, or if working with a preferred grading vector it is possible to specify the maximum degree.
Arguments:
min_point
(int, optional): Specifies the minimum number of points to find. The degree will be increased until this minimum number is achieved.max_deg
(int, optional): The maximum degree of the points to find. This is useful when working with a preferred grading.grading_vector
(array_like, optional): The grading vector that will be used. If it is not specified then it is computed.max_coord
(int, optional, default=1000): The maximum magnitude of the coordinates of the points.deg_window
(int, optional): If using min_points, search for lattice points with degrees in range [n(deg_window+1), n(deg_window+1)+deg_window] for 0<=nfilter_function
(function, optional): A function to use as a filter of the points that will be kept. It should return a boolean indicating whether to keep the point. Note thatmin_points
does not take the filtering into account.process_function
(function, optional): A function to process the points as they are found. This is useful to avoid first constructing a large list of points and then processing it.verbose
(boolean, optional): Whether to print extra diagnostic information (True) or not (False).
Returns: (numpy.ndarray) The list of points.
Example
We construct a cone and find at least 20 lattice points in it.
c = Cone([[3,2],[5,3]])
pts = c.find_lattice_points(min_points=20)
print(len(pts)) # We see that it found 21 points
# 21
Let's also give an example where we use a function to apply some filtering. This can be something very complicated, but here we just pick the points where all coordinates are odd.
def filter_function(pt):
return all(c%2 for c in pt)
c = Cone([[3,2],[5,3]])
pts = c.find_lattice_points(min_points=20, filter_function=filter_function)
print(len(pts)) # Now we get only 6 points instead of 21
# 6
Finally, let's give an example where we process the data as it comes instead of first constructing a list. In this simple example we just print each point with odd coordinates, but in general it can be a complex algorithm.
def process_function(pt):
if all(c%2 for c in pt):
print(f"Processing point {pt}")
c = Cone([[3,2],[5,3]])
c.find_lattice_points(min_points=20, process_function=process_function)
# Processing point (5, 3)
# Processing point (11, 7)
# Processing point (15, 9)
# Processing point (17, 11)
# Processing point (21, 13)
# Processing point (25, 15)
hilbert_basis
Description: Returns the Hilbert basis of the cone. Normaliz is used for the computation.
Arguments: None.
Returns: (numpy.ndarray) The list of vectors forming the Hilbert basis.
Example
We compute the Hilbert basis of a two-dimensional cone.
c = Cone([[1,3],[2,1]])
c.hilbert_basis()
# array([[1, 1],
# [1, 2],
# [1, 3],
# [2, 1]])
hyperplanes
Description: Returns the inward-pointing normals to the hyperplanes that define the cone.
Arguments: None.
Returns: (numpy.ndarray) The list of inward-pointing normals to the hyperplanes that define the cone.
Example
We construct two cones and find their hyperplane normals.
c1 = Cone([[0,1],[1,1]])
c2 = Cone(hyperplanes=[[0,1],[1,1]])
c1.hyperplanes()
# array([[ 1, 0],
# [-1, 1]])
c2.hyperplanes()
# array([[0, 1],
# [1, 1]])
intersection
Description: Computes the intersection with another cone, or with a list of cones.
Arguments:
other
(Cone or array_like): The other cone that is being intersected, or a list of cones to intersect with.
Returns: (Cone) The cone that results from the intersection.
Example
We construct two cones and find their intersection.
c1 = Cone([[1,0],[1,2]])
c2 = Cone([[0,1],[2,1]])
c3 = c1.intersection(c2)
c3.rays()
# array([[2, 1],
# [1, 2]])
is_pointed
Description: Returns True if the cone is pointed (i.e. strongly convex).
There are two available methods to perform the computation. Using NNLS it directly checks if it can find a linear subspace. Alternatively, it can check if the dual cone is solid. For extremely wide cones the second approach is more reliable, so that is the default one.
Arguments:
backend
(str, optional): Specifies which backend to use. Available options are "nnls", and any backends available for theis_solid
function. If not specified, it uses the default backend for theis_solid
function.tol
(float, optional, default=1e-7): The tolerance for determining when a linear subspace is found. This is only used for the NNLS backend.
Returns: (bool) The truth value of the cone being pointed.
Aliases:
is_strongly_convex
.
Example
We construct two cones and check if they are pointed.
c1 = Cone([[1,0],[0,1]])
c2 = Cone([[1,0],[0,1],[-1,0]])
c1.is_pointed()
# True
c2.is_pointed()
# False
is_simplicial
Description: Returns True if the cone is simplicial.
Arguments: None.
Returns: (bool) The truth value of the cone being simplicial.
Example
We construct two cones and check if they are simplicial.
c1 = Cone([[1,0,0],[0,1,0],[0,0,1]])
c2 = Cone([[1,0,0],[0,1,0],[0,0,1],[1,1,-1]])
c1.is_simplicial()
# True
c2.is_simplicial()
# False
is_smooth
Description: Returns True if the cone is smooth, i.e. its extremal rays either form a basis of the ambient lattice, or they can be extended into one.
Arguments: None.
Returns: (bool) The truth value of the cone being smooth.
Example
We construct two cones and check if they are smooth.
c1 = Cone([[1,0,0],[0,1,0],[0,0,1]])
c2 = Cone([[2,0,1],[0,1,0],[1,0,2]])
c1.is_smooth()
# True
c2.is_smooth()
# False
is_solid
Description: Returns True if the cone is solid, i.e. if it is full-dimensional.
If the generating rays are known then this can simply be checked by computing the dimension of the linear space that they span. However, when only the hyperplane inequalities are known this can be a difficult problem. When using PPL as the backend, the convex hull is explicitly constructed and checked. The other backends try to find a point in the strict interior of the cone, which fails if the cone is not solid. The latter approach is much faster, but there could be extremely narrow cones where the optimization fails and this function returns a false negative.
Arguments:
backend
(str, optional): Specifies which backend to use. Available options are "ppl", and any backends available for thefind_interior_point
function. If not specified, it uses the default backend of thefind_interior_point
function.
Returns: (bool) The truth value of the cone being solid.
Aliases:
is_full_dimensional
.
Example
We construct two cones and check if they are solid.
c1 = Cone([[1,0],[0,1]])
c2 = Cone([[1,0,0],[0,1,0]])
c1.is_solid()
# True
c2.is_solid()
# False
rays
Description: Returns the (not necessarily extremal) rays that generate the cone.
Arguments:
use_extremal_hyperplanes
: Whether to use extremal hyperplanes in this computation, or just any hyperplanes.verbosity
: The verbosity level.
Returns: (numpy.ndarray) The list of rays that generate the cone.
Example
We construct two cones and find their generating rays.
c1 = Cone([[0,1],[1,1]])
c2 = Cone(hyperplanes=[[0,1],[1,1]])
c1.rays()
# array([[0, 1],
# [1, 1]])
c2.rays()
# array([[ 1, 0],
# [-1, 1]])
tip_of_stretched_cone
Description:
Finds the tip of the stretched cone. The stretched cone is defined as
the convex polyhedral region inside the cone that is at least a
distance c
from any of its defining hyperplanes. Its tip is defined
as the point in this region with the smallest norm.
This is a problem that requires quadratic programming since the norm of a vector is being minimized. For dimensions up to around 50, this can easily be done with open-source solvers like OSQP or CVXOPT, however for higher dimensions this becomes a difficult task that only the Mosek optimizer is able to handle. However, Mosek is closed-source and requires a license. For this reason we preferentially use ORTools, which is open-source, to solve a linear problem and find a good approximation of the tip. Nevertheless, if Mosek is activated then it uses Mosek as it is faster and more accurate.
Arguments:
c
(float): A real positive number specifying the stretching of the cone (i.e. the minimum distance to the defining hyperplanes).backend
(str, optional, default=None): String that specifies the optimizer to use. Options are "mosek", "osqp", "cvxopt", and "glop". If it is not specified then for it uses "osqp" by default. For it uses "mosek" if it is activated, or "glop" otherwise.check
(bool, optional, default=True): Flag that specifies whether to check if the output of the optimizer is consistent and satisfiesconstraint_error_tol
.constraint_error_tol
(float, optional, default=1e-2): Error tolerance for the linear constraints.max_iter
(int, optional, default=10**6): The maximum number of iterations allowed for the non-GLOP backends. If this function is returning None, then increasing this parameter (maximum permissible value: 2**31-1) might resolve the issue. For backend=="glop", this does nothing.show_hints
: Whether to show hints about odd backend behavior.verbose
(boolean, optional): Whether to print extra diagnostic information (True) or not (False).
Returns: (numpy.ndarray) The vector specifying the location of the tip. If it could not be found then None is returned.
Example
We construct two cones and find the locations of the tips of the stretched cones.
c1 = Cone([[1,0],[0,1]])
c2 = Cone([[3,2],[5,3]])
c1.tip_of_stretched_cone(1)
# array([1., 1.])
c2.tip_of_stretched_cone(1)
# array([8., 5.])
Hidden Functions
__eq__
Description: Implements comparison of cones with ==.
The comparison of cones that are not pointed, and whose duals are also not pointed, is not supported.
Arguments:
other
(Cone): The other cone that is being compared.
Returns: (bool) The truth value of the cones being equal.
Example
We construct two cones and compare them.
c1 = Cone([[0,1],[1,1]])
c2 = Cone(hyperplanes=[[1,0],[-1,1]])
c1 == c2
# True
__hash__
Description: Implements the ability to obtain hash values from cones.
Cones that are not pointed, and whose duals are also not pointed, are not supported.
Arguments: None.
Returns: (int) The hash value of the cone.
Example
We compute the hash value of a cone. Also, we construct a set and a dictionary with a cone, which make use of the hash function.
c = Cone([[0,1],[1,1]])
h = hash(c) # Obtain hash value
d = {c: 1} # Create dictionary with cone keys
s = {c} # Create a set of cones
__init__
Description:
Initializes a Cone
object.
Arguments:
rays
(array_like, optional): A list of rays that generates the cone. If it is not specified then the hyperplane normals must be specified.hyperplanes
(array_like, optional): A list of inward-pointing hyperplane normals that define the cone. If it is not specified then the generating rays must be specified.check
(bool, optional, default=True): Whether to check the input. Recommended if constructing a cone directly.copy
(bool, optional, default=True): Whether to ensure we copy the input rays/hyperplanes. Recommended.
Exactly one of rays
or hyperplanes
must be specified. Otherwise, an
exception is raised.
Returns: Nothing.
Example
This is the function that is called when creating a new Cone
object.
We construct a cone in two different ways. First from a list of rays
then from a list of hyperplane normals. We verify that the two inputs
result in the same cone.
from cytools import Cone
c1 = Cone([[0,1],[1,1]]) # Create a cone using rays. It can also be done with Cone(rays=[[0,1],[1,1]])
c2 = Cone(hyperplanes=[[1,0],[-1,1]]) # Create a cone using hyperplane normals.
c1 == c2 # We verify that the two cones are the same.
# True
__ne__
Description: Implements comparison of cones with !=.
The comparison of cones that are not pointed, and whose duals are also not pointed, is not supported.
Arguments:
other
(Cone): The other cone that is being compared.
Returns: (bool) The truth value of the cones being different.
Example
We construct two cones and compare them.
c1 = Cone([[0,1],[1,1]])
c2 = Cone(hyperplanes=[[1,0],[-1,1]])
c1 != c2
# False
__repr__
Description: Returns a string describing the polytope.
Arguments: None.
Returns: (str) A string describing the polytope.
Example
This function can be used to convert the Cone to a string or to print information about the cone.
c = Cone([[1,0],[1,1],[0,1]])
cone_info = str(c) # Converts to string
print(c) # Prints cone info
# A 2-dimensional rational polyhedral cone in RR^2 generated by 3 rays