?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
laurent.cc: laurent_summator: member initializer list: put base classes first
[barvinok.git] / testlib.cc
blob77d907f20a8b56ce0fe24e55e04820315ac66392
1 #include <assert.h>
2 #include <stdlib.h>
3 #include <sstream>
4 #include <barvinok/barvinok.h>
5 #include <barvinok/set.h>
6 #include <barvinok/options.h>
7 #include <barvinok/basis_reduction.h>
8 #include <barvinok/evalue.h>
9 #include <barvinok/util.h>
10 #include "conversion.h"
11 #include "evalue_read.h"
12 #include "dpoly.h"
13 #include "lattice_point.h"
14 #include "counter.h"
15 #include "bernoulli.h"
16 #include "hilbert.h"
17 #include "hull.h"
18 #include "ilp.h"
19 #include "laurent.h"
20 #include "matrix_read.h"
21 #include "remove_equalities.h"
22 #include "config.h"
23
24 using std::cout;
25 using std::cerr;
26 using std::endl;
27
28 template <typename T>
29 void set_from_string(T& v, const char *s)
30 {
31 std::istringstream str(s);
32 str >> v;
33 }
34
35 static Matrix *matrix_read_from_str(const char *s)
36 {
37 std::istringstream str(s);
38 return Matrix_Read(str);
39 }
40
41 static int test_equalities(struct barvinok_options *options)
42 {
43 Matrix *M = matrix_read_from_str(
44 "11 11\n"
45 " 0 23 0 0 -10 0 0 0 7 -44 -8 \n"
46 " 0 0 23 0 5 0 0 0 -15 114 27 \n"
47 " 0 0 0 1 0 0 0 0 0 -1 2 \n"
48 " 0 0 0 0 0 1 0 0 -1 8 0 \n"
49 " 0 0 0 0 0 0 1 0 0 -1 2 \n"
50 " 0 0 0 0 0 0 0 1 0 -1 -1 \n"
51 " 1 0 0 0 10 0 0 0 -7 44 8 \n"
52 " 1 0 0 0 -5 0 0 0 15 -114 -27 \n"
53 " 1 0 0 0 1 0 0 0 0 0 0 \n"
54 " 1 0 0 0 0 0 0 0 1 -8 0 \n"
55 " 1 0 0 0 0 0 0 0 0 1 -2 \n");
56 Polyhedron *P = Constraints2Polyhedron(M, options->MaxRays);
57 Matrix_Free(M);
58 Polyhedron *Q = P;
59 remove_all_equalities(&P, NULL, NULL, NULL, 2, options->MaxRays);
60 assert(P->NbEq == 0);
61 if (P != Q)
62 Polyhedron_Free(Q);
63 Polyhedron_Free(P);
64 return 0;
65 }
66
67 int test_evalue_read(struct barvinok_options *options)
68 {
69 unsigned nvar, nparam;
70 const char **all_vars;
71 evalue *e1, *e2;
72
73 e1 = evalue_read_from_str("(1 * aa + 2 * a)",
74 NULL, &all_vars, &nvar, &nparam, options->MaxRays);
75 Free_ParamNames(all_vars, nvar+nparam);
76 e2 = evalue_read_from_str("(3 * aa)",
77 NULL, &all_vars, &nvar, &nparam, options->MaxRays);
78 Free_ParamNames(all_vars, nvar+nparam);
79 assert(!eequal(e1, e2));
80 evalue_free(e1);
81 evalue_free(e2);
82 return 0;
83 }
84
85 static void evalue_check_disjoint(evalue *e)
86 {
87 int i, j;
88
89 if (!e)
90 return;
91 if (EVALUE_IS_ZERO(*e))
92 return;
93 for (i = 0; i < e->x.p->size/2; ++i) {
94 Polyhedron *A = EVALUE_DOMAIN(e->x.p->arr[2*i]);
95 for (j = i+1; j < e->x.p->size/2; ++j) {
96 Polyhedron *B = EVALUE_DOMAIN(e->x.p->arr[2*j]);
97 Polyhedron *I = DomainIntersection(A, B, 0);
98 assert(emptyQ(I));
99 Polyhedron_Free(I);
100 }
101 }
102 }
103
104 static int test_eadd(struct barvinok_options *options)
105 {
106 unsigned nvar, nparam;
107 const char **all_vars;
108 evalue *e1, *e2;
109
110 e1 = evalue_read_from_str(" d -1 = 0\n"
111 " h -3 >= 0\n"
112 " - h + 100 >= 0\n"
113 "\n"
114 "(1)\n",
115 "d,h", &all_vars, &nvar, &nparam,
116 options->MaxRays);
117 Free_ParamNames(all_vars, nvar+nparam);
118 e2 = evalue_read_from_str(
119 " h -3 = 0\n"
120 " d -1 >= 0\n"
121 " - d + 3 >= 0\n"
122 "\n"
123 "(0)\n"
124 " d -1 >= 0\n"
125 " - d + 3 >= 0\n"
126 " h -4 >= 0\n"
127 " - h + 100 >= 0\n"
128 "\n"
129 "(1)\n",
130 "d,h", &all_vars, &nvar, &nparam,
131 options->MaxRays);
132 Free_ParamNames(all_vars, nvar+nparam);
133 eadd(e2, e1);
134 evalue_check_disjoint(e1);
135 evalue_free(e1);
136 evalue_free(e2);
137 return 0;
138 }
139
140 int test_evalue(struct barvinok_options *options)
141 {
142 unsigned nvar, nparam;
143 const char **all_vars;
144 evalue *poly1, poly2;
145
146 poly1 = evalue_read_from_str("(1/4 * n^4 + 1/2 * n^3 + 1/4 * n^2)",
147 NULL, &all_vars, &nvar, &nparam,
148 options->MaxRays);
149 Free_ParamNames(all_vars, nvar+nparam);
150
151 value_init(poly2.d);
152 evalue_copy(&poly2, poly1);
153 evalue_negate(poly1);
154 eadd(&poly2, poly1);
155 reduce_evalue(poly1);
156 assert(EVALUE_IS_ZERO(*poly1));
157 evalue_free(poly1);
158 free_evalue_refs(&poly2);
159 return 0;
160 }
161
162 int test_substitute(struct barvinok_options *options)
163 {
164 unsigned nvar, nparam;
165 const char **all_vars;
166 const char *vars = "a,b";
167 evalue *e1, *e2;
168 evalue *subs[2];
169
170 e1 = evalue_read_from_str("[ { 1/3 * a } = 0 ] * \n"
171 " ([ { 1/5 * b + 2/5 } = 0 ] * 5) + \n"
172 "[ { 1/3 * a } != 0 ] * 42",
173 vars, &all_vars, &nvar, &nparam,
174 options->MaxRays);
175 Free_ParamNames(all_vars, nvar+nparam);
176
177 subs[0] = evalue_read_from_str("(2 * b + 5)",
178 vars, &all_vars, &nvar, &nparam,
179 options->MaxRays);
180 Free_ParamNames(all_vars, nvar+nparam);
181 subs[1] = evalue_read_from_str("(a + 1)",
182 vars, &all_vars, &nvar, &nparam,
183 options->MaxRays);
184 Free_ParamNames(all_vars, nvar+nparam);
185
186 evalue_substitute(e1, subs);
187 evalue_free(subs[0]);
188 evalue_free(subs[1]);
189 reduce_evalue(e1);
190
191 e2 = evalue_read_from_str("[ { 2/3 * b + 2/3 } = 0 ] * \n"
192 " ([ { 1/5 * a + 3/5 } = 0 ] * 5) + \n"
193 "[ { 2/3 * b + 2/3 } != 0 ] * 42",
194 vars, &all_vars, &nvar, &nparam,
195 options->MaxRays);
196 Free_ParamNames(all_vars, nvar+nparam);
197 reduce_evalue(e2);
198
199 assert(eequal(e1, e2));
200
201 evalue_free(e1);
202 evalue_free(e2);
203 return 0;
204 }
205
206 int test_specialization(struct barvinok_options *options)
207 {
208 Value v;
209 value_init(v);
210 mpq_t count;
211 mpq_init(count);
212
213 value_set_si(v, 5);
214 dpoly n(2, v);
215 assert(value_cmp_si(n.coeff->p[0], 1) == 0);
216 assert(value_cmp_si(n.coeff->p[1], 5) == 0);
217 assert(value_cmp_si(n.coeff->p[2], 10) == 0);
218
219 value_set_si(v, 1);
220 dpoly d(2, v, 1);
221 value_set_si(v, 2);
222 dpoly d2(2, v, 1);
223 d *= d2;
224 assert(value_cmp_si(d.coeff->p[0], 2) == 0);
225 assert(value_cmp_si(d.coeff->p[1], 1) == 0);
226 assert(value_cmp_si(d.coeff->p[2], 0) == 0);
227
228 n.div(d, count, 1);
229 mpq_canonicalize(count);
230 assert(value_cmp_si(mpq_numref(count), 31) == 0);
231 assert(value_cmp_si(mpq_denref(count), 8) == 0);
232
233 value_set_si(v, -2);
234 dpoly n2(2, v);
235 assert(value_cmp_si(n2.coeff->p[0], 1) == 0);
236 assert(value_cmp_si(n2.coeff->p[1], -2) == 0);
237 assert(value_cmp_si(n2.coeff->p[2], 3) == 0);
238
239 n2.div(d, count, 1);
240 mpq_canonicalize(count);
241 assert(value_cmp_si(mpq_numref(count), 6) == 0);
242 assert(value_cmp_si(mpq_denref(count), 1) == 0);
243
244 value_clear(v);
245 mpq_clear(count);
246 return 0;
247 }
248
249 int test_lattice_points(struct barvinok_options *options)
250 {
251 Param_Vertices V;
252 mat_ZZ tmp;
253 set_from_string(tmp, "[[0 0 0 0 4][0 0 0 0 4][-1 0 1 0 4]]");
254 V.Vertex = zz2matrix(tmp);
255 vec_ZZ lambda;
256 set_from_string(lambda, "[3 5 7]");
257 mat_ZZ rays;
258 set_from_string(rays, "[[0 1 0][4 0 1][0 0 -1]]");
259 term_info num;
260 evalue *point[4];
261
262 unsigned nvar, nparam;
263 const char **all_vars;
264 point[0] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
265 "( 7 * {( 1/4 * a + ( 3/4 * c + 3/4 ) ) } + -21/4 ) ) )",
266 "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
267 Free_ParamNames(all_vars, nvar+nparam);
268 point[1] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
269 "( 7 * {( 1/4 * a + ( 3/4 * c + 1/2 ) ) } + -1/2 ) ) )",
270 "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
271 Free_ParamNames(all_vars, nvar+nparam);
272 point[2] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
273 "( 7 * {( 1/4 * a + ( 3/4 * c + 1/4 ) ) } + 17/4 ) ) )",
274 "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
275 Free_ParamNames(all_vars, nvar+nparam);
276 point[3] = evalue_read_from_str("( -7/4 * a + ( 7/4 * c + "
277 "( 7 * {( 1/4 * a + ( 3/4 * c + 0 ) ) } + 9 ) ) )",
278 "a,b,c", &all_vars, &nvar, &nparam, options->MaxRays);
279 Free_ParamNames(all_vars, nvar+nparam);
280
281 lattice_point(&V, rays, lambda, &num, 4, options);
282 Matrix_Free(V.Vertex);
283
284 for (int i = 0; i < 4; ++i) {
285 assert(eequal(num.E[i], point[i]));
286 evalue_free(point[i]);
287 evalue_free(num.E[i]);
288 }
289 delete [] num.E;
290 return 0;
291 }
292
293 static int test_icounter(struct barvinok_options *options)
294 {
295 icounter cnt(2);
296 vec_QQ n_coeff;
297 mat_ZZ n_power;
298 mat_ZZ d_power;
299 set_from_string(n_coeff, "[-2/1 1/1]");
300 set_from_string(n_power, "[[2 6][3 6]]");
301 d_power.SetDims(0, 2);
302 cnt.reduce(n_coeff, n_power, d_power);
303 assert(value_cmp_si(mpq_numref(cnt.count), -1) == 0);
304 assert(value_cmp_si(mpq_denref(cnt.count), 1) == 0);
305 return 0;
306 }
307
308 static int test_infinite_counter(struct barvinok_options *options)
309 {
310 Matrix *M = matrix_read_from_str("1 3\n 1 1 0\n");
311 Polyhedron *ctx = Constraints2Polyhedron(M, options->MaxRays);
312 Matrix_Free(M);
313
314 /* (1 -1/2 x^5 - 1/2 x^7)/(1-x) */
315 infinite_counter *cnt = new infinite_counter(1, 1);
316 cnt->init(ctx, 0);
317 vec_QQ n_coeff;
318 mat_ZZ n_power;
319 mat_ZZ d_power;
320 set_from_string(n_coeff, "[1/1 -1/2 -1/2]");
321 set_from_string(n_power, "[[0][5][7]]");
322 set_from_string(d_power, "[[1]]");
323 cnt->reduce(n_coeff, n_power, d_power);
324 assert(value_cmp_si(mpq_numref(cnt->count[0]), 6) == 0);
325 assert(value_cmp_si(mpq_denref(cnt->count[0]), 1) == 0);
326 assert(value_cmp_si(mpq_numref(cnt->count[1]), 0) == 0);
327 assert(value_cmp_si(mpq_denref(cnt->count[1]), 1) == 0);
328 delete cnt;
329 Polyhedron_Free(ctx);
330
331 M = matrix_read_from_str("2 4\n 1 1 0 0\n 1 0 1 0\n");
332 ctx = Constraints2Polyhedron(M, options->MaxRays);
333 Matrix_Free(M);
334
335 /* (1 - xy)/((1-x)(1-xy)) */
336 cnt = new infinite_counter(2, 3);
337 cnt->init(ctx, 0);
338 set_from_string(n_coeff, "[1/1 -1/1]");
339 set_from_string(n_power, "[[0 0][1 1]]");
340 set_from_string(d_power, "[[1 0][1 1]]");
341 cnt->reduce(n_coeff, n_power, d_power);
342 assert(value_cmp_si(mpq_numref(cnt->count[1]), 0) != 0);
343 assert(value_cmp_si(mpq_numref(cnt->count[2]), 0) == 0);
344 assert(value_cmp_si(mpq_denref(cnt->count[2]), 1) == 0);
345 assert(value_cmp_si(mpq_numref(cnt->count[3]), 0) == 0);
346 assert(value_cmp_si(mpq_denref(cnt->count[3]), 1) == 0);
347 delete cnt;
348
349 cnt = new infinite_counter(2, 2);
350 cnt->init(ctx, 0);
351 set_from_string(n_coeff, "[-1/2 1/1 -1/3]");
352 set_from_string(n_power, "[[2 6][3 6]]");
353 d_power.SetDims(0, 2);
354 cnt->reduce(n_coeff, n_power, d_power);
355 assert(value_cmp_si(mpq_numref(cnt->count[0]), 1) == 0);
356 assert(value_cmp_si(mpq_denref(cnt->count[0]), 6) == 0);
357 assert(value_cmp_si(mpq_numref(cnt->count[1]), 0) == 0);
358 assert(value_cmp_si(mpq_denref(cnt->count[1]), 1) == 0);
359 assert(value_cmp_si(mpq_numref(cnt->count[2]), 0) == 0);
360 assert(value_cmp_si(mpq_denref(cnt->count[2]), 1) == 0);
361 delete cnt;
362
363 cnt = new infinite_counter(2, 2);
364 cnt->init(ctx, 0);
365 set_from_string(n_coeff, "[1/1]");
366 set_from_string(n_power, "[[0 11]]");
367 set_from_string(d_power, "[[0 1]]");
368 cnt->reduce(n_coeff, n_power, d_power);
369 assert(value_cmp_si(mpq_numref(cnt->count[1]), 0) != 0);
370 assert(value_cmp_si(mpq_numref(cnt->count[2]), 0) == 0);
371 assert(value_cmp_si(mpq_denref(cnt->count[2]), 1) == 0);
372 delete cnt;
373
374 Polyhedron_Free(ctx);
375
376 return 0;
377 }
378
379 static int test_series(struct barvinok_options *options)
380 {
381 Matrix *M;
382 Polyhedron *P, *C;
383 gen_fun *gf;
384
385 M = matrix_read_from_str(
386 "2 3\n"
387 "1 1 0\n"
388 "1 -1 10\n");
389 P = Constraints2Polyhedron(M, options->MaxRays);
390 Matrix_Free(M);
391 C = Universe_Polyhedron(1);
392 gf = barvinok_series_with_options(P, C, options);
393 Polyhedron_Free(P);
394 Polyhedron_Free(C);
395 gen_fun *sum = gf->summate(1, options);
396 delete gf;
397 delete sum;
398
399 return 0;
400 }
401
402 int test_todd(struct barvinok_options *options)
403 {
404 tcounter t(2, options->max_index);
405 assert(value_cmp_si(t.todd.coeff->p[0], 1) == 0);
406 assert(value_cmp_si(t.todd.coeff->p[1], -3) == 0);
407 assert(value_cmp_si(t.todd.coeff->p[2], 3) == 0);
408 assert(value_cmp_si(t.todd_denom->p[0], 1) == 0);
409 assert(value_cmp_si(t.todd_denom->p[1], 6) == 0);
410 assert(value_cmp_si(t.todd_denom->p[2], 36) == 0);
411
412 vec_ZZ lambda;
413 set_from_string(lambda, "[1 -1]");
414 zz2values(lambda, t.lambda->p);
415
416 mat_ZZ rays;
417 set_from_string(rays, "[[-1 0][-1 1]]");
418
419 QQ one(1, 1);
420
421 vec_ZZ v;
422 set_from_string(v, "[2 0 1]");
423 Vector *vertex = Vector_Alloc(3);
424 zz2values(v, vertex->p);
425
426 t.handle(rays, vertex->p, one, 1, options);
427 assert(value_cmp_si(mpq_numref(t.count), 71) == 0);
428 assert(value_cmp_si(mpq_denref(t.count), 24) == 0);
429
430 set_from_string(rays, "[[0 -1][1 -1]]");
431 set_from_string(v, "[0 2 1]");
432 zz2values(v, vertex->p);
433
434 t.handle(rays, vertex->p, one, 1, options);
435 assert(value_cmp_si(mpq_numref(t.count), 71) == 0);
436 assert(value_cmp_si(mpq_denref(t.count), 12) == 0);
437
438 set_from_string(rays, "[[1 0][0 1]]");
439 set_from_string(v, "[0 0 1]");
440 zz2values(v, vertex->p);
441
442 t.handle(rays, vertex->p, one, 1, options);
443 assert(value_cmp_si(mpq_numref(t.count), 6) == 0);
444 assert(value_cmp_si(mpq_denref(t.count), 1) == 0);
445
446 Vector_Free(vertex);
447 return 0;
448 }
449
450 int test_bernoulli(struct barvinok_options *options)
451 {
452 struct bernoulli_coef *bernoulli_coef;
453 struct poly_list *faulhaber, *bernoulli;
454 bernoulli_coef = bernoulli_coef_compute(2);
455 faulhaber = faulhaber_compute(4);
456 bernoulli_coef = bernoulli_coef_compute(8);
457 assert(value_cmp_si(bernoulli_coef->num->p[6], 1) == 0);
458 assert(value_cmp_si(bernoulli_coef->den->p[6], 42) == 0);
459 assert(value_cmp_si(faulhaber->poly[3]->p[0], 0) == 0);
460 assert(value_cmp_si(faulhaber->poly[3]->p[1], 0) == 0);
461 assert(value_cmp_si(faulhaber->poly[3]->p[2], 1) == 0);
462 assert(value_cmp_si(faulhaber->poly[3]->p[3], -2) == 0);
463 assert(value_cmp_si(faulhaber->poly[3]->p[4], 1) == 0);
464 assert(value_cmp_si(faulhaber->poly[3]->p[5], 4) == 0);
465
466 bernoulli = bernoulli_compute(6);
467 assert(value_cmp_si(bernoulli->poly[6]->p[0], 1) == 0);
468 assert(value_cmp_si(bernoulli->poly[6]->p[1], 0) == 0);
469 assert(value_cmp_si(bernoulli->poly[6]->p[2], -21) == 0);
470 assert(value_cmp_si(bernoulli->poly[6]->p[3], 0) == 0);
471 assert(value_cmp_si(bernoulli->poly[6]->p[4], 105) == 0);
472 assert(value_cmp_si(bernoulli->poly[6]->p[5], -126) == 0);
473 assert(value_cmp_si(bernoulli->poly[6]->p[6], 42) == 0);
474 assert(value_cmp_si(bernoulli->poly[6]->p[7], 42) == 0);
475
476 unsigned nvar, nparam;
477 const char **all_vars;
478 evalue *base, *sum1, *sum2;
479 base = evalue_read_from_str("(1 * n + 1)", NULL, &all_vars, &nvar, &nparam,
480 options->MaxRays);
481
482 sum1 = evalue_polynomial(faulhaber->poly[3], base);
483 Free_ParamNames(all_vars, nvar+nparam);
484
485 sum2 = evalue_read_from_str("(1/4 * n^4 + 1/2 * n^3 + 1/4 * n^2)",
486 NULL, &all_vars, &nvar, &nparam,
487 options->MaxRays);
488 Free_ParamNames(all_vars, nvar+nparam);
489 assert(eequal(sum1, sum2));
490 evalue_free(base);
491 evalue_free(sum1);
492 evalue_free(sum2);
493 return 0;
494 }
495
496 int test_bernoulli_sum(struct barvinok_options *options)
497 {
498 int summation = options->summation;
499 options->summation = BV_SUM_BERNOULLI;
500
501 unsigned nvar, nparam;
502 const char **all_vars;
503 evalue *e, *sum1, *sum2;
504 e = evalue_read_from_str("i + -1 >= 0\n -i + n >= 0\n\n 1 + (-1 *i) + i^2",
505 "i", &all_vars, &nvar, &nparam,
506 options->MaxRays);
507 Free_ParamNames(all_vars, nvar+nparam);
508
509 sum1 = barvinok_summate(e, 1, options);
510 sum2 = evalue_read_from_str("n -1 >= 0\n\n (1/3 * n^3 + 2/3 * n)",
511 NULL, &all_vars, &nvar, &nparam,
512 options->MaxRays);
513 Free_ParamNames(all_vars, nvar+nparam);
514 evalue_negate(sum1);
515 eadd(sum2, sum1);
516 reduce_evalue(sum1);
517 assert(EVALUE_IS_ZERO(*sum1));
518 evalue_free(e);
519 evalue_free(sum1);
520
521 e = evalue_read_from_str("-i + -1 >= 0\n i + n >= 0\n\n 1 + i + i^2",
522 "i", &all_vars, &nvar, &nparam,
523 options->MaxRays);
524 Free_ParamNames(all_vars, nvar+nparam);
525 sum1 = barvinok_summate(e, 1, options);
526 evalue_negate(sum1);
527 eadd(sum2, sum1);
528 reduce_evalue(sum1);
529 assert(EVALUE_IS_ZERO(*sum1));
530 evalue_free(e);
531 evalue_free(sum1);
532
533 evalue_free(sum2);
534
535 e = evalue_read_from_str("i + 4 >= 0\n -i + n >= 0\n\n i",
536 "i", &all_vars, &nvar, &nparam,
537 options->MaxRays);
538 Free_ParamNames(all_vars, nvar+nparam);
539 sum1 = barvinok_summate(e, 1, options);
540 sum2 = evalue_read_from_str("n + 0 >= 0\n\n (1/2 * n^2 + 1/2 * n + (-10))\n"
541 "n + 4 >= 0\n -n -1 >= 0\n\n (1/2 * n^2 + 1/2 * n + (-10))",
542 NULL, &all_vars, &nvar, &nparam,
543 options->MaxRays);
544 Free_ParamNames(all_vars, nvar+nparam);
545 evalue_negate(sum1);
546 eadd(sum2, sum1);
547 reduce_evalue(sum1);
548 assert(EVALUE_IS_ZERO(*sum1));
549 evalue_free(e);
550 evalue_free(sum1);
551 evalue_free(sum2);
552
553 e = evalue_read_from_str("i -5 >= 0\n -i + n >= 0\n j -1 >= 0\n i -j >= 0\n"
554 "k -1 >= 0\n j -k >= 0\n\n1",
555 "i,j,k", &all_vars, &nvar, &nparam,
556 options->MaxRays);
557 Free_ParamNames(all_vars, nvar+nparam);
558 sum1 = barvinok_summate(e, 3, options);
559 sum2 = evalue_read_from_str("n -5 >= 0\n\n"
560 "1/6 * n^3 + 1/2 * n^2 + 1/3 * n + -20",
561 NULL, &all_vars, &nvar, &nparam,
562 options->MaxRays);
563 Free_ParamNames(all_vars, nvar+nparam);
564 evalue_negate(sum1);
565 eadd(sum2, sum1);
566 reduce_evalue(sum1);
567 assert(EVALUE_IS_ZERO(*sum1));
568 evalue_free(e);
569 evalue_free(sum1);
570 evalue_free(sum2);
571
572 options->summation = summation;
573 return 0;
574 }
575
576 int test_hilbert(struct barvinok_options *options)
577 {
578 #ifdef USE_ZSOLVE
579 Matrix *M = matrix_read_from_str(
580 "2 4\n"
581 " 1 4 -3 0 \n"
582 " 1 3 2 0 \n");
583 Polyhedron *P = Constraints2Polyhedron(M, options->MaxRays);
584 Matrix_Free(M);
585
586 M = Cone_Hilbert_Basis(P, options->MaxRays);
587 assert(M->NbRows == 5);
588 assert(M->NbColumns == 3);
589 Matrix_Free(M);
590
591 M = Cone_Integer_Hull(P, NULL, 0, options);
592 assert(M->NbRows == 4);
593 assert(M->NbColumns == 3);
594 Matrix_Free(M);
595
596 Polyhedron_Free(P);
597 #endif
598 return 0;
599 }
600
601 int test_ilp(struct barvinok_options *options)
602 {
603 Matrix *M = matrix_read_from_str(
604 "2 4\n"
605 " 1 4 -3 0 \n"
606 " 1 3 2 0 \n");
607 Polyhedron *P = Constraints2Polyhedron(M, options->MaxRays);
608 Matrix_Free(M);
609 Vector *obj = Vector_Alloc(2);
610 value_set_si(obj->p[0], 7);
611 value_set_si(obj->p[1], -1);
612 Value min, max;
613 value_init(min);
614 value_init(max);
615
616 value_set_si(min, 1);
617 value_set_si(max, 17);
618 Vector *opt = Polyhedron_Integer_Minimum(P, obj->p, min, max,
619 NULL, 0, options);
620 assert(opt);
621 assert(value_cmp_si(opt->p[0], 1) == 0);
622 assert(value_cmp_si(opt->p[1], 1) == 0);
623 assert(value_cmp_si(opt->p[2], 1) == 0);
624 Vector_Free(opt);
625
626 value_clear(min);
627 value_clear(max);
628 Vector_Free(obj);
629 Polyhedron_Free(P);
630 return 0;
631 }
632
633 int test_hull(struct barvinok_options *options)
634 {
635 Matrix *M = matrix_read_from_str(
636 "4 4\n"
637 "1 32 -20 7\n"
638 "1 8 -44 187\n"
639 "1 -48 -4 285\n"
640 "1 8 68 -199\n");
641 Polyhedron *P = Constraints2Polyhedron(M, options->MaxRays);
642 Matrix_Free(M);
643
644 Matrix *hull = Polyhedron_Integer_Hull(P, options);
645 Polyhedron_Free(P);
646 assert(hull->NbRows == 4);
647 M = Matrix_Alloc(hull->NbRows, 1+hull->NbColumns);
648 for (int i = 0; i < hull->NbRows; ++i) {
649 value_set_si(M->p[i][0], 1);
650 Vector_Copy(hull->p[i], M->p[i]+1, hull->NbColumns);
651 }
652 Matrix_Free(hull);
653 Polyhedron *H = Constraints2Polyhedron(M, options->MaxRays);
654 Matrix_Free(M);
655
656 M = matrix_read_from_str(
657 "4 4\n"
658 "1 2 3 1 \n"
659 "1 3 4 1 \n"
660 "1 5 3 1 \n"
661 "1 5 5 1 \n");
662 P = Constraints2Polyhedron(M, options->MaxRays);
663 Matrix_Free(M);
664 assert(PolyhedronIncludes(P, H) && PolyhedronIncludes(H, P));
665 Polyhedron_Free(P);
666 Polyhedron_Free(H);
667
668 M = matrix_read_from_str(
669 "3 4\n"
670 "1 2 6 -3 \n"
671 "1 2 -6 3 \n"
672 "1 -2 0 3 \n");
673 P = Constraints2Polyhedron(M, options->MaxRays);
674 Matrix_Free(M);
675 assert(!emptyQ(P));
676 hull = Polyhedron_Integer_Hull(P, options);
677 Polyhedron_Free(P);
678 assert(hull->NbRows == 0);
679 Matrix_Free(hull);
680 return 0;
681 }
682
683 static int test_laurent(struct barvinok_options *options)
684 {
685 unsigned nvar, nparam;
686 const char **all_vars;
687 evalue *e, *sum, *res;
688
689 e = evalue_read_from_str(" x1 >= 0\n"
690 " x2 >= 0\n"
691 " -x1 -x2 + 2 >= 0\n"
692 "\n"
693 "(N + M * x2)\n",
694 "x1,x2", &all_vars, &nvar, &nparam,
695 options->MaxRays);
696 Free_ParamNames(all_vars, nvar+nparam);
697
698 int summation = options->summation;
699 options->summation = BV_SUM_LAURENT;
700 sum = barvinok_summate(e, nvar, options);
701 options->summation = summation;
702
703 res = evalue_read_from_str("(6 * N + 4 * M)",
704 "", &all_vars, &nvar, &nparam,
705 options->MaxRays);
706 Free_ParamNames(all_vars, nvar+nparam);
707
708 assert(value_zero_p(sum->d));
709 assert(sum->x.p->type == partition);
710 assert(sum->x.p->size == 2);
711
712 assert(eequal(res, &sum->x.p->arr[1]));
713
714 evalue_free(e);
715 evalue_free(sum);
716 evalue_free(res);
717 return 0;
718 }
719
720 /* Check that Polyhedron_Reduced_Basis produces a result
721 * of the expected dimensions (without crashing).
722 */
723 static int test_basis_reduction(struct barvinok_options *options)
724 {
725 Matrix *M;
726 Polyhedron *P;
727
728 M = matrix_read_from_str(
729 "4 4\n"
730 "1 1 0 0 \n"
731 "1 0 1 0 \n"
732 "1 -1 0 1 \n"
733 "1 0 -1 1 \n");
734 P = Constraints2Polyhedron(M, options->MaxRays);
735 Matrix_Free(M);
736
737 M = Polyhedron_Reduced_Basis(P, options);
738
739 assert(M);
740 assert(M->NbRows == 2);
741 assert(M->NbColumns == 2);
742
743 Polyhedron_Free(P);
744 Matrix_Free(M);
745
746 return 0;
747 }
748
749 int main(int argc, char **argv)
750 {
751 struct barvinok_options *options = barvinok_options_new_with_defaults();
752 test_equalities(options);
753 test_evalue_read(options);
754 test_eadd(options);
755 test_evalue(options);
756 test_substitute(options);
757 test_specialization(options);
758 test_lattice_points(options);
759 test_icounter(options);
760 test_infinite_counter(options);
761 test_series(options);
762 test_todd(options);
763 test_bernoulli(options);
764 test_bernoulli_sum(options);
765 test_hilbert(options);
766 test_ilp(options);
767 test_hull(options);
768 test_laurent(options);
769 test_basis_reduction(options);
770 barvinok_options_free(options);
771
772 return EXIT_SUCCESS;
773 }
lattice point enumeration library