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