conversion.cc

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