isl_sample.c

   1 #include "isl_sample.h"
   2 #include "isl_sample_piplib.h"
   3 #include "isl_vec.h"
   4 #include "isl_mat.h"
   5 #include "isl_seq.h"
   6 #include "isl_map_private.h"
   7 #include "isl_equalities.h"
   8
   9 static struct isl_vec *empty_sample(struct isl_basic_set *bset)
  10 {
  11         struct isl_vec *vec;
  12
  13         vec = isl_vec_alloc(bset->ctx, 0);
  14         isl_basic_set_free(bset);
  15         return vec;
  16 }
  17
  18 /* Construct a zero sample of the same dimension as bset.
  19  * As a special case, if bset is zero-dimensional, this
  20  * function creates a zero-dimensional sample point.
  21  */
  22 static struct isl_vec *zero_sample(struct isl_basic_set *bset)
  23 {
  24         unsigned dim;
  25         struct isl_vec *sample;
  26
  27         dim = isl_basic_set_total_dim(bset);
  28         sample = isl_vec_alloc(bset->ctx, 1 + dim);
  29         if (sample) {
  30                 isl_int_set_si(sample->el[0], 1);
  31                 isl_seq_clr(sample->el + 1, dim);
  32         }
  33         isl_basic_set_free(bset);
  34         return sample;
  35 }
  36
  37 static struct isl_vec *interval_sample(struct isl_ctx *ctx,
  38         struct isl_basic_set *bset)
  39 {
  40         int i;
  41         isl_int t;
  42         struct isl_vec *sample;
  43
  44         bset = isl_basic_set_simplify(bset);
  45         if (!bset)
  46                 return NULL;
  47         if (isl_basic_set_fast_is_empty(bset))
  48                 return empty_sample(bset);
  49         if (bset->n_eq == 0 && bset->n_ineq == 0)
  50                 return zero_sample(bset);
  51
  52         sample = isl_vec_alloc(ctx, 2);
  53         isl_int_set_si(sample->block.data[0], 1);
  54
  55         if (bset->n_eq > 0) {
  56                 isl_assert(bset->ctx, bset->n_eq == 1, goto error);
  57                 isl_assert(bset->ctx, bset->n_ineq == 0, goto error);
  58                 if (isl_int_is_one(bset->eq[0][1]))
  59                         isl_int_neg(sample->el[1], bset->eq[0][0]);
  60                 else {
  61                         isl_assert(bset->ctx, isl_int_is_negone(bset->eq[0][1]),
  62                                    goto error);
  63                         isl_int_set(sample->el[1], bset->eq[0][0]);
  64                 }
  65                 isl_basic_set_free(bset);
  66                 return sample;
  67         }
  68
  69         isl_int_init(t);
  70         if (isl_int_is_one(bset->ineq[0][1]))
  71                 isl_int_neg(sample->block.data[1], bset->ineq[0][0]);
  72         else
  73                 isl_int_set(sample->block.data[1], bset->ineq[0][0]);
  74         for (i = 1; i < bset->n_ineq; ++i) {
  75                 isl_seq_inner_product(sample->block.data,
  76                                         bset->ineq[i], 2, &t);
  77                 if (isl_int_is_neg(t))
  78                         break;
  79         }
  80         isl_int_clear(t);
  81         if (i < bset->n_ineq) {
  82                 isl_vec_free(sample);
  83                 return empty_sample(bset);
  84         }
  85
  86         isl_basic_set_free(bset);
  87         return sample;
  88 error:
  89         isl_basic_set_free(bset);
  90         isl_vec_free(sample);
  91         return NULL;
  92 }
  93
  94 static struct isl_mat *independent_bounds(struct isl_ctx *ctx,
  95         struct isl_basic_set *bset)
  96 {
  97         int i, j, n;
  98         struct isl_mat *dirs = NULL;
  99         unsigned dim;
 100
 101         if (!bset)
 102                 return NULL;
 103
 104         dim = isl_basic_set_n_dim(bset);
 105         if (bset->n_ineq == 0)
 106                 return isl_mat_alloc(ctx, 0, dim);
 107
 108         dirs = isl_mat_alloc(ctx, dim, dim);
 109         if (!dirs)
 110                 return NULL;
 111         isl_seq_cpy(dirs->row[0], bset->ineq[0]+1, dirs->n_col);
 112         for (j = 1, n = 1; n < dim && j < bset->n_ineq; ++j) {
 113                 int pos;
 114
 115                 isl_seq_cpy(dirs->row[n], bset->ineq[j]+1, dirs->n_col);
 116
 117                 pos = isl_seq_first_non_zero(dirs->row[n], dirs->n_col);
 118                 if (pos < 0)
 119                         continue;
 120                 for (i = 0; i < n; ++i) {
 121                         int pos_i;
 122                         pos_i = isl_seq_first_non_zero(dirs->row[i], dirs->n_col);
 123                         if (pos_i < pos)
 124                                 continue;
 125                         if (pos_i > pos)
 126                                 break;
 127                         isl_seq_elim(dirs->row[n], dirs->row[i], pos,
 128                                         dirs->n_col, NULL);
 129                         pos = isl_seq_first_non_zero(dirs->row[n], dirs->n_col);
 130                         if (pos < 0)
 131                                 break;
 132                 }
 133                 if (pos < 0)
 134                         continue;
 135                 if (i < n) {
 136                         int k;
 137                         isl_int *t = dirs->row[n];
 138                         for (k = n; k > i; --k)
 139                                 dirs->row[k] = dirs->row[k-1];
 140                         dirs->row[i] = t;
 141                 }
 142                 ++n;
 143         }
 144         dirs->n_row = n;
 145         return dirs;
 146 }
 147
 148 static struct isl_basic_set *remove_lineality(struct isl_ctx *ctx,
 149         struct isl_basic_set *bset, struct isl_mat *bounds, struct isl_mat **T)
 150 {
 151         struct isl_mat *U = NULL;
 152         unsigned old_dim, new_dim;
 153
 154         old_dim = isl_basic_set_n_dim(bset);
 155         new_dim = bounds->n_row;
 156         *T = NULL;
 157         bounds = isl_mat_left_hermite(ctx, bounds, 0, &U, NULL);
 158         if (!bounds)
 159                 goto error;
 160         U = isl_mat_lin_to_aff(ctx, U);
 161         U = isl_mat_drop_cols(ctx, U, 1 + new_dim, old_dim - new_dim);
 162         bset = isl_basic_set_preimage(bset, isl_mat_copy(ctx, U));
 163         if (!bset)
 164                 goto error;
 165         *T = U;
 166         isl_mat_free(ctx, bounds);
 167         return bset;
 168 error:
 169         isl_mat_free(ctx, bounds);
 170         isl_mat_free(ctx, U);
 171         isl_basic_set_free(bset);
 172         return NULL;
 173 }
 174
 175 struct isl_vec *isl_basic_set_sample(struct isl_basic_set *bset)
 176 {
 177         struct isl_ctx *ctx;
 178         struct isl_mat *bounds;
 179         unsigned dim;
 180         if (!bset)
 181                 return NULL;
 182
 183         ctx = bset->ctx;
 184         if (isl_basic_set_fast_is_empty(bset))
 185                 return empty_sample(bset);
 186
 187         dim = isl_basic_set_n_dim(bset);
 188         isl_assert(ctx, isl_basic_set_n_param(bset) == 0, goto error);
 189         isl_assert(ctx, bset->n_div == 0, goto error);
 190
 191         if (bset->n_eq > 0) {
 192                 struct isl_mat *T;
 193                 struct isl_vec *sample;
 194
 195                 bset = isl_basic_set_remove_equalities(bset, &T, NULL);
 196                 sample = isl_basic_set_sample(bset);
 197                 if (sample && sample->size != 0)
 198                         sample = isl_mat_vec_product(ctx, T, sample);
 199                 else
 200                         isl_mat_free(ctx, T);
 201                 return sample;
 202         }
 203         if (dim == 0)
 204                 return zero_sample(bset);
 205         if (dim == 1)
 206                 return interval_sample(ctx, bset);
 207         bounds = independent_bounds(ctx, bset);
 208         if (!bounds)
 209                 goto error;
 210         if (bounds->n_row == dim)
 211                 isl_mat_free(ctx, bounds);
 212         else {
 213                 struct isl_mat *T;
 214                 struct isl_vec *sample;
 215
 216                 bset = remove_lineality(ctx, bset, bounds, &T);
 217                 sample = isl_basic_set_sample(bset);
 218                 if (sample && sample->size != 0)
 219                         sample = isl_mat_vec_product(ctx, T, sample);
 220                 else
 221                         isl_mat_free(ctx, T);
 222                 return sample;
 223         }
 224         return isl_pip_basic_set_sample(bset);
 225 error:
 226         isl_basic_set_free(bset);
 227         return NULL;
 228 }