conversion.cc

   1 #include <gmp.h>
   2 #include <NTL/mat_ZZ.h>
   3 extern "C" {
   4 #include <polylib/polylibgmp.h>
   5 }
   6 #include <barvinok/util.h>
   7 #include "conversion.h"
   8
   9 #define SIZE(p) (((long *) (p))[1])
  10 #define DATA(p) ((mp_limb_t *) (((long *) (p)) + 2))
  11
  12 void value2zz(Value v, ZZ& z)
  13 {
  14     int sa = v[0]._mp_size;
  15     int abs_sa = sa < 0 ? -sa : sa;
  16
  17     _ntl_gsetlength(&z.rep, abs_sa);
  18     mp_limb_t * adata = DATA(z.rep);
  19     for (int i = 0; i < abs_sa; ++i)
  20         adata[i] = v[0]._mp_d[i];
  21     SIZE(z.rep) = sa;
  22 }
  23
  24 void zz2value(const ZZ& z, Value& v)
  25 {
  26     if (!z.rep) {
  27         value_set_si(v, 0);
  28         return;
  29     }
  30
  31     int sa = SIZE(z.rep);
  32     int abs_sa = sa < 0 ? -sa : sa;
  33
  34     mp_limb_t * adata = DATA(z.rep);
  35     _mpz_realloc(v, abs_sa);
  36     for (int i = 0; i < abs_sa; ++i)
  37         v[0]._mp_d[i] = adata[i];
  38     v[0]._mp_size = sa;
  39 }
  40
  41 void values2zz(Value *p, vec_ZZ& v, int len)
  42 {
  43     v.SetLength(len);
  44
  45     for (int i = 0; i < len; ++i) {
  46         value2zz(p[i], v[i]);
  47     }
  48 }
  49
  50 /*
  51  */
  52 void zz2values(const vec_ZZ& v, Value *p)
  53 {
  54     for (int i = 0; i < v.length(); ++i)
  55         zz2value(v[i], p[i]);
  56 }
  57
  58 /*
  59  * We just ignore the last column and row
  60  * If the final element is not equal to one
  61  * then the result will actually be a multiple of the input
  62  */
  63 void matrix2zz(Matrix *M, mat_ZZ& m, unsigned nr, unsigned nc)
  64 {
  65     m.SetDims(nr, nc);
  66
  67     for (int i = 0; i < nr; ++i) {
  68 //      assert(value_one_p(M->p[i][M->NbColumns - 1]));
  69         for (int j = 0; j < nc; ++j) {
  70             value2zz(M->p[i][j], m[i][j]);
  71         }
  72     }
  73 }
  74
  75 Matrix *rays(Polyhedron *C)
  76 {
  77     unsigned dim = C->NbRays - 1; /* don't count zero vertex */
  78     assert(C->NbRays - 1 == C->Dimension);
  79
  80     Matrix *M = Matrix_Alloc(dim+1, dim+1);
  81     assert(M);
  82
  83     int i, c;
  84     for (i = 0, c = 0; i <= dim && c < dim; ++i)
  85         if (value_zero_p(C->Ray[i][dim+1])) {
  86             Vector_Copy(C->Ray[i] + 1, M->p[c], dim);
  87             value_set_si(M->p[c++][dim], 0);
  88         }
  89     assert(c == dim);
  90     value_set_si(M->p[dim][dim], 1);
  91
  92     return M;
  93 }
  94
  95 Matrix *rays2(Polyhedron *C)
  96 {
  97     unsigned dim = C->NbRays - 1; /* don't count zero vertex */
  98     assert(C->NbRays - 1 == C->Dimension);
  99
 100     Matrix *M = Matrix_Alloc(dim, dim);
 101     assert(M);
 102
 103     int i, c;
 104     for (i = 0, c = 0; i <= dim && c < dim; ++i)
 105         if (value_zero_p(C->Ray[i][dim+1]))
 106             Vector_Copy(C->Ray[i] + 1, M->p[c++], dim);
 107     assert(c == dim);
 108
 109     return M;
 110 }
 111
 112 void randomvector(Polyhedron *P, vec_ZZ& lambda, int nvar)
 113 {
 114     Value tmp;
 115     int max = 10 * 16;
 116     unsigned int dim = P->Dimension;
 117     value_init(tmp);
 118
 119     for (int i = 0; i < P->NbRays; ++i) {
 120         for (int j = 1; j <= dim; ++j) {
 121             value_absolute(tmp, P->Ray[i][j]);
 122             int t = VALUE_TO_LONG(tmp) * 16;
 123             if (t > max)
 124                 max = t;
 125         }
 126     }
 127     for (int i = 0; i < P->NbConstraints; ++i) {
 128         for (int j = 1; j <= dim; ++j) {
 129             value_absolute(tmp, P->Constraint[i][j]);
 130             int t = VALUE_TO_LONG(tmp) * 16;
 131             if (t > max)
 132                 max = t;
 133         }
 134     }
 135     value_clear(tmp);
 136
 137     lambda.SetLength(nvar);
 138     for (int k = 0; k < nvar; ++k) {
 139         int r = random_int(max*dim)+2;
 140         int v = (2*(r%2)-1) * (max/2*dim + (r >> 1));
 141         lambda[k] = v;
 142     }
 143 }