update pet for explicitly linking in gmp
[barvinok.git] / verify.c
blob285fe3a4be36deee4b23a2679d2eb2bd56ed70e1
1 #include <assert.h>
2 #include <stdlib.h>
3 #include <barvinok/options.h>
4 #include <barvinok/util.h>
5 #include "verify.h"
7 #define ALLOCN(type,n) (type*)malloc((n) * sizeof(type))
9 /* RANGE : normal range for evalutations (-RANGE -> RANGE) */
10 #define RANGE 50
12 /* SRANGE : small range for evalutations */
13 #define SRANGE 15
15 /* if dimension >= BIDDIM, use SRANGE */
16 #define BIGDIM 5
18 /* VSRANGE : very small range for evalutations */
19 #define VSRANGE 5
21 /* if dimension >= VBIDDIM, use VSRANGE */
22 #define VBIGDIM 8
24 static int set_m(void *opt, long val)
26 struct verify_options *options = (struct verify_options *)opt;
27 options->m = val;
28 options->verify = 1;
29 return 0;
32 static int set_M(void *opt, long val)
34 struct verify_options *options = (struct verify_options *)opt;
35 options->M = val;
36 options->verify = 1;
37 return 0;
40 static int set_r(void *opt, long val)
42 struct verify_options *options = (struct verify_options *)opt;
43 options->r = val;
44 options->m = -val;
45 options->M = val;
46 options->verify = 1;
47 return 0;
50 struct isl_arg verify_options_arg[] = {
51 ISL_ARG_CHILD(struct verify_options, barvinok, NULL, barvinok_options_arg, NULL)
52 ISL_ARG_BOOL(struct verify_options, verify, 'T', "verify", 0, NULL)
53 ISL_ARG_BOOL(struct verify_options, exact, 'E', "exact", 0, NULL)
54 ISL_ARG_BOOL(struct verify_options, print_all, 'A', "print-all", 0, NULL)
55 ISL_ARG_BOOL(struct verify_options, continue_on_error, 'C',
56 "continue-on-error", 0, NULL)
57 ISL_ARG_USER_LONG(struct verify_options, m, 'm', "min", set_m, INT_MAX, NULL)
58 ISL_ARG_USER_LONG(struct verify_options, M, 'M', "max", set_M, INT_MIN, NULL)
59 ISL_ARG_USER_LONG(struct verify_options, r, 'r', "range", set_r, -1, NULL)
60 ISL_ARG_END
63 void verify_options_set_range(struct verify_options *options, int dim)
65 int r;
67 if (dim >= VBIGDIM)
68 r = VSRANGE;
69 else if (dim >= BIGDIM)
70 r = SRANGE;
71 else
72 r = RANGE;
73 /* If the user didn't set m or M, then we try to adjust the window
74 * to the context in check_poly_context_scan.
76 if (options->m == INT_MAX && options->M == INT_MIN)
77 options->r = r;
78 if (options->M == INT_MIN)
79 options->M = r;
80 if (options->m == INT_MAX)
81 options->m = -r;
83 if (options->verify && options->m > options->M) {
84 fprintf(stderr,"Nothing to do: min > max !\n");
85 exit(0);
89 static Polyhedron *project_on(Polyhedron *P, int i)
91 unsigned dim = P->Dimension;
92 Matrix *T;
93 Polyhedron *I;
95 if (dim == 1)
96 return Polyhedron_Copy(P);
98 T = Matrix_Alloc(2, dim+1);
99 value_set_si(T->p[0][i], 1);
100 value_set_si(T->p[1][dim], 1);
101 I = Polyhedron_Image(P, T, P->NbConstraints);
102 Matrix_Free(T);
103 return I;
106 static void set_bounds(Polyhedron *C, Value **rows, int i, unsigned nparam,
107 const struct verify_options *options)
109 Value min;
110 Value max;
112 value_init(min);
113 value_init(max);
114 value_set_si(min, options->m);
115 value_set_si(max, options->M);
117 if (options->r > 0) {
118 Polyhedron *I = project_on(C, i);
119 line_minmax(I, &min, &max);
120 if (value_cmp_si(min, options->M) >= 0)
121 value_add_int(max, min, options->r);
122 else if (value_cmp_si(max, options->m) <= 0)
123 value_sub_int(min, max, options->r);
124 else {
125 value_set_si(min, options->m);
126 value_set_si(max, options->M);
130 value_set_si(rows[0][0], 1);
131 value_set_si(rows[0][1+i], 1);
132 value_oppose(rows[0][1+nparam], min);
133 value_set_si(rows[1][0], 1);
134 value_set_si(rows[1][1+i], -1);
135 value_assign(rows[1][1+nparam], max);
137 value_clear(min);
138 value_clear(max);
141 Polyhedron *check_poly_context_scan(Polyhedron *P, Polyhedron **C,
142 unsigned nparam,
143 const struct verify_options *options)
145 int i;
146 Matrix *MM;
147 Polyhedron *CC, *CC2, *CS, *U;
148 unsigned MaxRays = options->barvinok->MaxRays;
150 if (nparam <= 0)
151 return NULL;
153 if (!P)
154 CC = *C;
155 else {
156 CC = Polyhedron_Project(P, nparam);
157 if (*C) {
158 CC2 = DomainIntersection(*C, CC, MaxRays);
159 Domain_Free(CC);
160 CC = CC2;
164 /* Intersect context with range */
165 MM = Matrix_Alloc(2*nparam, nparam+2);
166 for (i = 0; i < nparam; ++i)
167 set_bounds(CC, &MM->p[2*i], i, nparam, options);
168 CC2 = AddConstraints(MM->p[0], 2*nparam, CC, options->barvinok->MaxRays);
169 if (CC != *C)
170 Domain_Free(CC);
171 CC = CC2;
172 U = Universe_Polyhedron(0);
173 CS = Polyhedron_Scan(CC, U, MaxRays & POL_NO_DUAL ? 0 : MaxRays);
174 Polyhedron_Free(U);
175 *C = CC;
176 Matrix_Free(MM);
177 return CS;
180 void check_poly_init(Polyhedron *C, struct verify_options *options)
182 int d, i;
184 if (options->print_all)
185 return;
186 if (C->Dimension <= 0)
187 return;
189 d = options->M - options->m;
190 if (d > 80)
191 options->st = 1+d/80;
192 else
193 options->st = 1;
194 for (i = options->m; i <= options->M; i += options->st)
195 printf(".");
196 printf( "\r" );
197 fflush(stdout);
200 static void print_rational(FILE *out, Value n, Value d)
202 value_print(out, VALUE_FMT, n);
203 if (value_notone_p(d)) {
204 fprintf(out, "/");
205 value_print(out, VALUE_FMT, d);
209 void check_poly_print(int ok, int nparam, Value *z,
210 Value want_n, Value want_d,
211 Value got_n, Value got_d,
212 const char *op, const char *check,
213 const char *long_op,
214 const struct verify_options *options)
216 int k;
218 if (options->print_all) {
219 printf("%s(", op);
220 value_print(stdout, VALUE_FMT, z[0]);
221 for (k = 1; k < nparam; ++k) {
222 printf(", ");
223 value_print(stdout, VALUE_FMT, z[k]);
225 printf(") = ");
226 print_rational(stdout, got_n, got_d);
227 printf(", %s = ", check);
228 print_rational(stdout, want_n, want_d);
229 printf(". ");
232 if (!ok) {
233 printf("\n");
234 fflush(stdout);
235 fprintf(stderr, "Error !\n");
236 fprintf(stderr, "%s(", op);
237 value_print(stderr, VALUE_FMT, z[0]);
238 for (k = 1; k < nparam; ++k) {
239 fprintf(stderr,", ");
240 value_print(stderr, VALUE_FMT, z[k]);
242 fprintf(stderr, ") should be ");
243 print_rational(stderr, want_n, want_d);
244 fprintf(stderr, ", while %s gives ", long_op);
245 print_rational(stderr, got_n, got_d);
246 fprintf(stderr, ".\n");
247 } else if (options->print_all)
248 printf("OK.\n");
251 /****************************************************/
252 /* function check_poly : */
253 /* scans the parameter space from Min to Max (all */
254 /* directions). Computes the number of points in */
255 /* the polytope using both methods, and compare them*/
256 /* returns 1 on success */
257 /****************************************************/
259 int check_poly(Polyhedron *CS, const struct check_poly_data *data,
260 int nparam, int pos, Value *z,
261 const struct verify_options *options)
263 if (pos == nparam) {
264 if (!data->check(data, nparam, z, options) && !options->continue_on_error)
265 return 0;
266 } else {
267 Value LB, UB;
268 int ok;
269 value_init(LB);
270 value_init(UB);
271 ok = !(lower_upper_bounds(1+pos, CS, z-1, &LB, &UB));
272 assert(ok);
273 for (; value_le(LB, UB); value_increment(LB, LB)) {
274 if (!options->print_all) {
275 int k = VALUE_TO_INT(LB);
276 if (!pos && !(k % options->st)) {
277 printf("o");
278 fflush(stdout);
282 value_assign(z[pos], LB);
283 if (!check_poly(CS->next, data, nparam, pos+1, z, options)) {
284 value_clear(LB);
285 value_clear(UB);
286 return 0;
289 value_set_si(z[pos], 0);
290 value_clear(LB);
291 value_clear(UB);
293 return 1;
294 } /* check_poly */
296 void check_EP_set_scan(struct check_EP_data *data, Polyhedron *C,
297 unsigned MaxRays)
299 const evalue *EP = data->EP;
300 int i;
301 int n_S = 0;
303 for (i = 0; i < EP->x.p->size/2; ++i) {
304 Polyhedron *A = EVALUE_DOMAIN(EP->x.p->arr[2*i]);
305 for (; A; A = A->next)
306 ++n_S;
309 data->n_S = n_S;
310 data->S = ALLOCN(Polyhedron *, n_S);
311 n_S = 0;
312 for (i = 0; i < EP->x.p->size/2; ++i) {
313 Polyhedron *A = EVALUE_DOMAIN(EP->x.p->arr[2*i]);
314 for (; A; A = A->next) {
315 Polyhedron *next = A->next;
316 A->next = NULL;
317 data->S[n_S++] = Polyhedron_Scan(A, C,
318 MaxRays & POL_NO_DUAL ? 0 : MaxRays);
319 A->next = next;
324 void check_EP_clear_scan(struct check_EP_data *data)
326 int i;
328 for (i = 0; i < data->n_S; ++i)
329 Domain_Free(data->S[i]);
330 free(data->S);
333 static int check_EP_on_poly(Polyhedron *P,
334 struct check_EP_data *data,
335 unsigned nvar, unsigned nparam,
336 struct verify_options *options)
338 Polyhedron *CS;
339 unsigned MaxRays = options->barvinok->MaxRays;
340 int ok = 1;
341 const evalue *EP = data->EP;
343 CS = check_poly_context_scan(NULL, &P, P->Dimension, options);
345 check_poly_init(P, options);
347 if (!(CS && emptyQ2(CS))) {
348 check_EP_set_scan(data, P, MaxRays);
349 ok = check_poly(CS, &data->cp, nparam, 0, data->cp.z+1+nvar, options);
350 check_EP_clear_scan(data);
353 if (!options->print_all)
354 printf("\n");
356 if (CS) {
357 Domain_Free(CS);
358 Domain_Free(P);
361 return ok;
365 * Project on final dim dimensions
367 Polyhedron *DomainProject(Polyhedron *D, unsigned dim, unsigned MaxRays)
369 Polyhedron *P;
370 Polyhedron *R;
372 R = Polyhedron_Project(D, dim);
373 for (P = D->next; P; P = P->next) {
374 Polyhedron *R2 = Polyhedron_Project(P, dim);
375 Polyhedron *R3 = DomainUnion(R, R2, MaxRays);
376 Polyhedron_Free(R2);
377 Domain_Free(R);
378 R = R3;
380 return R;
383 static Polyhedron *evalue_parameter_domain(const evalue *e, unsigned nparam,
384 unsigned MaxRays)
386 int i;
387 Polyhedron *U = NULL;
389 if (EVALUE_IS_ZERO(*e))
390 return Universe_Polyhedron(0);
392 assert(value_zero_p(e->d));
393 assert(e->x.p->type == partition);
394 assert(e->x.p->size >= 2);
396 for (i = 0; i < e->x.p->size/2; ++i) {
397 Polyhedron *D = EVALUE_DOMAIN(e->x.p->arr[2*i]);
398 Polyhedron *P = DomainProject(D, nparam, MaxRays);
399 if (!U) {
400 U = P;
401 } else {
402 Polyhedron *T = U;
403 U = DomainUnion(U, P, MaxRays);
404 Domain_Free(P);
405 Domain_Free(T);
408 return U;
411 int check_EP(struct check_EP_data *data, unsigned nvar, unsigned nparam,
412 struct verify_options *options)
414 Vector *p;
415 int i;
416 int ok = 1;
417 Polyhedron *D, *P;
419 p = Vector_Alloc(nvar+nparam+2);
420 value_set_si(p->p[nvar+nparam+1], 1);
422 data->cp.z = p->p;
424 D = evalue_parameter_domain(data->EP, nparam, options->barvinok->MaxRays);
426 for (P = D; P; P = P->next) {
427 ok = check_EP_on_poly(P, data, nvar, nparam, options);
428 if (!ok && !options->continue_on_error)
429 break;
432 Domain_Free(D);
433 Vector_Free(p);
435 return ok;
438 __isl_give isl_set *verify_context_set_bounds(__isl_take isl_set *set,
439 const struct verify_options *options)
441 int i;
442 unsigned nparam;
443 isl_point *pt, *pt2;
444 isl_set *box;
446 nparam = isl_set_dim(set, isl_dim_param);
448 if (options->r > 0) {
449 pt = isl_set_sample_point(isl_set_copy(set));
450 pt2 = isl_point_copy(pt);
452 for (i = 0; i < nparam; ++i) {
453 pt = isl_point_add_ui(pt, isl_dim_param, i, options->r);
454 pt2 = isl_point_sub_ui(pt2, isl_dim_param, i, options->r);
456 } else {
457 isl_int v;
459 isl_int_init(v);
460 pt = isl_point_zero(isl_set_get_space(set));
461 isl_int_set_si(v, options->m);
462 for (i = 0; i < nparam; ++i)
463 pt = isl_point_set_coordinate(pt, isl_dim_param, i, v);
465 pt2 = isl_point_zero(isl_set_get_space(set));
466 isl_int_set_si(v, options->M);
467 for (i = 0; i < nparam; ++i)
468 pt2 = isl_point_set_coordinate(pt2, isl_dim_param, i, v);
470 isl_int_clear(v);
473 box = isl_set_box_from_points(pt, pt2);
475 return isl_set_intersect(set, box);
478 int verify_point_data_init(struct verify_point_data *vpd,
479 __isl_keep isl_set *context)
481 isl_int v;
482 int i;
483 int r;
485 isl_int_init(v);
486 r = isl_set_count(context, &v);
487 vpd->n = isl_int_cmp_si(v, 200) < 0 ? isl_int_get_si(v) : 200;
488 isl_int_clear(v);
490 if (!vpd->options->print_all) {
491 vpd->s = vpd->n < 80 ? 1 : 1 + vpd->n/80;
492 for (i = 0; i < vpd->n; i += vpd->s)
493 printf(".");
494 printf("\r");
495 fflush(stdout);
498 vpd->error = r < 0 ? -1 : 0;
500 return r;
503 void verify_point_data_fini(struct verify_point_data *vpd)
505 if (!vpd->options->print_all)
506 printf("\n");
508 if (vpd->error)
509 fprintf(stderr, "Check failed !\n");