Skip to main content

Miscellaneous Functions

There are various other functions in CYTools that don't belong to any particular class. They are defined in different places according to where they most closely belong. Here we list the location of the definitions of these functions as well as their documentation.

Functions in cytools.polytope

is_reflexive_barebones

Description: Minimal code to check if conv(points) is reflexive.

Arguments:

  • points: The points defining the hull.
  • backend: The backend to use. See poly_v_to_h.

Returns: Whether conv(points) is reflexive


poly_v_to_h

Description: Generate the H-representation of a polytope, given the V-representation. I.e., map points/vertices to hyperplanes inequalities.

The o inequalities are in the form c0 * x_0 + ... + c{d-1} * x_{d-1} + c_d >= 0

Arguments:

  • pts: The input points. Each row is a point.
  • backend: The backend to use. Currently, support "ppl", "qhull", and "palp".

Returns: The hyperplane inequalities in the form c0 * x_0 + ... + c{d-1} * x_{d-1} + c_d >= 0 and, depending on backend/dimension, the formal convex hull of the points.


saturating_lattice_pts

Description: Computes the lattice points contained in conv(pts), along with the indices of the hyperplane inequalities that they saturate.

Arguments:

  • pts: A list of points spanning the hull.
  • ineq: Hyperplane inqualities defining the hull. Same format as output by poly_v_to_h
  • dim: The dimension of the hull.
  • backend: The backend to use. Either "palp" or defaults to native.

Returns: An array of all lattice points (the rows). A list of sets of all inequalities each lattice point saturates.


Functions in cytools.triangulation

_cgal_triangulate

Description: Computes a regular triangulation using CGAL.

note

This function is not intended to be called by the end user. Instead, it is used by the Triangulation class when using CGAL as the backend.

Arguments:

  • points: A list of points.
  • heights: A list of heights defining the regular triangulation.

Returns: A list of simplices defining a regular triangulation.

Example

This function is not intended to be directly used, but it is used in the following example. We construct a triangulation using CGAL.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-6,-9]])
p.triangulate(backend="cgal")
# A fine, star triangulation of a 4-dimensional point configuration with 7
# points in ZZ^4

_qhull_triangulate

Description: Computes a regular triangulation using QHull.

note

This function is not intended to be called by the end user. Instead, it is used by the Triangulation class when using QHull as the backend.

Arguments:

  • points: A list of points.
  • heights: A list of heights defining the regular triangulation.

Returns: A list of simplices defining a regular triangulation.

Example

This function is not intended to be directly used, but it is used in the following example. We construct a triangulation using QHull.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-6,-9]])
p.triangulate(backend="qhull")
# A fine, star triangulation of a 4-dimensional point configuration with 7
# points in ZZ^4

_to_star

Description: Turns a triangulation into a star triangulation by deleting internal lines and connecting all points to the origin.

note

This function is not intended to be called by the end user. Instead, it is used by the Triangulation class when needed.

info

This function is only reliable for triangulations of reflexive polytopes and may produce invalid triangulations for other polytopes.

Arguments:

  • triang: The triangulation to convert to star.

Returns: Nothing.


_topcom_triangulate

Description: Computes the placing/pushing triangulation using TOPCOM.

note

This function is not intended to be called by the end user. Instead, it is used by the Triangulation class when using TOPCOM as the backend.

Arguments:

  • points: A list of points.

Returns: A list of simplices defining a triangulation.

Example

This function is not intended to be directly used, but it is used in the following example. We construct a triangulation using TOMCOM.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-6,-9]])
p.triangulate(backend="topcom")
# A fine, star triangulation of a 4-dimensional point configuration with 7
# points in ZZ^4

all_triangulations

Description: Computes all triangulations of the input point configuration using TOPCOM. There is the option to only compute fine, regular or fine triangulations.

note

This function is not intended to be called by the end user. Instead, it is used by the all_triangulations function of the Polytope class.

Arguments:

  • poly: The ambient polytope.
  • pts: The list of points to be triangulated. Specified by labels.
  • only_fine: Restricts to only fine triangulations.
  • only_regular: Restricts to only regular triangulations.
  • only_star: Restricts to only star triangulations.
  • star_origin: The index of the point that will be used as the star origin. If the polytope is reflexive this is set to 0, but otherwise it must be specified.
  • backend: The optimizer used to check regularity computation. The available options are "topcom" and the backends of the is_solid function of the Cone class. If not specified, it will be picked automatically. Note that using TOPCOM to check regularity is slower.
  • raw_output: Return the triangulations as lists of simplices instead of as Triangulation objects.

Returns: A generator of Triangulation objects with the specified properties. If raw_output is set to True then numpy arrays are used instead.

Example

This function is not intended to be directly used, but it is used in the following example. We construct a polytope and find all of its triangulations. We try picking different restrictions and see how the number of triangulations changes.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-2,-1,-1],[-2,-1,-1,-1]])
g = p.all_triangulations()
next(g) # Takes some time while TOPCOM finds all the triangulations
# A fine, regular, star triangulation of a 4-dimensional point
# configuration with 7 points in ZZ^4
next(g) # Produces the next triangulation immediately
# A fine, regular, star triangulation of a 4-dimensional point
# configuration with 7 points in ZZ^4
len(p.all_triangulations(as_list=True)) # Number of fine, regular, star triangulations
# 2
len(p.all_triangulations(only_regular=False, only_star=False, only_fine=False, as_list=True) )# Number of triangularions, no matter if fine, regular, or star
# 6

random_triangulations_fair_generator

Description: Constructs pseudorandom regular (optionally star) triangulations of a given point set. Implements Algorithm #3 from the paper Bounding the Kreuzer-Skarke Landscape by Mehmet Demirtas, Liam McAllister, and Andres Rios-Tascon. arXiv:2008.01730

This is a Markov chain Monte Carlo algorithm that involves taking random walks inside the subset of the secondary fan corresponding to fine triangulations and performing random flips. For details, please see Section 4.1 in the paper.

notes
  • This function is not intended to be called by the end user. Instead, it is used by the random_triangulations_fair function of the Polytope class.
  • This function is designed mainly for large polytopes where sampling triangulations is challenging. When small polytopes are used it is likely to get stuck.

Arguments:

  • poly: The ambient polytope.
  • pts: The list of points to be triangulated. Specified by labels.
  • N: Number of desired unique triangulations. If not specified, it will generate as many triangulations as it can find until it has to retry more than max_retries times to obtain a new triangulation.
  • n_walk: Number of hit-and-run steps per triangulation.
  • n_flip: Number of random flips performed per triangulation.
  • initial_walk_steps: Number of hit-and-run steps to take before starting to record triangulations. Small values may result in a bias towards Delaunay-like triangulations.
  • walk_step_size: Determines size of random steps taken in the secondary fan. Algorithm may stall if too small.
  • max_steps_to_wall: Maximum number of steps to take towards a wall of the subset of the secondary fan that correspond to fine triangulations. If a wall is not found, a new random direction is selected. Setting this to be very large (>100) reduces performance. If this, or walk_step_size, is set to be too low, the algorithm may stall.
  • fine_tune_steps: Number of steps to determine the location of a wall. Decreasing improves performance, but might result in biased samples.
  • max_retries: Maximum number of attempts to obtain a new triangulation before the process is terminated.
  • make_star: Converts the obtained triangulations into star triangulations.
  • backend: Specifies the backend used to compute the triangulation. The available options are "cgal" and "qhull".
  • seed: A seed for the random number generator. This can be used to obtain reproducible results.

Returns: A generator of Triangulation objects with the specified properties.

Example

This function is not intended to be directly used, but it is used in the following example. We construct a polytope and find some random triangulations. The computation takes considerable time, but they should be a fair sample from the full set of triangulations (if the parameters are chosen correctly). For (some) machine learning purposes or when the fairness of the sample is not crucial, the random_triangulations_fast function should be used instead.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-1,-1]]).dual()
g = p.random_triangulations_fast()
next(g) # Takes a long time (around a minute)
# A fine, regular, star triangulation of a 4-dimensional point configuration with 106 points in ZZ^4
next(g) # Takes slightly shorter (still around a minute)
# A fine, regular, star triangulation of a 4-dimensional point
# configuration with 106 points in ZZ^4
rand_triangs = p.random_triangulations_fair(N=10, as_list=True) # Produces the list of 10 triangulations, but takes a long time (around 10 minutes)

It is worth noting that the time it takes to obtain each triangulation varies very significantly on the parameters used. The function tries to guess reasonable parameters, but it is better to adjust them to your desired balance between speed and fairness of the sampling.


random_triangulations_fast_generator

Constructs pseudorandom regular (optionally fine and star) triangulations of a given point set. This is done by picking random heights around the Delaunay heights from a Gaussian distribution.

note

This function is not intended to be called by the end user. Instead, it is used by the random_triangulations_fast function of the Polytope class.

important

This function produces random triangulations very quickly, but it does not produce a fair sample. When a fair sampling is required the random_triangulations_fair function should be used.

Arguments:

  • poly: The ambient polytope.
  • pts: The list of points to be triangulated. Specified by labels.
  • N: Number of desired unique triangulations. If not specified, it will generate as many triangulations as it can find until it has to retry more than max_retries times to obtain a new triangulation.
  • c: A contant used as the standard deviation of the Gaussian distribution used to pick the heights. A larger c results in a wider range of possible triangulations, but with a larger fraction of them being non-fine, which slows down the process when using only_fine=True.
  • max_retries: Maximum number of attempts to obtain a new triangulation before the process is terminated.
  • make_star: Converts the obtained triangulations into star triangulations.
  • only_fine: Restricts to fine triangulations.
  • backend: Specifies the backend used to compute the triangulation. The available options are "cgal" and "qhull".
  • seed: A seed for the random number generator. This can be used to obtain reproducible results.
  • verbosity: The verbosity level.

Returns: A generator of Triangulation objects with the specified properties.

Example

This function is not intended to be directly used, but it is used in the following example. We construct a polytope and find some random triangulations. The triangulations are obtained very quicly, but they are not a fair sample of the space of triangulations. For a fair sample, the random_triangulations_fair function should be used.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-1,-1]]).dual()
g = p.random_triangulations_fast()
next(g) # Runs very quickly
# A fine, regular, star triangulation of a 4-dimensional point
# configuration with 106 points in ZZ^4
next(g) # Keeps producing triangulations until it has trouble finding more
# A fine, regular, star triangulation of a 4-dimensional point
# configuration with 106 points in ZZ^4
rand_triangs = p.random_triangulations_fast(N=10, as_list=True) # Produces the list of 10 triangulations very quickly

Functions in cytools.calabiyau

_group_by_deg

Description: Organize the charges by their degrees.

Arguments:

  • charges: The charges.
  • grading_vec: The grading vector.
  • as_np_arr: Whether to map charges to np.array (True) or leave as set (False).

Returns: (dictionary) Dictionary mapping degree to charges.


def __init__


def charges


def cone


def cutoff


def grading_vec


def gvs


def gws


def invariant


def size


Functions in cytools.cone

feasibility

Description: Solve a feasibility problem Ax>=c.

Arguments:

  • hyperplanes: The constraining hyperplanes, A.
  • c: The 'stretching'.
  • ambient_dim: The ambient dimension... A.shape[1].
  • backend: The backend to use. Options are "glop", "scip", or "cpsat".
  • verbose: Whether to print extra diagnostic info.

Returns: A feasible point, if it exists. Else, None.


is_extremal

Description: Auxiliary function that is used to find the extremal rays of cones. Returns True if the ray is extremal and False otherwise. It has additional parameters that are used when parallelizing.

Arguments:

  • A (array_like): A matrix where the columns are rays (excluding b).
  • b (array_like): The vector that will be checked if it can be expressed as a positive linear combination of the columns of A.
  • i (int, optional): An id number that is used when parallelizing.
  • q (multiprocessing.Queue, optional): A queue that is used when parallelizing.
  • tol (float, optional, default=1e-4): The tolerance for determining whether a ray is extremal.

Returns: (bool or None) The truth value of the ray being extremal. If the process fails then it returns nothing.

Example

This function is not meant to be directly used by the end user. Instead it is used by the extremal_rays function. 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]])

Functions in cytools.utils

array_from_flint

Description: Converts a numpy array with fmpz/fmpq (Flint's integer/float number class) entries to 64-bit integer/float entries.

Arguments:

  • arr: A numpy array with fmpz/fmpq entries.

Returns: A numpy array with 64-bit integer/float entries.


array_to_flint

Description: Converts a numpy array with either: 1) 64-bit integer entries or 2) float entries to Flint type (fmpz or fmpq for integer or rational numbers, respectively).

See https://flintlib.org/doc/fmpz.html and https://flintlib.org/doc/fmpq.html

Arguments:

  • arr: A numpy array with either 64-bit integer or float entries.

Returns: A numpy array with either fmpz or fmpq entries.

Example

We convert an integer array to an fmpz array.

from cytools.utils import array_int_to_fmpz
arr = [[1,0,0],[0,2,0],[0,0,3]]
array_int_to_fmpz(arr)
# array([[1, 0, 0],
# [0, 2, 0],
# [0, 0, 3]], dtype=object)

fetch_polytopes

Description: Fetches reflexive polytopes from the Kreuzer-Skarke database or from the Schöller-Skarke database. The data is fetched from the websites http://hep.itp.tuwien.ac.at/~kreuzer/CY/ and http://rgc.itp.tuwien.ac.at/fourfolds/ respectively.

note

The Kreuzer-Skarke database does not store favorability data. Thus, when setting favorable to True or False it fetches additional polytopes so that after filtering by favorability it can saturate the requested limit. However, it may happen that fewer polytopes than requested are returned even though more exist. To verify that no more polytopes with the requested conditions exist one can increase the limit significantly and check if more polytopes are returned.

Arguments:

  • h11: The Hodge number h1,1h^{1,1} of the Calabi-Yau hypersurface.
  • h12: The Hodge number h1,2h^{1,2} of the Calabi-Yau hypersurface.
  • h13: The Hodge number h1,3h^{1,3} of the Calabi-Yau hypersurface.
  • h21: The Hodge number h2,1h^{2,1} of the Calabi-Yau hypersurface. This is equivalent to the h12 parameter.
  • h22: The Hodge number h2,2h^{2,2} of the Calabi-Yau hypersurface.
  • h31: The Hodge number h3,1h^{3,1} of the Calabi-Yau hypersurface. This is equivalent to the h13 parameter.
  • chi: The Euler characteristic of the Calabi-Yau hypersurface.
  • lattice: The lattice on which the polytope is defined. Options are "N" and "M". Has to be specified if the Hodge numbers or the Euler characteristic is specified.
  • dim: The dimension of the polytope. Only available options are 4 and 5.
  • n_points: The number of lattice points of the desired polytopes.
  • n_vertices: The number of vertices of the desired polytopes.
  • n_dual_points: The number of points of the dual polytopes of the desired polytopes.
  • n_facets: The number of facets of the desired polytopes.
  • limit: The maximum number of fetched polytopes.
  • timeout: The maximum number of seconds to wait for the server to return the data.
  • as_list: Return the list of polytopes instead of a generator.
  • backend: A string that specifies the backend used for the Polytope class.
  • dualize: Flag that indicates whether to dualize all the polytopes before yielding them.
  • favorable: Yield or return only polytopes that are favorable when set to True, or non-favorable when set to False. If not specified then it yields both favorable and non-favorable polytopes.

Returns: A generator of Polytope objects, or the full list when as_list is set to True.

Example

We fetch polytopes from the Kreuzer-Skarke and Schöller-Skarke databases with a few different parameters.

from cytools import fetch_polytopes # Note that it can directly be imported from the root
g = fetch_polytopes(h11=27, lattice="N") # Constructs a generator of polytopes
next(g)
# A 4-dimensional reflexive lattice polytope in ZZ^4
l = fetch_polytopes(h11=27, lattice="N", as_list=True, limit=100) # Constructs a generator of polytopes
print(f"Fetched {len(l)} polytopes")
# Fetched 100 polytopes
g_5d = fetch_polytopes(h11=1000, lattice="N", dim=5, limit=100) # Generator of 5D polytopes
next(g_5d)
# A 5-dimensional reflexive lattice polytope in ZZ^5

filter_tensor_indices

Description: Selects a specific subset of indices from a tensor.

The tensor is reindexed so that indices are in the range 0..len(indices) with the ordering specified by the input indices. This function can be used to convert the tensor of intersection numbers to a given basis.

Arguments:

  • tensor: The input symmetric sparse tensor of the form of a dictionary {(a,b,...,c):M_ab...c, ...}.
  • indices: The list of indices that will be preserved.

Returns: A dictionary describing a tensor in the same format as the input, but only with the desired indices.

Example

We construct a simple tensor and then filter some of the indices. We also give a concrete example of when this is used for intersection numbers.

from cytools.utils import filter_tensor_indices
tensor = {(0,1):0, (1,1):1, (1,2):2, (1,3):3, (2,3):4}
filter_tensor_indices(tensor, [1,3])
# {(0, 0): 1, (0, 1): 3}
p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-1,-1]])
t = p.triangulate()
v = t.get_toric_variety()
intnums_nobasis = v.intersection_numbers(in_basis=False)
intnums_inbasis = v.intersection_numbers(in_basis=True)
intnums_inbasis == filter_tensor_indices(intnums_nobasis, v.divisor_basis())
# True

find_new_affinely_independent_points

Description: Finds new points that are affinely independent to the input list of points.

This is useful when one wants to turn a polytope that is not full-dimensional into one that is, without affecting the structure of the triangulations.

Arguments:

  • pts: A list of points.

Returns: A list of affinely independent points with respect to the ones inputted.

Example

We construct a list of points and then find a set of affinely independent points.

pts = [[1,0,1],[0,0,1],[0,1,1]]
find_new_affinely_independent_points(pts)
array([[1, 0, 2]])

float_to_fmpq

Description: Converts a float to an fmpq (Flint's rational number class).

See https://flintlib.org/doc/fmpq.html

Arguments:

  • c: The input number.

Returns: The rational number that most reasonably approximates the input.

Example

We convert a few floats to rational numbers.

from cytools.utils import float_to_fmpq
float_to_fmpq(0.1), float_to_fmpq(0.333333333333), float_to_fmpq(2.45)
# (1/10, 1/3, 49/20)

fmpq_to_float

Description: Converts an fmpq (Flint's rational number class) to a float.

See https://flintlib.org/doc/fmpq.html

Arguments:

  • c: The input rational number.

Returns: The number as a float.

Example

We convert a few rational numbers to floats.

from cytools.utils import fmpq_to_float
from flint import fmpq
fmpq_to_float(fmpq(1,2)), fmpq_to_float(fmpq(1,3)),\
fmpq_to_float(fmpq(49,20))
# (0.5, 0.3333333333333333, 2.45)

gcd_float

Description: Compute the greatest common (floating-point) divisor of a and b. This is simply the largest floating point number that divides a and b. Uses the Euclidean algorithm.

Warning - unexpected/buggy behavior can occur if b starts tiny. E.g., gcd_float(100,0.1,0.2) returns 100.

This only seems to be a risk if b starts below tol.

Arguments:

  • a: The first number.
  • b: The second number.
  • tol: The tolerance for rounding.

Returns: The gcd of a and b.

Example

We compute the gcd of two floats. This function is useful since the standard gcd functions raise an error for non-integer values.

from cytools.utils import gcd_float
gcd_float(0.2, 0.5)
# Should be 0.1, but there are small rounding errors
# 0.09999999999999998

lll_reduce

Apply lll-reduction to the input points (the rows).

Arguments:

  • pts: A list of points.

Returns: The reduced points (pts_red; as rows of a numpy array). If transform==True, also return the transformation matrix/inverse (A, Ainv) s.t. pts_red.T = A*pts_in.T. As numpy arrays.


polytope_generator

Description: Reads polytopes from a file or a string. The polytopes can be specified with their vertices, as used in the Kreuzer-Skarke database, or from a weight system.

note

This function is not intended to be called by the end user. Instead, it is used by the read_polytopes and fetch_polytopes functions.

Arguments:

  • input: Specifies the name of the file to read or the string containing the polytopes.
  • input_type: Specifies whether to read from a file or from the input string. Options are "file" or "str".
  • format: Specifies the format to read. The options are "ks", which is the format used in the KS database, and "ws", if the polytopes should be constructed from weight systems.
  • backend: A string that specifies the backend used for the Polytope class.
  • dualize: Flag that indicates whether to dualize all the polytopes before yielding them.
  • favorable: Yield only polytopes that are favorable when set to True, or non-favorable when set to False. If not specified then it yields both favorable and non-favorable polytopes.
  • lattice: The lattice to use when checking favorability. This parameter is only required when favorable is set. Options are "M" and "N".
  • limit: Sets a maximum numbers of polytopes to yield.

Returns: A generator of Polytope objects.

Example

Since this function should not be used directly, we show an example of it being used with the read_polytopes function. We take a string obtained from the KS database and read the polytope it specifies.

from cytools import read_polytopes # Note - it cannot be imported from root
poly_data = '''4 5 M:10 5 N:376 5 H:272,2 [540]
1 0 0 0 -9
0 1 0 0 -6
0 0 1 0 -1
0 0 0 1 -1
'''
read_polytopes(poly_data, input_type="str", as_list=True)
# [A 4-dimensional reflexive lattice polytope in ZZ^4]

read_polytopes

Description: Reads polytopes from a file or a string. The polytopes can be specified with their vertices, as used in the Kreuzer-Skarke database, or from a weight system.

Arguments:

  • input: Specifies the name of the file to read or the string containing the polytopes.
  • input_type: Specifies whether to read from a file or from the input string. Options are "file" or "str".
  • format: Specifies the format to read. The options are "ks", which is the format used in the KS database, and "ws", if the polytopes should be constructed from weight systems.
  • backend: A string that specifies the backend used for the Polytope class.
  • as_list: Return the list of polytopes instead of a generator.
  • dualize: Flag that indicates whether to dualize all the polytopes before yielding them.
  • favorable: Yield or return only polytopes that are favorable when set to True, or non-favorable when set to False. If not specified then it yields both favorable and non-favorable polytopes.
  • lattice: The lattice to use when checking favorability. This parameter is only required when favorable is set. Options are "M" and "N".
  • limit: Sets a maximum numbers of polytopes to yield.

Returns: A generator of Polytope objects, or the full list when as_list is set to True.

Example

We take a string obtained from the KS database and read the polytope it specifies.

from cytools import read_polytopes # Note that it can directly be imported from the root
poly_data = '''4 5 M:10 5 N:376 5 H:272,2 [540]
1 0 0 0 -9
0 1 0 0 -6
0 0 1 0 -1
0 0 0 1 -1
'''
read_polytopes(poly_data, input_type="str", as_list=True)
# [A 4-dimensional reflexive lattice polytope in ZZ^4]

set_curve_basis

Description: Specifies a basis of curves of the toric variety, which in turn specifies a dual basis of divisors. This can be done with a vector specifying the indices of the dual prime toric divisors or as a matrix with the rows being the basis curves, and the entries are the intersection numbers with the prime toric divisors. Note that when using a vector it is equivalent to using the same vector in the set_divisor_basis function.

info

This function should generally not be called directly by the user. Instead, it is called by the set_curve_basis function of the ToricVariety class, or the set_curve_basis function of the CalabiYau class.

note

Only integral bases are supported by CYTools, meaning that all toric curves must be able to be written as an integral linear combination of the basis curves.

Arguments:

  • tv_or_cy: The toric variety or Calabi-Yau whose basis will be set.
  • basis: Vector or matrix specifying a basis. When a vector is used, the entries will be taken as indices of the standard basis of the dual to the lattice of prime toric divisors. When a matrix is used, the rows are taken as linear combinations of the aforementioned elements.
  • include_origin: Whether to interpret the indexing specified by the input vector as including the origin.

Returns: Nothing.

Example

This function is not intended to be directly used, but it is used in the following example. We consider a simple toric variety with two independent curves. We first find the default basis of curves it picks and then set a basis of our choice.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-6,-9]])
t = p.triangulate()
v = t.get_toric_variety()
v.curve_basis() # We haven't set any basis
# array([1, 6])
v.set_curve_basis([5,6]) # Here we set a basis
v.curve_basis() # We get the basis we set
# array([5, 6])
v.curve_basis(as_matrix=True) # We get the basis in matrix form
# array([[-18, 1, 9, 6, 1, 1, 0],
# [ -6, 0, 3, 2, 0, 0, 1]])

Note that when setting a curve basis in this way, the function behaves exactly the same as set_divisor_basis. For a more advanced example involving generic bases these two functions differ. An example can be found in the experimental features section.


set_divisor_basis

Description: Specifies a basis of divisors for the toric variety or Calabi-Yau manifold, which in turn specifies a dual basis of curves. This can be done with a vector specifying the indices of the prime toric divisors, or as a matrix where each row is a linear combination of prime toric divisors.

info

This function should generally not be called directly by the user. Instead, it is called by the set_divisor_basis function of the ToricVariety class, or the set_divisor_basis function of the CalabiYau class.

note

Only integral bases are supported by CYTools, meaning that all prime toric divisors must be able to be written as an integral linear combination of the basis divisors.

Arguments:

  • tv_or_cy: The toric variety or Calabi-Yau whose basis will be set.
  • basis: Vector or matrix specifying a basis. When a vector is used, the entries will be taken as the indices of points of the polytope or prime divisors of the toric variety. When a matrix is used, the rows are taken as linear combinations of the aforementioned divisors.
  • include_origin: Whether to interpret the indexing specified by the input vector as including the origin.

Returns: Nothing.

Example

This function is not intended to be directly used, but it is used in the following example. We consider a simple toric variety with two independent divisors. We first find the default basis it picks and then we set a basis of our choice.

p = Polytope([[1,0,0,0],[0,1,0,0],[0,0,1,0],[0,0,0,1],[-1,-1,-6,-9]])
t = p.triangulate()
v = t.get_toric_variety()
v.divisor_basis() # We haven't set any basis
# array([1, 6])
v.set_divisor_basis([5,6]) # Here we set a basis
v.divisor_basis() # We get the basis we set
# array([5, 6])
v.divisor_basis(as_matrix=True) # We get the basis in matrix form
# array([[0, 0, 0, 0, 0, 1, 0],
# [0, 0, 0, 0, 0, 0, 1]])

An example for more generic basis choices can be found in the experimental features section.


solve_linear_system

Description: Solves the sparse linear system M*x + C = 0.

Arguments:

  • M: The matrix.
  • C: The constant term.
  • backend: The solver to use. Options are "all", "sksparse" and "scipy". When set to "all" it tries all available backends.
  • check: Whether to explicitly check the solution to the linear system.
  • backend_error_tol: Error tolerance for the solution.
  • verbosity: The verbosity level.
    • verbosity = 0: Do not print anything.
    • verbosity = 1: Print warnings when backends fail.

Returns: Floating-point solution to the linear system.

Example

We solve a very simple linear equation.

from cytools.utils import to_sparse, solve_linear_system
id_array = [[0,0,1],[1,1,1],[2,2,1],[3,3,1],[4,4,1]]
M = to_sparse(id_array, sparse_type="csr")
C = [1,1,1,1,1]
solve_linear_system(M, C)
# array([-1., -1., -1., -1., -1.])

symmetric_dense_to_sparse

Description: Converts a dense symmetric tensor to a sparse tensor of the form {(a,b,...,c):M_ab...c, ...}.

The upper triangular indices are used. That is, a<=b<=...<=c.

Optionally, it applies a basis transformation.

Arguments:

  • tensor: A dense symmetric tensor.
  • basis: A matrix where the rows are the basis elements.

Returns: A sparse tensor of the form {(a,b,...,c):M_ab...c, ...}.

Example

We construct a simple tensor and then convert it to a dense array. We consider the same example as for the filter_tensor_indices function, but now we have to specify the basis in matrix form.

from cytools.utils import symmetric_dense_to_sparse
tensor = [[1,2],[2,3]]
symmetric_dense_to_sparse(tensor)
# {(0, 0): 1, (0, 1): 2, (1, 1): 3}

symmetric_sparse_to_dense

Description: Converts a symmetric sparse tensor of the form {(a,b,...,c): M_ab...c, ...} to a dense tensor.

Optionally, it then applies a basis transformation.

Arguments:

  • tensor: A sparse tensor of the form {(a,b,...,c):M_ab...c, ...}.
  • basis: A matrix where the rows are the basis elements.

Returns: A dense tensor.

Example

We construct a simple tensor and then convert it to a dense array. We consider the same example as for the filter_tensor_indices function, but now we have to specify the basis in matrix form.

from cytools.utils import symmetric_sparse_to_dense_in_basis
tensor = {(0,1):0, (1,1):1, (1,2):2, (1,3):3, (2,3):4}
basis = [[0,1,0,0],[0,0,0,1]]
symmetric_sparse_to_dense(tensor, basis)
# array([[1, 3],
# [3, 0]])

to_sparse

Description: Converts a (manually implemented) sparse matrix of the form [[a,b,M_ab], ...] or a dictionary of the form {(a,b):M_ab, ...} to a formal dok_matrix or to a csr_matrix.

Arguments:

  • arr: A list of the form [[a,b,M_ab],...] or a dictionary of the form [(a,b):M_ab,...].
  • sparse_type: The type of sparse matrix to return. The options are "dok" and "csr".

Returns: The sparse matrix in dok_matrix or csr_matrix format.

Example

We convert the 5x5 identity matrix into a sparse matrix.

from cytools.utils import to_sparse
id_array = [[0,0,1],[1,1,1],[2,2,1],[3,3,1],[4,4,1]]
to_sparse(id_array)
# <5x5 sparse matrix of type '<class 'numpy.float64'>'
# with 5 stored elements in Dictionary Of Keys format>
to_sparse(id_array, sparse_type="csr")
# <5x5 sparse matrix of type '<class 'numpy.float64'>'
# with 5 stored elements in Compressed Sparse Row format>

Functions in cytools.config

check_mosek_license

Description: Checks if the Mosek license is valid. If it is not, it prints the reason.

Arguments: None.

Returns: Nothing.

Example

The Mosek license should be automatically checked, but it can also be checked as follows.

import cytools
cytools.config.check_mosek_license()
# It will print an error if it is not working, and if nothing is printed
# then it is working correctly

enable_experimental_features

Description: Enables the experimental features of CYTools. For more information read the experimental features page.

Arguments: None.

Returns: Nothing.

Example

We enable the experimental features.

import cytools
cytools.config.enable_experimental_features()

set_mosek_path

Description: Sets a custom path to the Mosek license. This is useful if the Docker image was built without the license, but it is stored somewhere in your computer. The license will be checked after the new path is set.

Arguments:

  • path (str): The path to the Mosek license. Note that the mounted directory on the Docker container is /home/cytools/mounted_volume/.

Returns: Nothing.

Example

We set the path to the original one for illustation purposes, and show how to set the path to a directory on the host computer.

import cytools
cytools.config.set_mosek_path("/opt/cytools/external/mosek/mosek.lic") # Original path
cytools.config.set_mosek_path("/home/cytools/mounted_volume/[path-to-license]") # If not in the Docker image

Functions in cytools.__init__

check_for_updates

Description: Checks for updates of CYTools. It prints a message if a new version is avaiable, and displays a warning if the current version has a serious bug.

Arguments: None.

Returns: Nothing.

Example

We check for updates of CYTools. This is done automatically, so there is usually no need to do this.

import cytools
cytools.check_for_updates()