param_util.c

   1 #include <barvinok/options.h>
   2 #include "param_util.h"
   3
   4 #define ALLOC(type) (type*)malloc(sizeof(type))
   5 #define ALLOCN(type,n) (type*)malloc((n) * sizeof(type))
   6
   7 void Param_Vertex_Common_Denominator(Param_Vertices *V)
   8 {
   9     unsigned dim;
  10     Value lcm;
  11     int i;
  12
  13     assert(V->Vertex->NbRows > 0);
  14     dim = V->Vertex->NbColumns-2;
  15
  16     value_init(lcm);
  17
  18     value_assign(lcm, V->Vertex->p[0][dim+1]);
  19     for (i = 1; i < V->Vertex->NbRows; ++i)
  20         value_lcm(V->Vertex->p[i][dim+1], lcm, &lcm);
  21
  22     for (i = 0; i < V->Vertex->NbRows; ++i) {
  23         if (value_eq(V->Vertex->p[i][dim+1], lcm))
  24             continue;
  25         value_division(V->Vertex->p[i][dim+1], lcm, V->Vertex->p[i][dim+1]);
  26         Vector_Scale(V->Vertex->p[i], V->Vertex->p[i],
  27                      V->Vertex->p[i][dim+1], dim+1);
  28         value_assign(V->Vertex->p[i][dim+1], lcm);
  29     }
  30
  31     value_clear(lcm);
  32 }
  33
  34 /* Plug in the parametric vertex Vertex (nvar x (nparam + 2))
  35  * in the constraint constraint (1 + nvar + nparam + 1).
  36  * The result is stored in row (1 + nparam + 1),
  37  * with the denominator in position 0.
  38  */
  39 void Param_Inner_Product(Value *constraint, Matrix *Vertex, Value *row)
  40 {
  41     unsigned nparam = Vertex->NbColumns - 2;
  42     unsigned nvar = Vertex->NbRows;
  43     int j;
  44     Value tmp, tmp2;
  45
  46     value_set_si(row[0], 1);
  47     Vector_Set(row+1, 0, nparam+1);
  48
  49     value_init(tmp);
  50     value_init(tmp2);
  51
  52     for (j = 0 ; j < nvar; ++j) {
  53         value_set_si(tmp, 1);
  54         value_assign(tmp2,  constraint[1+j]);
  55         if (value_ne(row[0], Vertex->p[j][nparam+1])) {
  56             value_assign(tmp, row[0]);
  57             value_lcm(row[0], Vertex->p[j][nparam+1], &row[0]);
  58             value_division(tmp, row[0], tmp);
  59             value_multiply(tmp2, tmp2, row[0]);
  60             value_division(tmp2, tmp2, Vertex->p[j][nparam+1]);
  61         }
  62         Vector_Combine(row+1, Vertex->p[j], row+1, tmp, tmp2, nparam+1);
  63     }
  64     value_set_si(tmp, 1);
  65     Vector_Combine(row+1, constraint+1+nvar, row+1, tmp, row[0], nparam+1);
  66
  67     value_clear(tmp);
  68     value_clear(tmp2);
  69 }
  70
  71 Param_Polyhedron *Polyhedron2Param_Polyhedron(Polyhedron *Din, Polyhedron *Cin,
  72                                               struct barvinok_options *options)
  73 {
  74     unsigned MaxRays = options->MaxRays;
  75     if (MaxRays & POL_NO_DUAL)
  76         MaxRays = 0;
  77     return Polyhedron2Param_Domain(Din, Cin, MaxRays);
  78 }
  79
  80 /* Compute a dummy Param_Domain that contains all vertices of Param_Domain D
  81  * (which contains the vertices of P) that lie on the facet obtained by
  82  * saturating constraint c of P
  83  */
  84 Param_Domain *Param_Polyhedron_Facet(Param_Polyhedron *PP, Param_Domain *D,
  85                                      Polyhedron *P, int c)
  86 {
  87     int nv;
  88     Param_Vertices *V;
  89     unsigned nparam = PP->V->Vertex->NbColumns-2;
  90     Vector *row = Vector_Alloc(1+nparam+1);
  91     Param_Domain *FD = ALLOC(Param_Domain);
  92     FD->Domain = 0;
  93     FD->next = 0;
  94
  95     nv = (PP->nbV - 1)/(8*sizeof(int)) + 1;
  96     FD->F = ALLOCN(unsigned, nv);
  97     memset(FD->F, 0, nv * sizeof(unsigned));
  98
  99     FORALL_PVertex_in_ParamPolyhedron(V, D, PP) /* _i, _ix, _bx internal counters */
 100         int n;
 101         Param_Inner_Product(P->Constraint[c], V->Vertex, row->p);
 102         if (First_Non_Zero(row->p+1, nparam+1) == -1)
 103             FD->F[_ix] |= _bx;
 104     END_FORALL_PVertex_in_ParamPolyhedron;
 105
 106     Vector_Free(row);
 107
 108     return FD;
 109 }