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