/* * Copyright 2008-2009 Katholieke Universiteit Leuven * * Use of this software is governed by the MIT license * * Written by Sven Verdoolaege, K.U.Leuven, Departement * Computerwetenschappen, Celestijnenlaan 200A, B-3001 Leuven, Belgium */ #include <isl_ctx_private.h> #include <isl_map_private.h> #include "isl_sample.h" #include <isl/vec.h> #include <isl/mat.h> #include <isl_seq.h> #include "isl_equalities.h" #include "isl_tab.h" #include "isl_basis_reduction.h" #include <isl_factorization.h> #include <isl_point_private.h> #include <isl_options_private.h> #include <isl_vec_private.h> #include <bset_from_bmap.c> #include <set_to_map.c> static __isl_give isl_vec *isl_basic_set_sample_bounded( __isl_take isl_basic_set *bset); static __isl_give isl_vec *empty_sample(__isl_take isl_basic_set *bset) { … } /* Construct a zero sample of the same dimension as bset. * As a special case, if bset is zero-dimensional, this * function creates a zero-dimensional sample point. */ static __isl_give isl_vec *zero_sample(__isl_take isl_basic_set *bset) { … } static __isl_give isl_vec *interval_sample(__isl_take isl_basic_set *bset) { … } /* Find a sample integer point, if any, in bset, which is known * to have equalities. If bset contains no integer points, then * return a zero-length vector. * We simply remove the known equalities, compute a sample * in the resulting bset, using the specified recurse function, * and then transform the sample back to the original space. */ static __isl_give isl_vec *sample_eq(__isl_take isl_basic_set *bset, __isl_give isl_vec *(*recurse)(__isl_take isl_basic_set *)) { … } /* Return a matrix containing the equalities of the tableau * in constraint form. The tableau is assumed to have * an associated bset that has been kept up-to-date. */ static struct isl_mat *tab_equalities(struct isl_tab *tab) { … } /* Compute and return an initial basis for the bounded tableau "tab". * * If the tableau is either full-dimensional or zero-dimensional, * the we simply return an identity matrix. * Otherwise, we construct a basis whose first directions correspond * to equalities. */ static struct isl_mat *initial_basis(struct isl_tab *tab) { … } /* Compute the minimum of the current ("level") basis row over "tab" * and store the result in position "level" of "min". * * This function assumes that at least one more row and at least * one more element in the constraint array are available in the tableau. */ static enum isl_lp_result compute_min(isl_ctx *ctx, struct isl_tab *tab, __isl_keep isl_vec *min, int level) { … } /* Compute the maximum of the current ("level") basis row over "tab" * and store the result in position "level" of "max". * * This function assumes that at least one more row and at least * one more element in the constraint array are available in the tableau. */ static enum isl_lp_result compute_max(isl_ctx *ctx, struct isl_tab *tab, __isl_keep isl_vec *max, int level) { … } /* Perform a greedy search for an integer point in the set represented * by "tab", given that the minimal rational value (rounded up to the * nearest integer) at "level" is smaller than the maximal rational * value (rounded down to the nearest integer). * * Return 1 if we have found an integer point (if tab->n_unbounded > 0 * then we may have only found integer values for the bounded dimensions * and it is the responsibility of the caller to extend this solution * to the unbounded dimensions). * Return 0 if greedy search did not result in a solution. * Return -1 if some error occurred. * * We assign a value half-way between the minimum and the maximum * to the current dimension and check if the minimal value of the * next dimension is still smaller than (or equal) to the maximal value. * We continue this process until either * - the minimal value (rounded up) is greater than the maximal value * (rounded down). In this case, greedy search has failed. * - we have exhausted all bounded dimensions, meaning that we have * found a solution. * - the sample value of the tableau is integral. * - some error has occurred. */ static int greedy_search(isl_ctx *ctx, struct isl_tab *tab, __isl_keep isl_vec *min, __isl_keep isl_vec *max, int level) { … } /* Given a tableau representing a set, find and return * an integer point in the set, if there is any. * * We perform a depth first search * for an integer point, by scanning all possible values in the range * attained by a basis vector, where an initial basis may have been set * by the calling function. Otherwise an initial basis that exploits * the equalities in the tableau is created. * tab->n_zero is currently ignored and is clobbered by this function. * * The tableau is allowed to have unbounded direction, but then * the calling function needs to set an initial basis, with the * unbounded directions last and with tab->n_unbounded set * to the number of unbounded directions. * Furthermore, the calling functions needs to add shifted copies * of all constraints involving unbounded directions to ensure * that any feasible rational value in these directions can be rounded * up to yield a feasible integer value. * In particular, let B define the given basis x' = B x * and let T be the inverse of B, i.e., X = T x'. * Let a x + c >= 0 be a constraint of the set represented by the tableau, * or a T x' + c >= 0 in terms of the given basis. Assume that * the bounded directions have an integer value, then we can safely * round up the values for the unbounded directions if we make sure * that x' not only satisfies the original constraint, but also * the constraint "a T x' + c + s >= 0" with s the sum of all * negative values in the last n_unbounded entries of "a T". * The calling function therefore needs to add the constraint * a x + c + s >= 0. The current function then scans the first * directions for an integer value and once those have been found, * it can compute "T ceil(B x)" to yield an integer point in the set. * Note that during the search, the first rows of B may be changed * by a basis reduction, but the last n_unbounded rows of B remain * unaltered and are also not mixed into the first rows. * * The search is implemented iteratively. "level" identifies the current * basis vector. "init" is true if we want the first value at the current * level and false if we want the next value. * * At the start of each level, we first check if we can find a solution * using greedy search. If not, we continue with the exhaustive search. * * The initial basis is the identity matrix. If the range in some direction * contains more than one integer value, we perform basis reduction based * on the value of ctx->opt->gbr * - ISL_GBR_NEVER: never perform basis reduction * - ISL_GBR_ONCE: only perform basis reduction the first * time such a range is encountered * - ISL_GBR_ALWAYS: always perform basis reduction when * such a range is encountered * * When ctx->opt->gbr is set to ISL_GBR_ALWAYS, then we allow the basis * reduction computation to return early. That is, as soon as it * finds a reasonable first direction. */ __isl_give isl_vec *isl_tab_sample(struct isl_tab *tab) { … } static __isl_give isl_vec *sample_bounded(__isl_take isl_basic_set *bset); /* Internal data for factored_sample. * "sample" collects the sample and may get reset to a zero-length vector * signaling the absence of a sample vector. * "pos" is the position of the contribution of the next factor. */ struct isl_factored_sample_data { … }; /* isl_factorizer_every_factor_basic_set callback that extends * the sample in data->sample with the contribution * of the factor "bset". * If "bset" turns out to be empty, then the product is empty too and * no further factors need to be considered. */ static isl_bool factor_sample(__isl_keep isl_basic_set *bset, void *user) { … } /* Compute a sample point of the given basic set, based on the given, * non-trivial factorization. */ static __isl_give isl_vec *factored_sample(__isl_take isl_basic_set *bset, __isl_take isl_factorizer *f) { … } /* Given a basic set that is known to be bounded, find and return * an integer point in the basic set, if there is any. * * After handling some trivial cases, we construct a tableau * and then use isl_tab_sample to find a sample, passing it * the identity matrix as initial basis. */ static __isl_give isl_vec *sample_bounded(__isl_take isl_basic_set *bset) { … } /* Given a basic set "bset" and a value "sample" for the first coordinates * of bset, plug in these values and drop the corresponding coordinates. * * We do this by computing the preimage of the transformation * * [ 1 0 ] * x = [ s 0 ] x' * [ 0 I ] * * where [1 s] is the sample value and I is the identity matrix of the * appropriate dimension. */ static __isl_give isl_basic_set *plug_in(__isl_take isl_basic_set *bset, __isl_take isl_vec *sample) { … } /* Given a basic set "bset", return any (possibly non-integer) point * in the basic set. */ static __isl_give isl_vec *rational_sample(__isl_take isl_basic_set *bset) { … } /* Given a linear cone "cone" and a rational point "vec", * construct a polyhedron with shifted copies of the constraints in "cone", * i.e., a polyhedron with "cone" as its recession cone, such that each * point x in this polyhedron is such that the unit box positioned at x * lies entirely inside the affine cone 'vec + cone'. * Any rational point in this polyhedron may therefore be rounded up * to yield an integer point that lies inside said affine cone. * * Denote the constraints of cone by "<a_i, x> >= 0" and the rational * point "vec" by v/d. * Let b_i = <a_i, v>. Then the affine cone 'vec + cone' is given * by <a_i, x> - b/d >= 0. * The polyhedron <a_i, x> - ceil{b/d} >= 0 is a subset of this affine cone. * We prefer this polyhedron over the actual affine cone because it doesn't * require a scaling of the constraints. * If each of the vertices of the unit cube positioned at x lies inside * this polyhedron, then the whole unit cube at x lies inside the affine cone. * We therefore impose that x' = x + \sum e_i, for any selection of unit * vectors lies inside the polyhedron, i.e., * * <a_i, x'> - ceil{b/d} = <a_i, x> + sum a_i - ceil{b/d} >= 0 * * The most stringent of these constraints is the one that selects * all negative a_i, so the polyhedron we are looking for has constraints * * <a_i, x> + sum_{a_i < 0} a_i - ceil{b/d} >= 0 * * Note that if cone were known to have only non-negative rays * (which can be accomplished by a unimodular transformation), * then we would only have to check the points x' = x + e_i * and we only have to add the smallest negative a_i (if any) * instead of the sum of all negative a_i. */ static __isl_give isl_basic_set *shift_cone(__isl_take isl_basic_set *cone, __isl_take isl_vec *vec) { … } /* Given a rational point vec in a (transformed) basic set, * such that cone is the recession cone of the original basic set, * "round up" the rational point to an integer point. * * We first check if the rational point just happens to be integer. * If not, we transform the cone in the same way as the basic set, * pick a point x in this cone shifted to the rational point such that * the whole unit cube at x is also inside this affine cone. * Then we simply round up the coordinates of x and return the * resulting integer point. */ static __isl_give isl_vec *round_up_in_cone(__isl_take isl_vec *vec, __isl_take isl_basic_set *cone, __isl_take isl_mat *U) { … } /* Concatenate two integer vectors, i.e., two vectors with denominator * (stored in element 0) equal to 1. */ static __isl_give isl_vec *vec_concat(__isl_take isl_vec *vec1, __isl_take isl_vec *vec2) { … } /* Give a basic set "bset" with recession cone "cone", compute and * return an integer point in bset, if any. * * If the recession cone is full-dimensional, then we know that * bset contains an infinite number of integer points and it is * fairly easy to pick one of them. * If the recession cone is not full-dimensional, then we first * transform bset such that the bounded directions appear as * the first dimensions of the transformed basic set. * We do this by using a unimodular transformation that transforms * the equalities in the recession cone to equalities on the first * dimensions. * * The transformed set is then projected onto its bounded dimensions. * Note that to compute this projection, we can simply drop all constraints * involving any of the unbounded dimensions since these constraints * cannot be combined to produce a constraint on the bounded dimensions. * To see this, assume that there is such a combination of constraints * that produces a constraint on the bounded dimensions. This means * that some combination of the unbounded dimensions has both an upper * bound and a lower bound in terms of the bounded dimensions, but then * this combination would be a bounded direction too and would have been * transformed into a bounded dimensions. * * We then compute a sample value in the bounded dimensions. * If no such value can be found, then the original set did not contain * any integer points and we are done. * Otherwise, we plug in the value we found in the bounded dimensions, * project out these bounded dimensions and end up with a set with * a full-dimensional recession cone. * A sample point in this set is computed by "rounding up" any * rational point in the set. * * The sample points in the bounded and unbounded dimensions are * then combined into a single sample point and transformed back * to the original space. */ __isl_give isl_vec *isl_basic_set_sample_with_cone( __isl_take isl_basic_set *bset, __isl_take isl_basic_set *cone) { … } static void vec_sum_of_neg(__isl_keep isl_vec *v, isl_int *s) { … } /* Given a tableau "tab", a tableau "tab_cone" that corresponds * to the recession cone and the inverse of a new basis U = inv(B), * with the unbounded directions in B last, * add constraints to "tab" that ensure any rational value * in the unbounded directions can be rounded up to an integer value. * * The new basis is given by x' = B x, i.e., x = U x'. * For any rational value of the last tab->n_unbounded coordinates * in the update tableau, the value that is obtained by rounding * up this value should be contained in the original tableau. * For any constraint "a x + c >= 0", we therefore need to add * a constraint "a x + c + s >= 0", with s the sum of all negative * entries in the last elements of "a U". * * Since we are not interested in the first entries of any of the "a U", * we first drop the columns of U that correpond to bounded directions. */ static int tab_shift_cone(struct isl_tab *tab, struct isl_tab *tab_cone, struct isl_mat *U) { … } /* Compute and return an initial basis for the possibly * unbounded tableau "tab". "tab_cone" is a tableau * for the corresponding recession cone. * Additionally, add constraints to "tab" that ensure * that any rational value for the unbounded directions * can be rounded up to an integer value. * * If the tableau is bounded, i.e., if the recession cone * is zero-dimensional, then we just use inital_basis. * Otherwise, we construct a basis whose first directions * correspond to equalities, followed by bounded directions, * i.e., equalities in the recession cone. * The remaining directions are then unbounded. */ int isl_tab_set_initial_basis_with_cone(struct isl_tab *tab, struct isl_tab *tab_cone) { … } /* Compute and return a sample point in bset using generalized basis * reduction. We first check if the input set has a non-trivial * recession cone. If so, we perform some extra preprocessing in * sample_with_cone. Otherwise, we directly perform generalized basis * reduction. */ static __isl_give isl_vec *gbr_sample(__isl_take isl_basic_set *bset) { … } static __isl_give isl_vec *basic_set_sample(__isl_take isl_basic_set *bset, int bounded) { … } __isl_give isl_vec *isl_basic_set_sample_vec(__isl_take isl_basic_set *bset) { … } /* Compute an integer sample in "bset", where the caller guarantees * that "bset" is bounded. */ __isl_give isl_vec *isl_basic_set_sample_bounded(__isl_take isl_basic_set *bset) { … } __isl_give isl_basic_set *isl_basic_set_from_vec(__isl_take isl_vec *vec) { … } __isl_give isl_basic_map *isl_basic_map_sample(__isl_take isl_basic_map *bmap) { … } __isl_give isl_basic_set *isl_basic_set_sample(__isl_take isl_basic_set *bset) { … } __isl_give isl_basic_map *isl_map_sample(__isl_take isl_map *map) { … } __isl_give isl_basic_set *isl_set_sample(__isl_take isl_set *set) { … } __isl_give isl_point *isl_basic_set_sample_point(__isl_take isl_basic_set *bset) { … } __isl_give isl_point *isl_set_sample_point(__isl_take isl_set *set) { … }
Codebase Browser
llvm