testlib.cc

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