use isl for argument parsing
[barvinok.git] / test_approx.c
blob4e37f0931ce3968496eda77e4a73c8a9b77ecbb2
1 #include <assert.h>
2 #include <limits.h>
3 #include <barvinok/barvinok.h>
4 #include "verify.h"
6 #ifdef HAVE_SYS_TIMES_H
8 #include <sys/times.h>
10 typedef clock_t my_clock_t;
12 static my_clock_t time_diff(struct tms *before, struct tms *after)
14 return after->tms_utime - before->tms_utime;
17 #else
19 typedef int my_clock_t;
21 struct tms { int dummy; };
22 static void times(struct tms* time)
25 static my_clock_t time_diff(struct tms *before, struct tms *after)
27 return 0;
30 #endif
32 struct {
33 int sign;
34 int method;
35 int flags;
36 } methods[] = {
37 { BV_APPROX_SIGN_NONE, BV_APPROX_NONE, 0 },
38 { BV_APPROX_SIGN_APPROX, BV_APPROX_BERNOULLI, 0 },
39 { BV_APPROX_SIGN_APPROX, BV_APPROX_DROP, 0 },
40 { BV_APPROX_SIGN_APPROX, BV_APPROX_VOLUME, BV_VOL_LIFT },
41 { BV_APPROX_SIGN_APPROX, BV_APPROX_VOLUME, BV_VOL_VERTEX },
42 { BV_APPROX_SIGN_APPROX, BV_APPROX_VOLUME, BV_VOL_BARYCENTER },
43 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, 0 },
44 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, BV_APPROX_SCALE_CHAMBER },
45 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST },
46 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
47 { BV_APPROX_SIGN_LOWER, BV_APPROX_DROP, 0 },
48 { BV_APPROX_SIGN_LOWER, BV_APPROX_VOLUME, BV_VOL_LIFT },
49 { BV_APPROX_SIGN_LOWER, BV_APPROX_VOLUME, BV_VOL_VERTEX },
50 { BV_APPROX_SIGN_LOWER, BV_APPROX_VOLUME, BV_VOL_BARYCENTER },
51 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, 0 },
52 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_CHAMBER },
53 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST },
54 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
55 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW },
56 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER},
57 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW },
58 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
59 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 },
60 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
61 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 },
62 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
63 { BV_APPROX_SIGN_UPPER, BV_APPROX_DROP, 0 },
64 { BV_APPROX_SIGN_UPPER, BV_APPROX_VOLUME, BV_VOL_LIFT },
65 { BV_APPROX_SIGN_UPPER, BV_APPROX_VOLUME, BV_VOL_VERTEX },
66 { BV_APPROX_SIGN_UPPER, BV_APPROX_VOLUME, BV_VOL_BARYCENTER },
67 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, 0 },
68 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_CHAMBER },
69 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST },
70 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
71 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW },
72 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
73 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW },
74 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
75 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 },
76 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
77 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 },
78 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
81 #define nr_methods (sizeof(methods)/sizeof(*methods))
83 struct options {
84 int quiet;
85 struct verify_options *verify;
88 struct isl_arg options_arg[] = {
89 ISL_ARG_CHILD(struct options, verify, NULL, verify_options_arg, NULL)
90 ISL_ARG_BOOL(struct options, quiet, 'q', "quiet", 0, NULL)
91 ISL_ARG_END
94 ISL_ARG_DEF(options, struct options, options_arg)
96 struct result_data {
97 Value n;
98 double RE_sum[nr_methods];
100 my_clock_t ticks[nr_methods];
101 size_t size[nr_methods];
104 void result_data_init(struct result_data *result)
106 int i;
107 for (i = 0; i < nr_methods; ++i) {
108 result->RE_sum[i] = 0;
109 result->ticks[i] = 0;
110 result->size[i] = 0;
112 value_init(result->n);
115 void result_data_clear(struct result_data *result)
117 int i;
118 value_clear(result->n);
121 void result_data_print(struct result_data *result, int n)
123 int i;
125 fprintf(stderr, "%d", result->ticks[0]/n);
126 for (i = 1; i < nr_methods; ++i)
127 fprintf(stderr, ", %d", result->ticks[i]/n);
128 fprintf(stderr, "\n");
130 fprintf(stderr, "%d", result->size[0]/n);
131 for (i = 1; i < nr_methods; ++i)
132 fprintf(stderr, ", %d", result->size[i]/n);
133 fprintf(stderr, "\n");
135 fprintf(stderr, "%g\n", VALUE_TO_DOUBLE(result->n));
136 fprintf(stderr, "%g", result->RE_sum[0]/VALUE_TO_DOUBLE(result->n));
137 for (i = 1; i < nr_methods; ++i)
138 fprintf(stderr, ", %g", result->RE_sum[i]/VALUE_TO_DOUBLE(result->n));
139 fprintf(stderr, "\n");
142 struct test_approx_data {
143 struct check_poly_data cp;
144 evalue **EP;
145 struct result_data *result;
148 static void eval(const evalue *EP, Value *z, int sign, Value *v)
150 evalue *res;
152 res = evalue_eval(EP, z);
153 if (sign == BV_APPROX_SIGN_LOWER)
154 mpz_cdiv_q(*v, res->x.n, res->d);
155 else if (sign == BV_APPROX_SIGN_UPPER)
156 mpz_fdiv_q(*v, res->x.n, res->d);
157 else if (sign == BV_APPROX_SIGN_APPROX)
158 mpz_tdiv_q(*v, res->x.n, res->d);
159 else {
160 assert(value_one_p(res->d));
161 value_assign(*v, res->x.n);
163 evalue_free(res);
166 static int test_approx(const struct check_poly_data *data, int nparam, Value *z,
167 const struct verify_options *options)
169 struct test_approx_data* ta_data = (struct test_approx_data*) data;
170 Value exact, approx;
171 int i;
173 value_init(exact);
174 value_init(approx);
176 eval(ta_data->EP[0], z, BV_APPROX_SIGN_NONE, &exact);
179 value_print(stderr, VALUE_FMT, exact);
182 value_increment(ta_data->result->n, ta_data->result->n);
183 for (i = 1; i < nr_methods; ++i) {
184 double error;
185 eval(ta_data->EP[i], z, methods[i].sign, &approx);
187 fprintf(stderr, ", ");
188 value_print(stderr, VALUE_FMT, approx);
190 if (methods[i].sign == BV_APPROX_SIGN_LOWER)
191 assert(value_le(approx, exact));
192 if (methods[i].sign == BV_APPROX_SIGN_UPPER)
193 assert(value_ge(approx, exact));
194 value_subtract(approx, approx, exact);
195 if (value_zero_p(exact))
196 error = abs(VALUE_TO_DOUBLE(approx));
197 else
198 error = abs(VALUE_TO_DOUBLE(approx)) / VALUE_TO_DOUBLE(exact);
199 ta_data->result->RE_sum[i] += error;
203 fprintf(stderr, "\n");
206 value_clear(exact);
207 value_clear(approx);
208 return 1;
211 static void test(Polyhedron *P, Polyhedron *C, evalue **EP,
212 struct result_data *result,
213 struct verify_options *options)
215 Polyhedron *CS;
216 Vector *p;
217 unsigned nparam = C->Dimension;
218 struct test_approx_data data;
220 CS = check_poly_context_scan(P, &C, C->Dimension, options);
222 p = Vector_Alloc(P->Dimension+2);
223 value_set_si(p->p[P->Dimension+1], 1);
225 check_poly_init(C, options);
227 data.cp.z = p->p;
228 data.cp.check = test_approx;
229 data.EP = EP;
230 data.result = result;
231 check_poly(CS, &data.cp, nparam, 0, p->p+P->Dimension-nparam+1,
232 options);
233 if (!options->print_all)
234 printf("\n");
236 Vector_Free(p);
237 if (CS) {
238 Domain_Free(CS);
239 Domain_Free(C);
243 void Matrix_File_Read_Input(FILE *in, Matrix *Mat)
245 Value *p;
246 int i,j,n;
247 char *c, s[1024],str[1024];
249 p = Mat->p_Init;
250 for (i=0;i<Mat->NbRows;i++) {
251 do {
252 c = fgets(s, 1024, in);
253 while(isspace(*c) && *c!='\n')
254 ++c;
255 } while(c && (*c =='#' || *c== '\n'));
257 if (!c) {
258 errormsg1( "Matrix_Read", "baddim", "not enough rows" );
259 break;
261 for (j=0;j<Mat->NbColumns;j++) {
262 if(!c || *c=='\n' || *c=='#') {
263 errormsg1("Matrix_Read", "baddim", "not enough columns");
264 break;
266 if (sscanf(c,"%s%n",str,&n) == 0) {
267 errormsg1( "Matrix_Read", "baddim", "not enough columns" );
268 break;
270 value_read(*(p++),str);
271 c += n;
274 } /* Matrix_Read_Input */
277 * Read the contents of the matrix 'Mat' from standard input.
278 * A '#' in the first column is a comment line
280 Matrix *Matrix_File_Read(FILE *in)
282 Matrix *Mat;
283 unsigned NbRows, NbColumns;
284 char s[1024];
286 fgets(s, 1024, in);
287 while ((*s=='#' || *s=='\n') ||
288 (sscanf(s, "%d %d", &NbRows, &NbColumns)<2))
289 fgets(s, 1024, in);
290 Mat = Matrix_Alloc(NbRows,NbColumns);
291 if(!Mat) {
292 errormsg1("Matrix_Read", "outofmem", "out of memory space");
293 return(NULL);
295 Matrix_File_Read_Input(in, Mat);
296 return Mat;
297 } /* Matrix_Read */
299 void handle(FILE *in, struct result_data *result, struct verify_options *options)
301 int i;
302 Polyhedron *A, *C;
303 Matrix *M;
304 const char **param_name;
305 evalue *EP[nr_methods];
307 M = Matrix_File_Read(in);
308 A = Constraints2Polyhedron(M, options->barvinok->MaxRays);
309 Matrix_Free(M);
310 M = Matrix_File_Read(in);
311 C = Constraints2Polyhedron(M, options->barvinok->MaxRays);
312 Matrix_Free(M);
313 param_name = Read_ParamNames(in, C->Dimension);
315 for (i = 0; i < nr_methods; ++i) {
316 struct tms st_cpu;
317 struct tms en_cpu;
318 options->barvinok->approx->approximation = methods[i].sign;
319 options->barvinok->approx->method = methods[i].method;
320 if (methods[i].method == BV_APPROX_SCALE)
321 options->barvinok->approx->scale_flags = methods[i].flags;
322 else if (methods[i].method == BV_APPROX_VOLUME)
323 options->barvinok->approx->volume_triangulate = methods[i].flags;
325 times(&st_cpu);
326 EP[i] = barvinok_enumerate_with_options(A, C, options->barvinok);
327 times(&en_cpu);
328 result->ticks[i] = time_diff(&en_cpu, &st_cpu);
330 print_evalue(stdout, EP[i], param_name);
333 for (i = 0; i < nr_methods; ++i)
334 result->size[i] = evalue_size(EP[i])/4;
335 test(A, C, EP, result, options);
336 for (i = 0; i < nr_methods; ++i)
337 evalue_free(EP[i]);
339 Free_ParamNames(param_name, C->Dimension);
340 Polyhedron_Free(A);
341 Polyhedron_Free(C);
344 int main(int argc, char **argv)
346 char path[PATH_MAX+1];
347 struct result_data all_result;
348 int n = 0;
349 struct options *options = options_new_with_defaults();
351 argc = options_parse(options, argc, argv, ISL_ARG_ALL);
353 if (options->verify->M == INT_MIN)
354 options->verify->M = 100;
355 if (options->verify->m == INT_MAX)
356 options->verify->m = -100;
358 result_data_init(&all_result);
360 while (fgets(path, sizeof(path), stdin)) {
361 struct result_data result;
362 FILE *in;
363 int i;
365 ++n;
366 result_data_init(&result);
367 fprintf(stderr, "%s", path);
368 *strchr(path, '\n') = '\0';
369 in = fopen(path, "r");
370 assert(in);
371 handle(in, &result, options->verify);
372 fclose(in);
374 if (!options->quiet)
375 result_data_print(&result, 1);
377 value_addto(all_result.n, all_result.n, result.n);
378 for (i = 0; i < nr_methods; ++i) {
379 all_result.RE_sum[i] += result.RE_sum[i];
380 all_result.ticks[i] += result.ticks[i];
381 all_result.size[i] += result.size[i];
384 result_data_clear(&result);
386 if (!options->quiet) {
387 fprintf(stderr, "average\n");
388 result_data_print(&all_result, n);
392 result_data_clear(&all_result);
394 options_free(options);
396 return 0;