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