polysign_polylib.c

   1 #include <assert.h>
   2 #include <barvinok/util.h>
   3 #include <barvinok/options.h>
   4 #include "polysign.h"
   5
   6 static enum order_sign interval_minmax(Polyhedron *I)
   7 {
   8     int i;
   9     int min = 1;
  10     int max = -1;
  11     assert(I->Dimension == 1);
  12     POL_ENSURE_VERTICES(I);
  13     for (i = 0; i < I->NbRays; ++i) {
  14         if (value_pos_p(I->Ray[i][1]))
  15             max = 1;
  16         else if (value_neg_p(I->Ray[i][1]))
  17             min = -1;
  18         else {
  19             if (max < 0)
  20                 max = 0;
  21             if (min > 0)
  22                 min = 0;
  23         }
  24     }
  25     if (min > 0)
  26         return order_gt;
  27     if (max < 0)
  28         return order_lt;
  29     if (min == max)
  30         return order_eq;
  31     if (max == 0)
  32         return order_le;
  33     if (min == 0)
  34         return order_ge;
  35     return order_unknown;
  36 }
  37
  38 /* Returns the sign of the affine function specified by T on the polyhedron D */
  39 enum order_sign PL_polyhedron_affine_sign(Polyhedron *D, Matrix *T,
  40                                             struct barvinok_options *options)
  41 {
  42     enum order_sign sign;
  43     Polyhedron *I = Polyhedron_Image(D, T, options->MaxRays);
  44     if (POL_ISSET(options->MaxRays, POL_INTEGER))
  45         I = DomainConstraintSimplify(I, options->MaxRays);
  46     if (emptyQ2(I)) {
  47         Polyhedron_Free(I);
  48         I = Polyhedron_Image(D, T, options->MaxRays);
  49     }
  50     sign = interval_minmax(I);
  51     Polyhedron_Free(I);
  52     return sign;
  53 }
  54
  55 enum lp_result PL_constraints_opt(Matrix *C, Value *obj, Value denom,
  56                                     enum lp_dir dir, Value *opt,
  57                                     unsigned MaxRays)
  58 {
  59     enum lp_result res;
  60     Polyhedron *P = Constraints2Polyhedron(C, MaxRays);
  61     res = PL_polyhedron_opt(P, obj, denom, dir, opt);
  62     Polyhedron_Free(P);
  63     return res;
  64 }
  65
  66 enum lp_result PL_polyhedron_opt(Polyhedron *P, Value *obj, Value denom,
  67                                 enum lp_dir dir, Value *opt)
  68 {
  69     int i;
  70     int first = 1;
  71     Value val, d;
  72     enum lp_result res = lp_empty;
  73
  74     POL_ENSURE_VERTICES(P);
  75     if (emptyQ(P))
  76         return res;
  77
  78     value_init(val);
  79     value_init(d);
  80     for (i = 0; i < P->NbRays; ++ i) {
  81         Inner_Product(P->Ray[i]+1, obj, P->Dimension+1, &val);
  82         if (value_zero_p(P->Ray[i][0]) && value_notzero_p(val)) {
  83             res = lp_unbounded;
  84             break;
  85         }
  86         if (value_zero_p(P->Ray[i][1+P->Dimension])) {
  87             if ((dir == lp_min && value_neg_p(val)) ||
  88                 (dir == lp_max && value_pos_p(val))) {
  89                 res = lp_unbounded;
  90                 break;
  91             }
  92         } else {
  93             res = lp_ok;
  94             value_multiply(d, denom, P->Ray[i][1+P->Dimension]);
  95             if (dir == lp_min)
  96                 mpz_cdiv_q(val, val, d);
  97             else
  98                 mpz_fdiv_q(val, val, d);
  99             if (first || (dir == lp_min ? value_lt(val, *opt) :
 100                                           value_gt(val, *opt)))
 101                 value_assign(*opt, val);
 102             first = 0;
 103         }
 104     }
 105     value_clear(d);
 106     value_clear(val);
 107
 108     return res;
 109 }
 110
 111 enum lp_result PL_polyhedron_range(Polyhedron *D, Value *obj, Value denom,
 112                                 Value *min, Value *max,
 113                                 struct barvinok_options *options)
 114 {
 115     Polyhedron *I;
 116     int bounded;
 117
 118     Matrix *T = Matrix_Alloc(2, D->Dimension+1);
 119     Vector_Copy(obj, T->p[0], D->Dimension+1);
 120     value_assign(T->p[1][D->Dimension], denom);
 121
 122     I = Polyhedron_Image(D, T, options->MaxRays);
 123     Matrix_Free(T);
 124     POL_ENSURE_VERTICES(I);
 125
 126     if (emptyQ(I)) {
 127         Polyhedron_Free(I);
 128         return lp_empty;
 129     }
 130     bounded = line_minmax(I, min, max);
 131
 132     if (!bounded)
 133         return lp_unbounded;
 134     if (value_gt(*min, *max))
 135         return lp_empty;
 136     return lp_ok;
 137 }