polysign_isl.c

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