testlib.cc

   1 #include <sstream>
   2 #include <barvinok/set.h>
   3 #include <barvinok/options.h>
   4 #include <barvinok/evalue.h>
   5 #include <barvinok/util.h>
   6 #include "conversion.h"
   7 #include "evalue_read.h"
   8 #include "dpoly.h"
   9 #include "lattice_point.h"
  10 #include "tcounter.h"
  11 #include "bernoulli.h"
  12
  13 using std::cout;
  14 using std::cerr;
  15 using std::endl;
  16
  17 template <typename T>
  18 void set_from_string(T& v, const char *s)
  19 {
  20     std::istringstream str(s);
  21     str >> v;
  22 }
  23
  24 int test_evalue(struct barvinok_options *options)
  25 {
  26     unsigned nvar, nparam;
  27     char **all_vars;
  28     evalue *poly1, poly2;
  29
  30     poly1 = evalue_read_from_str("(1/4 * n^4 + 1/2 * n^3 + 1/4 * n^2)",
  31                                  NULL, &all_vars, &nvar, &nparam,
  32                                  options->MaxRays);
  33     Free_ParamNames(all_vars, nvar+nparam);
  34
  35     value_init(poly2.d);
  36     evalue_copy(&poly2, poly1);
  37     evalue_negate(poly1);
  38     eadd(&poly2, poly1);
  39     reduce_evalue(poly1);
  40     assert(EVALUE_IS_ZERO(*poly1));
  41     free_evalue_refs(poly1);
  42     free_evalue_refs(&poly2);
  43     free(poly1);
  44 }
  45
  46 int test_split_periods(struct barvinok_options *options)
  47 {
  48     unsigned nvar, nparam;
  49     char **all_vars;
  50     evalue *e;
  51
  52     e = evalue_read_from_str("U + 2V + 3 >= 0\n- U -2V  >= 0\n- U  10 >= 0\n"
  53                              "U  >= 0\n\n({( 1/3 * U + ( 2/3 * V + 0 ))})",
  54                              NULL, &all_vars, &nvar, &nparam,
  55                              options->MaxRays);
  56     Free_ParamNames(all_vars, nvar+nparam);
  57
  58     evalue_split_periods(e, 2, options->MaxRays);
  59     assert(value_zero_p(e->d));
  60     assert(e->x.p->type == partition);
  61     assert(e->x.p->size == 4);
  62     assert(value_zero_p(e->x.p->arr[1].d));
  63     assert(e->x.p->arr[1].x.p->type == polynomial);
  64     assert(value_zero_p(e->x.p->arr[3].d));
  65     assert(e->x.p->arr[3].x.p->type == polynomial);
  66     free_evalue_refs(e);
  67     free(e);
  68 }
  69
  70 int test_specialization(struct barvinok_options *options)
  71 {
  72     Value v;
  73     value_init(v);
  74     mpq_t count;
  75     mpq_init(count);
  76     ZZ sign;
  77
  78     value_set_si(v, 5);
  79     dpoly n(2, v);
  80     assert(value_cmp_si(n.coeff->p[0], 1) == 0);
  81     assert(value_cmp_si(n.coeff->p[1], 5) == 0);
  82     assert(value_cmp_si(n.coeff->p[2], 10) == 0);
  83
  84     value_set_si(v, 1);
  85     dpoly d(2, v, 1);
  86     value_set_si(v, 2);
  87     dpoly d2(2, v, 1);
  88     d *= d2;
  89     assert(value_cmp_si(d.coeff->p[0], 2) == 0);
  90     assert(value_cmp_si(d.coeff->p[1], 1) == 0);
  91     assert(value_cmp_si(d.coeff->p[2], 0) == 0);
  92
  93     sign = 1;
  94     n.div(d, count, sign);
  95     mpq_canonicalize(count);
  96     assert(value_cmp_si(mpq_numref(count), 31) == 0);
  97     assert(value_cmp_si(mpq_denref(count), 8) == 0);
  98
  99     value_set_si(v, -2);
 100     dpoly n2(2, v);
 101     assert(value_cmp_si(n2.coeff->p[0], 1) == 0);
 102     assert(value_cmp_si(n2.coeff->p[1], -2) == 0);
 103     assert(value_cmp_si(n2.coeff->p[2], 3) == 0);
 104
 105     n2.div(d, count, sign);
 106     mpq_canonicalize(count);
 107     assert(value_cmp_si(mpq_numref(count), 6) == 0);
 108     assert(value_cmp_si(mpq_denref(count), 1) == 0);
 109
 110     value_clear(v);
 111     mpq_clear(count);
 112 }
 113
 114 int test_lattice_points(struct barvinok_options *options)
 115 {
 116     Param_Vertices V;
 117     mat_ZZ tmp;
 118     set_from_string(tmp, "[[0 0 0 0 4][0 0 0 0 4][-1 0 1 0 4]]");
 119     V.Vertex = zz2matrix(tmp);
 120     vec_ZZ lambda;
 121     set_from_string(lambda, "[3 5 7]");
 122     mat_ZZ rays;
 123     set_from_string(rays, "[[0 1 0][4 0 1][0 0 -1]]");
 124     term_info num;
 125     evalue *point[4];
 126
 127     unsigned nvar, nparam;
 128     char **all_vars;
 129     point[0] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
 130                     "( 7 * {( 1/4 * a + ( 3/4 * c + 3/4 ) ) } + -21/4 ) ) )",
 131                     "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
 132     Free_ParamNames(all_vars, nvar+nparam);
 133     point[1] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
 134                     "( 7 * {( 1/4 * a + ( 3/4 * c + 1/2 ) ) } + -1/2 ) ) )",
 135                     "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
 136     Free_ParamNames(all_vars, nvar+nparam);
 137     point[2] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
 138                     "( 7 * {( 1/4 * a + ( 3/4 * c + 1/4 ) ) } + 17/4 ) ) )",
 139                     "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
 140     Free_ParamNames(all_vars, nvar+nparam);
 141     point[3] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
 142                     "( 7 * {( 1/4 * a + ( 3/4 * c + 0 ) ) } + 9 ) ) )",
 143                     "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
 144     Free_ParamNames(all_vars, nvar+nparam);
 145
 146     lattice_point(&V, rays, lambda, &num, 4, NULL, options);
 147     Matrix_Free(V.Vertex);
 148
 149     for (int i = 0; i < 4; ++i) {
 150         assert(eequal(num.E[i], point[i]));
 151         free_evalue_refs(point[i]);
 152         free(point[i]);
 153         free_evalue_refs(num.E[i]);
 154         delete num.E[i];
 155     }
 156     delete [] num.E;
 157 }
 158
 159 int test_todd(struct barvinok_options *options)
 160 {
 161     tcounter t(2, options->max_index);
 162     assert(value_cmp_si(t.todd.coeff->p[0], 1) == 0);
 163     assert(value_cmp_si(t.todd.coeff->p[1], -3) == 0);
 164     assert(value_cmp_si(t.todd.coeff->p[2], 3) == 0);
 165     assert(value_cmp_si(t.todd_denom->p[0], 1) == 0);
 166     assert(value_cmp_si(t.todd_denom->p[1], 6) == 0);
 167     assert(value_cmp_si(t.todd_denom->p[2], 36) == 0);
 168
 169     vec_ZZ lambda;
 170     set_from_string(lambda, "[1 -1]");
 171     zz2values(lambda, t.lambda->p);
 172
 173     mat_ZZ rays;
 174     set_from_string(rays, "[[-1 0][-1 1]]");
 175
 176     QQ one(1, 1);
 177
 178     vec_ZZ v;
 179     set_from_string(v, "[2 0 1]");
 180     Vector *vertex = Vector_Alloc(3);
 181     zz2values(v, vertex->p);
 182
 183     t.handle(rays, vertex->p, one, 1, NULL, options);
 184     assert(value_cmp_si(mpq_numref(t.count), 71) == 0);
 185     assert(value_cmp_si(mpq_denref(t.count), 24) == 0);
 186
 187     set_from_string(rays, "[[0 -1][1 -1]]");
 188     set_from_string(v, "[0 2 1]");
 189     zz2values(v, vertex->p);
 190
 191     t.handle(rays, vertex->p, one, 1, NULL, options);
 192     assert(value_cmp_si(mpq_numref(t.count), 71) == 0);
 193     assert(value_cmp_si(mpq_denref(t.count), 12) == 0);
 194
 195     set_from_string(rays, "[[1 0][0 1]]");
 196     set_from_string(v, "[0 0 1]");
 197     zz2values(v, vertex->p);
 198
 199     t.handle(rays, vertex->p, one, 1, NULL, options);
 200     assert(value_cmp_si(mpq_numref(t.count), 6) == 0);
 201     assert(value_cmp_si(mpq_denref(t.count), 1) == 0);
 202
 203     Vector_Free(vertex);
 204 }
 205
 206 int test_bernoulli(struct barvinok_options *options)
 207 {
 208     struct bernoulli *bernoulli;
 209     bernoulli = bernoulli_compute(2);
 210     bernoulli = bernoulli_compute(4);
 211     bernoulli = bernoulli_compute(8);
 212     assert(value_cmp_si(bernoulli->b_num->p[6], 1) == 0);
 213     assert(value_cmp_si(bernoulli->b_den->p[6], 42) == 0);
 214     assert(value_cmp_si(bernoulli->sum[3]->p[0], 0) == 0);
 215     assert(value_cmp_si(bernoulli->sum[3]->p[1], 0) == 0);
 216     assert(value_cmp_si(bernoulli->sum[3]->p[2], 1) == 0);
 217     assert(value_cmp_si(bernoulli->sum[3]->p[3], -2) == 0);
 218     assert(value_cmp_si(bernoulli->sum[3]->p[4], 1) == 0);
 219     assert(value_cmp_si(bernoulli->sum[3]->p[5], 4) == 0);
 220
 221     unsigned nvar, nparam;
 222     char **all_vars;
 223     evalue *base, *sum1, *sum2;
 224     base = evalue_read_from_str("(1 * n + 1)", NULL, &all_vars, &nvar, &nparam,
 225                                 options->MaxRays);
 226
 227     sum1 = evalue_polynomial(bernoulli->sum[3], base);
 228     Free_ParamNames(all_vars, nvar+nparam);
 229
 230     sum2 = evalue_read_from_str("(1/4 * n^4 + 1/2 * n^3 + 1/4 * n^2)",
 231                                 NULL, &all_vars, &nvar, &nparam,
 232                                 options->MaxRays);
 233     Free_ParamNames(all_vars, nvar+nparam);
 234     assert(eequal(sum1, sum2));
 235     free_evalue_refs(base);
 236     free_evalue_refs(sum1);
 237     free_evalue_refs(sum2);
 238     free(base);
 239     free(sum1);
 240     free(sum2);
 241 }
 242
 243 int test_bernoulli_sum(struct barvinok_options *options)
 244 {
 245     unsigned nvar, nparam;
 246     char **all_vars;
 247     evalue *e, *sum1, *sum2;
 248     e = evalue_read_from_str("i + -1 >= 0\n -i + n >= 0\n\n 1 + (-1 *i) + i^2",
 249                              "i", &all_vars, &nvar, &nparam,
 250                                 options->MaxRays);
 251     Free_ParamNames(all_vars, nvar+nparam);
 252
 253     sum1 = Bernoulli_sum_evalue(e, 1, options);
 254     sum2 = evalue_read_from_str("n -1 >= 0\n\n (1/3 * n^3 + 2/3 * n)",
 255                                 NULL, &all_vars, &nvar, &nparam,
 256                                 options->MaxRays);
 257     Free_ParamNames(all_vars, nvar+nparam);
 258     evalue_negate(sum1);
 259     eadd(sum2, sum1);
 260     reduce_evalue(sum1);
 261     assert(EVALUE_IS_ZERO(*sum1));
 262     free_evalue_refs(e);
 263     free_evalue_refs(sum1);
 264     free(e);
 265     free(sum1);
 266
 267     e = evalue_read_from_str("-i + -1 >= 0\n i + n >= 0\n\n 1 + i + i^2",
 268                              "i", &all_vars, &nvar, &nparam,
 269                                 options->MaxRays);
 270     Free_ParamNames(all_vars, nvar+nparam);
 271     sum1 = Bernoulli_sum_evalue(e, 1, options);
 272     evalue_negate(sum1);
 273     eadd(sum2, sum1);
 274     reduce_evalue(sum1);
 275     assert(EVALUE_IS_ZERO(*sum1));
 276     free_evalue_refs(e);
 277     free_evalue_refs(sum1);
 278     free(e);
 279     free(sum1);
 280
 281     free_evalue_refs(sum2);
 282     free(sum2);
 283
 284     e = evalue_read_from_str("i + 4 >= 0\n -i + n >= 0\n\n i",
 285                              "i", &all_vars, &nvar, &nparam,
 286                                 options->MaxRays);
 287     Free_ParamNames(all_vars, nvar+nparam);
 288     sum1 = Bernoulli_sum_evalue(e, 1, options);
 289     sum2 = evalue_read_from_str("n + 0 >= 0\n\n (1/2 * n^2 + 1/2 * n + (-10))\n"
 290                     "n + 4 >= 0\n -n -1 >= 0\n\n (1/2 * n^2 + 1/2 * n + (-10))",
 291                                 NULL, &all_vars, &nvar, &nparam,
 292                                 options->MaxRays);
 293     Free_ParamNames(all_vars, nvar+nparam);
 294     evalue_negate(sum1);
 295     eadd(sum2, sum1);
 296     reduce_evalue(sum1);
 297     assert(EVALUE_IS_ZERO(*sum1));
 298     free_evalue_refs(e);
 299     free_evalue_refs(sum1);
 300     free(e);
 301     free(sum1);
 302     free_evalue_refs(sum2);
 303     free(sum2);
 304
 305     e = evalue_read_from_str("i -5 >= 0\n -i + n >= 0\n j -1 >= 0\n i -j >= 0\n"
 306                              "k -1 >= 0\n j -k >= 0\n\n1",
 307                              "i,j,k", &all_vars, &nvar, &nparam,
 308                                 options->MaxRays);
 309     Free_ParamNames(all_vars, nvar+nparam);
 310     sum1 = Bernoulli_sum_evalue(e, 3, options);
 311     sum2 = evalue_read_from_str("n -5 >= 0\n\n"
 312                                 "1/6 * n^3 + 1/2 * n^2 + 1/3 * n + -20",
 313                                 NULL, &all_vars, &nvar, &nparam,
 314                                 options->MaxRays);
 315     Free_ParamNames(all_vars, nvar+nparam);
 316     evalue_negate(sum1);
 317     eadd(sum2, sum1);
 318     reduce_evalue(sum1);
 319     assert(EVALUE_IS_ZERO(*sum1));
 320     free_evalue_refs(e);
 321     free_evalue_refs(sum1);
 322     free(e);
 323     free(sum1);
 324     free_evalue_refs(sum2);
 325     free(sum2);
 326 }
 327
 328 int main(int argc, char **argv)
 329 {
 330     struct barvinok_options *options = barvinok_options_new_with_defaults();
 331     test_evalue(options);
 332     test_split_periods(options);
 333     test_specialization(options);
 334     test_lattice_points(options);
 335     test_todd(options);
 336     test_bernoulli(options);
 337     test_bernoulli_sum(options);
 338     barvinok_options_free(options);
 339 }