polysign_isl.c

   1 #include <assert.h>
   2 #include <isl/val_gmp.h>
   3 #include <isl/lp.h>
   4 #include <isl_set_polylib.h>
   5 #include "polysign.h"
   6
   7 static __isl_give isl_mat *extract_equalities(isl_ctx *ctx, Matrix *M)
   8 {
   9         int i, j;
  10         int n;
  11         isl_val *v;
  12         isl_mat *eq;
  13
  14         n = 0;
  15         for (i = 0; i < M->NbRows; ++i)
  16                 if (value_zero_p(M->p[i][0]))
  17                         ++n;
  18
  19         eq = isl_mat_alloc(ctx, n, M->NbColumns - 1);
  20         for (i = 0; i < M->NbRows; ++i) {
  21                 if (!value_zero_p(M->p[i][0]))
  22                         continue;
  23                 for (j = 0; j < M->NbColumns - 1; ++j) {
  24                         v = isl_val_int_from_gmp(ctx, M->p[i][1 + j]);
  25                         eq = isl_mat_set_element_val(eq, i, j, v);
  26                 }
  27         }
  28
  29         return eq;
  30 }
  31
  32 static __isl_give isl_mat *extract_inequalities(isl_ctx *ctx, Matrix *M)
  33 {
  34         int i, j;
  35         int n;
  36         isl_val *v;
  37         isl_mat *ineq;
  38
  39         n = 0;
  40         for (i = 0; i < M->NbRows; ++i)
  41                 if (!value_zero_p(M->p[i][0]))
  42                         ++n;
  43
  44         ineq = isl_mat_alloc(ctx, n, M->NbColumns - 1);
  45         for (i = 0; i < M->NbRows; ++i) {
  46                 if (value_zero_p(M->p[i][0]))
  47                         continue;
  48                 for (j = 0; j < M->NbColumns - 1; ++j) {
  49                         v = isl_val_int_from_gmp(ctx, M->p[i][1 + j]);
  50                         ineq = isl_mat_set_element_val(ineq, i, j, v);
  51                 }
  52         }
  53
  54         return ineq;
  55 }
  56
  57 enum order_sign isl_polyhedron_affine_sign(Polyhedron *D, Matrix *T,
  58                                             struct barvinok_options *options)
  59 {
  60         int i;
  61         isl_ctx *ctx = isl_ctx_alloc();
  62         isl_space *dim;
  63         isl_local_space *ls;
  64         isl_aff *aff;
  65         isl_basic_set *bset;
  66         isl_val *min, *max = NULL;
  67         isl_val *v;
  68         enum order_sign sign = order_undefined;
  69
  70         assert(D->Dimension == T->NbColumns - 1);
  71
  72         dim = isl_space_set_alloc(ctx, 0, D->Dimension);
  73         ls = isl_local_space_from_space(isl_space_copy(dim));
  74         bset = isl_basic_set_new_from_polylib(D, dim);
  75         aff = isl_aff_zero_on_domain(ls);
  76         for (i = 0; i < D->Dimension; ++i) {
  77                 v = isl_val_int_from_gmp(ctx, T->p[0][i]);
  78                 aff = isl_aff_set_coefficient_val(aff, isl_dim_in, i, v);
  79         }
  80         v = isl_val_int_from_gmp(ctx, T->p[0][D->Dimension]);
  81         aff = isl_aff_set_constant_val(aff, v);
  82         v = isl_val_int_from_gmp(ctx, T->p[1][D->Dimension]);
  83         aff = isl_aff_scale_down_val(aff, v);
  84
  85         min = isl_basic_set_min_lp_val(bset, aff);
  86         min = isl_val_ceil(min);
  87         assert(min);
  88
  89         if (isl_val_is_nan(min))
  90                 sign = order_undefined;
  91         else if (isl_val_is_pos(min))
  92                 sign = order_gt;
  93         else {
  94                 max = isl_basic_set_max_lp_val(bset, aff);
  95                 max = isl_val_floor(max);
  96                 assert(max);
  97
  98                 if (isl_val_is_neg(max))
  99                         sign = order_lt;
 100                 else if (isl_val_is_zero(min) && isl_val_is_zero(max))
 101                         sign = order_eq;
 102                 else if (isl_val_is_zero(min))
 103                         sign = order_ge;
 104                 else if (isl_val_is_zero(max))
 105                         sign = order_le;
 106                 else
 107                         sign = order_unknown;
 108         }
 109
 110         isl_basic_set_free(bset);
 111         isl_aff_free(aff);
 112         isl_val_free(min);
 113         isl_val_free(max);
 114         isl_ctx_free(ctx);
 115
 116         return sign;
 117 }
 118
 119 static enum lp_result isl_lp_result2lp_result(enum isl_lp_result res)
 120 {
 121         switch (res) {
 122         case isl_lp_ok:
 123                 return lp_ok;
 124         case isl_lp_unbounded:
 125                 return lp_unbounded;
 126         case isl_lp_empty:
 127                 return lp_empty;
 128         default:
 129                 assert(0);
 130         }
 131 }
 132
 133 enum lp_result isl_constraints_opt(Matrix *C, Value *obj, Value denom,
 134                                     enum lp_dir dir, Value *opt)
 135 {
 136         int i;
 137         isl_ctx *ctx = isl_ctx_alloc();
 138         isl_space *dim;
 139         isl_local_space *ls;
 140         isl_mat *eq, *ineq;
 141         isl_basic_set *bset;
 142         isl_aff *aff;
 143         isl_val *v;
 144         enum isl_lp_result res;
 145         int max = dir == lp_max;
 146
 147         eq = extract_equalities(ctx, C);
 148         ineq = extract_inequalities(ctx, C);
 149         dim = isl_space_set_alloc(ctx, 0, C->NbColumns - 2);
 150         ls = isl_local_space_from_space(isl_space_copy(dim));
 151         bset = isl_basic_set_from_constraint_matrices(dim, eq, ineq,
 152                         isl_dim_set, isl_dim_div, isl_dim_param, isl_dim_cst);
 153         aff = isl_aff_zero_on_domain(ls);
 154         for (i = 0; i < C->NbColumns - 2; ++i) {
 155                 v = isl_val_int_from_gmp(ctx, obj[i]);
 156                 aff = isl_aff_set_coefficient_val(aff, isl_dim_in, i, v);
 157         }
 158         v = isl_val_int_from_gmp(ctx, obj[C->NbColumns - 2]);
 159         aff = isl_aff_set_constant_val(aff, v);
 160         v = isl_val_int_from_gmp(ctx, denom);
 161         aff = isl_aff_scale_down_val(aff, v);
 162
 163         if (max)
 164                 v = isl_val_floor(isl_basic_set_max_lp_val(bset, aff));
 165         else
 166                 v = isl_val_ceil(isl_basic_set_min_lp_val(bset, aff));
 167         if (!v)
 168                 res = isl_lp_error;
 169         else if (isl_val_is_nan(v))
 170                 res = isl_lp_empty;
 171         else if (!isl_val_is_rat(v))
 172                 res = isl_lp_unbounded;
 173         else {
 174                 res = isl_lp_ok;
 175                 isl_val_get_num_gmp(v, *opt);
 176         }
 177
 178         isl_val_free(v);
 179         isl_aff_free(aff);
 180         isl_basic_set_free(bset);
 181         isl_ctx_free(ctx);
 182
 183         return isl_lp_result2lp_result(res);
 184 }