add isl_aff_mod_val
[platform/upstream/isl.git] / isl_sample.c
index ec63745..c6a0f31 100644 (file)
@@ -1,7 +1,7 @@
 /*
  * Copyright 2008-2009 Katholieke Universiteit Leuven
  *
- * Use of this software is governed by the GNU LGPLv2.1 license
+ * 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
@@ -347,6 +347,105 @@ static struct isl_mat *initial_basis(struct isl_tab *tab)
        return Q;
 }
 
+/* Compute the minimum of the current ("level") basis row over "tab"
+ * and store the result in position "level" of "min".
+ */
+static enum isl_lp_result compute_min(isl_ctx *ctx, struct isl_tab *tab,
+       __isl_keep isl_vec *min, int level)
+{
+       return isl_tab_min(tab, tab->basis->row[1 + level],
+                           ctx->one, &min->el[level], NULL, 0);
+}
+
+/* Compute the maximum of the current ("level") basis row over "tab"
+ * and store the result in position "level" of "max".
+ */
+static enum isl_lp_result compute_max(isl_ctx *ctx, struct isl_tab *tab,
+       __isl_keep isl_vec *max, int level)
+{
+       enum isl_lp_result res;
+       unsigned dim = tab->n_var;
+
+       isl_seq_neg(tab->basis->row[1 + level] + 1,
+                   tab->basis->row[1 + level] + 1, dim);
+       res = isl_tab_min(tab, tab->basis->row[1 + level],
+                   ctx->one, &max->el[level], NULL, 0);
+       isl_seq_neg(tab->basis->row[1 + level] + 1,
+                   tab->basis->row[1 + level] + 1, dim);
+       isl_int_neg(max->el[level], max->el[level]);
+
+       return res;
+}
+
+/* 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)
+{
+       struct isl_tab_undo *snap;
+       enum isl_lp_result res;
+
+       snap = isl_tab_snap(tab);
+
+       do {
+               isl_int_add(tab->basis->row[1 + level][0],
+                           min->el[level], max->el[level]);
+               isl_int_fdiv_q_ui(tab->basis->row[1 + level][0],
+                           tab->basis->row[1 + level][0], 2);
+               isl_int_neg(tab->basis->row[1 + level][0],
+                           tab->basis->row[1 + level][0]);
+               if (isl_tab_add_valid_eq(tab, tab->basis->row[1 + level]) < 0)
+                       return -1;
+               isl_int_set_si(tab->basis->row[1 + level][0], 0);
+
+               if (++level >= tab->n_var - tab->n_unbounded)
+                       return 1;
+               if (isl_tab_sample_is_integer(tab))
+                       return 1;
+
+               res = compute_min(ctx, tab, min, level);
+               if (res == isl_lp_error)
+                       return -1;
+               if (res != isl_lp_ok)
+                       isl_die(ctx, isl_error_internal,
+                               "expecting bounded rational solution",
+                               return -1);
+               res = compute_max(ctx, tab, max, level);
+               if (res == isl_lp_error)
+                       return -1;
+               if (res != isl_lp_ok)
+                       isl_die(ctx, isl_error_internal,
+                               "expecting bounded rational solution",
+                               return -1);
+       } while (isl_int_le(min->el[level], max->el[level]));
+
+       if (isl_tab_rollback(tab, snap) < 0)
+               return -1;
+
+       return 0;
+}
+
 /* Given a tableau representing a set, find and return
  * an integer point in the set, if there is any.
  *
@@ -386,6 +485,9 @@ static struct isl_mat *initial_basis(struct isl_tab *tab)
  * 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
@@ -455,34 +557,38 @@ struct isl_vec *isl_tab_sample(struct isl_tab *tab)
        reduced = 0;
 
        while (level >= 0) {
-               int empty = 0;
                if (init) {
-                       res = isl_tab_min(tab, tab->basis->row[1 + level],
-                                   ctx->one, &min->el[level], NULL, 0);
-                       if (res == isl_lp_empty)
-                               empty = 1;
-                       isl_assert(ctx, res != isl_lp_unbounded, goto error);
+                       int choice;
+
+                       res = compute_min(ctx, tab, min, level);
                        if (res == isl_lp_error)
                                goto error;
-                       if (!empty && isl_tab_sample_is_integer(tab))
+                       if (res != isl_lp_ok)
+                               isl_die(ctx, isl_error_internal,
+                                       "expecting bounded rational solution",
+                                       goto error);
+                       if (isl_tab_sample_is_integer(tab))
                                break;
-                       isl_seq_neg(tab->basis->row[1 + level] + 1,
-                                   tab->basis->row[1 + level] + 1, dim);
-                       res = isl_tab_min(tab, tab->basis->row[1 + level],
-                                   ctx->one, &max->el[level], NULL, 0);
-                       isl_seq_neg(tab->basis->row[1 + level] + 1,
-                                   tab->basis->row[1 + level] + 1, dim);
-                       isl_int_neg(max->el[level], max->el[level]);
-                       if (res == isl_lp_empty)
-                               empty = 1;
-                       isl_assert(ctx, res != isl_lp_unbounded, goto error);
+                       res = compute_max(ctx, tab, max, level);
                        if (res == isl_lp_error)
                                goto error;
-                       if (!empty && isl_tab_sample_is_integer(tab))
+                       if (res != isl_lp_ok)
+                               isl_die(ctx, isl_error_internal,
+                                       "expecting bounded rational solution",
+                                       goto error);
+                       if (isl_tab_sample_is_integer(tab))
                                break;
-                       if (!empty && !reduced &&
-                           ctx->opt->gbr != ISL_GBR_NEVER &&
-                           isl_int_lt(min->el[level], max->el[level])) {
+                       choice = isl_int_lt(min->el[level], max->el[level]);
+                       if (choice) {
+                               int g;
+                               g = greedy_search(ctx, tab, min, max, level);
+                               if (g < 0)
+                                       goto error;
+                               if (g)
+                                       break;
+                       }
+                       if (!reduced && choice &&
+                           ctx->opt->gbr != ISL_GBR_NEVER) {
                                unsigned gbr_only_first;
                                if (ctx->opt->gbr == ISL_GBR_ONCE)
                                        ctx->opt->gbr = ISL_GBR_NEVER;
@@ -502,7 +608,7 @@ struct isl_vec *isl_tab_sample(struct isl_tab *tab)
                } else
                        isl_int_add_ui(min->el[level], min->el[level], 1);
 
-               if (empty || isl_int_gt(min->el[level], max->el[level])) {
+               if (isl_int_gt(min->el[level], max->el[level])) {
                        level--;
                        init = 0;
                        if (level >= 0)
@@ -657,7 +763,7 @@ static struct isl_vec *sample_bounded(struct isl_basic_set *bset)
                
        ctx = bset->ctx;
 
-       tab = isl_tab_from_basic_set(bset);
+       tab = isl_tab_from_basic_set(bset, 1);
        if (tab && tab->empty) {
                isl_tab_free(tab);
                ISL_F_SET(bset, ISL_BASIC_SET_EMPTY);
@@ -666,8 +772,6 @@ static struct isl_vec *sample_bounded(struct isl_basic_set *bset)
                return sample;
        }
 
-       if (isl_tab_track_bset(tab, isl_basic_set_copy(bset)) < 0)
-               goto error;
        if (!ISL_F_ISSET(bset, ISL_BASIC_SET_NO_IMPLICIT))
                if (isl_tab_detect_implicit_equalities(tab) < 0)
                        goto error;
@@ -746,7 +850,7 @@ static struct isl_vec *rational_sample(struct isl_basic_set *bset)
        if (!bset)
                return NULL;
 
-       tab = isl_tab_from_basic_set(bset);
+       tab = isl_tab_from_basic_set(bset, 0);
        sample = isl_tab_get_sample_value(tab);
        isl_tab_free(tab);
 
@@ -1289,6 +1393,11 @@ error:
        return NULL;
 }
 
+__isl_give isl_basic_set *isl_basic_set_sample(__isl_take isl_basic_set *bset)
+{
+       return isl_basic_map_sample(bset);
+}
+
 __isl_give isl_basic_map *isl_map_sample(__isl_take isl_map *map)
 {
        int i;