?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
evalue.c: reorder_terms: fix typo
[barvinok.git] / test_approx.c
blob93f32d6a9325065c8638cd331a8f3fb3a7ad0e07
1 #include <limits.h>
2 #include <sys/times.h>
3 #include <barvinok/barvinok.h>
4 #include "verify.h"
5 #include "argp.h"
6 #include "progname.h"
7
8 struct {
9 int sign;
10 int method;
11 int flags;
12 } methods[] = {
13 { BV_APPROX_SIGN_NONE, BV_APPROX_NONE, 0 },
14 { BV_APPROX_SIGN_APPROX, BV_APPROX_BERNOULLI, 0 },
15 { BV_APPROX_SIGN_APPROX, BV_APPROX_DROP, 0 },
16 { BV_APPROX_SIGN_APPROX, BV_APPROX_VOLUME, BV_VOL_LIFT },
17 { BV_APPROX_SIGN_APPROX, BV_APPROX_VOLUME, BV_VOL_VERTEX },
18 { BV_APPROX_SIGN_APPROX, BV_APPROX_VOLUME, BV_VOL_BARYCENTER },
19 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, 0 },
20 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, BV_APPROX_SCALE_CHAMBER },
21 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST },
22 { BV_APPROX_SIGN_APPROX, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
23 { BV_APPROX_SIGN_LOWER, BV_APPROX_DROP, 0 },
24 { BV_APPROX_SIGN_LOWER, BV_APPROX_VOLUME, BV_VOL_LIFT },
25 { BV_APPROX_SIGN_LOWER, BV_APPROX_VOLUME, BV_VOL_VERTEX },
26 { BV_APPROX_SIGN_LOWER, BV_APPROX_VOLUME, BV_VOL_BARYCENTER },
27 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, 0 },
28 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_CHAMBER },
29 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST },
30 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
31 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW },
32 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER},
33 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW },
34 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
35 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 },
36 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
37 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 },
38 { BV_APPROX_SIGN_LOWER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
39 { BV_APPROX_SIGN_UPPER, BV_APPROX_DROP, 0 },
40 { BV_APPROX_SIGN_UPPER, BV_APPROX_VOLUME, BV_VOL_LIFT },
41 { BV_APPROX_SIGN_UPPER, BV_APPROX_VOLUME, BV_VOL_VERTEX },
42 { BV_APPROX_SIGN_UPPER, BV_APPROX_VOLUME, BV_VOL_BARYCENTER },
43 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, 0 },
44 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_CHAMBER },
45 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST },
46 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
47 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW },
48 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
49 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW },
50 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
51 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 },
52 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
53 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 },
54 { BV_APPROX_SIGN_UPPER, BV_APPROX_SCALE, BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
55 };
56
57 #define nr_methods (sizeof(methods)/sizeof(*methods))
58
59 struct argp_option argp_options[] = {
60 { "quiet", 'q' },
61 { 0 }
62 };
63
64 struct options {
65 int quiet;
66 struct verify_options verify;
67 };
68
69 static error_t parse_opt(int key, char *arg, struct argp_state *state)
70 {
71 struct options *options = (struct options*) state->input;
72
73 switch (key) {
74 case ARGP_KEY_INIT:
75 state->child_inputs[0] = options->verify.barvinok;
76 state->child_inputs[1] = &options->verify;
77 options->quiet = 0;
78 break;
79 case 'q':
80 options->quiet = 1;
81 break;
82 default:
83 return ARGP_ERR_UNKNOWN;
84 }
85 return 0;
86 }
87
88 struct result_data {
89 Value n;
90 double RE_sum[nr_methods];
91
92 clock_t ticks[nr_methods];
93 size_t size[nr_methods];
94 };
95
96 void result_data_init(struct result_data *result)
97 {
98 int i;
99 for (i = 0; i < nr_methods; ++i) {
100 result->RE_sum[i] = 0;
101 result->ticks[i] = 0;
102 result->size[i] = 0;
103 }
104 value_init(result->n);
105 }
106
107 void result_data_clear(struct result_data *result)
108 {
109 int i;
110 value_clear(result->n);
111 }
112
113 void result_data_print(struct result_data *result, int n)
114 {
115 int i;
116
117 fprintf(stderr, "%d", result->ticks[0]/n);
118 for (i = 1; i < nr_methods; ++i)
119 fprintf(stderr, ", %d", result->ticks[i]/n);
120 fprintf(stderr, "\n");
121
122 fprintf(stderr, "%d", result->size[0]/n);
123 for (i = 1; i < nr_methods; ++i)
124 fprintf(stderr, ", %d", result->size[i]/n);
125 fprintf(stderr, "\n");
126
127 fprintf(stderr, "%g\n", VALUE_TO_DOUBLE(result->n));
128 fprintf(stderr, "%g", result->RE_sum[0]/VALUE_TO_DOUBLE(result->n));
129 for (i = 1; i < nr_methods; ++i)
130 fprintf(stderr, ", %g", result->RE_sum[i]/VALUE_TO_DOUBLE(result->n));
131 fprintf(stderr, "\n");
132 }
133
134 struct test_approx_data {
135 struct check_poly_data cp;
136 evalue **EP;
137 struct result_data *result;
138 };
139
140 static void eval(const evalue *EP, Value *z, int sign, Value *v)
141 {
142 evalue *res;
143
144 res = evalue_eval(EP, z);
145 if (sign == BV_APPROX_SIGN_LOWER)
146 mpz_cdiv_q(*v, res->x.n, res->d);
147 else if (sign == BV_APPROX_SIGN_UPPER)
148 mpz_fdiv_q(*v, res->x.n, res->d);
149 else if (sign == BV_APPROX_SIGN_APPROX)
150 mpz_tdiv_q(*v, res->x.n, res->d);
151 else {
152 assert(value_one_p(res->d));
153 value_assign(*v, res->x.n);
154 }
155 free_evalue_refs(res);
156 free(res);
157 }
158
159 static int test_approx(const struct check_poly_data *data, int nparam, Value *z,
160 const struct verify_options *options)
161 {
162 struct test_approx_data* ta_data = (struct test_approx_data*) data;
163 Value exact, approx;
164 int i;
165
166 value_init(exact);
167 value_init(approx);
168
169 eval(ta_data->EP[0], z, BV_APPROX_SIGN_NONE, &exact);
170
171 /*
172 value_print(stderr, VALUE_FMT, exact);
173 */
174
175 value_increment(ta_data->result->n, ta_data->result->n);
176 for (i = 1; i < nr_methods; ++i) {
177 double error;
178 eval(ta_data->EP[i], z, methods[i].sign, &approx);
179 /*
180 fprintf(stderr, ", ");
181 value_print(stderr, VALUE_FMT, approx);
182 */
183 if (methods[i].sign == BV_APPROX_SIGN_LOWER)
184 assert(value_le(approx, exact));
185 if (methods[i].sign == BV_APPROX_SIGN_UPPER)
186 assert(value_ge(approx, exact));
187 value_subtract(approx, approx, exact);
188 if (value_zero_p(exact))
189 error = abs(VALUE_TO_DOUBLE(approx));
190 else
191 error = abs(VALUE_TO_DOUBLE(approx)) / VALUE_TO_DOUBLE(exact);
192 ta_data->result->RE_sum[i] += error;
193 }
194
195 /*
196 fprintf(stderr, "\n");
197 */
198
199 value_clear(exact);
200 value_clear(approx);
201 return 1;
202 }
203
204 static void test(Polyhedron *P, Polyhedron *C, evalue **EP,
205 struct result_data *result,
206 struct verify_options *options)
207 {
208 Polyhedron *CS;
209 Vector *p;
210 unsigned nparam = C->Dimension;
211 struct test_approx_data data;
212
213 CS = check_poly_context_scan(P, &C, C->Dimension, options);
214
215 p = Vector_Alloc(P->Dimension+2);
216 value_set_si(p->p[P->Dimension+1], 1);
217
218 check_poly_init(C, options);
219
220 data.cp.z = p->p;
221 data.cp.check = test_approx;
222 data.EP = EP;
223 data.result = result;
224 check_poly(CS, &data.cp, nparam, 0, p->p+P->Dimension-nparam+1,
225 options);
226 if (!options->print_all)
227 printf("\n");
228
229 Vector_Free(p);
230 if (CS) {
231 Domain_Free(CS);
232 Domain_Free(C);
233 }
234 }
235
236 void Matrix_File_Read_Input(FILE *in, Matrix *Mat)
237 {
238 Value *p;
239 int i,j,n;
240 char *c, s[1024],str[1024];
241
242 p = Mat->p_Init;
243 for (i=0;i<Mat->NbRows;i++) {
244 do {
245 c = fgets(s, 1024, in);
246 while(isspace(*c) && *c!='\n')
247 ++c;
248 } while(c && (*c =='#' || *c== '\n'));
249
250 if (!c) {
251 errormsg1( "Matrix_Read", "baddim", "not enough rows" );
252 break;
253 }
254 for (j=0;j<Mat->NbColumns;j++) {
255 if(!c || *c=='\n' || *c=='#') {
256 errormsg1("Matrix_Read", "baddim", "not enough columns");
257 break;
258 }
259 if (sscanf(c,"%s%n",str,&n) == 0) {
260 errormsg1( "Matrix_Read", "baddim", "not enough columns" );
261 break;
262 }
263 value_read(*(p++),str);
264 c += n;
265 }
266 }
267 } /* Matrix_Read_Input */
268
269 /*
270 * Read the contents of the matrix 'Mat' from standard input.
271 * A '#' in the first column is a comment line
272 */
273 Matrix *Matrix_File_Read(FILE *in)
274 {
275 Matrix *Mat;
276 unsigned NbRows, NbColumns;
277 char s[1024];
278
279 fgets(s, 1024, in);
280 while ((*s=='#' || *s=='\n') ||
281 (sscanf(s, "%d %d", &NbRows, &NbColumns)<2))
282 fgets(s, 1024, in);
283 Mat = Matrix_Alloc(NbRows,NbColumns);
284 if(!Mat) {
285 errormsg1("Matrix_Read", "outofmem", "out of memory space");
286 return(NULL);
287 }
288 Matrix_File_Read_Input(in, Mat);
289 return Mat;
290 } /* Matrix_Read */
291
292 void handle(FILE *in, struct result_data *result, struct verify_options *options)
293 {
294 int i;
295 Polyhedron *A, *C;
296 Matrix *M;
297 char **param_name;
298 evalue *EP[nr_methods];
299
300 M = Matrix_File_Read(in);
301 A = Constraints2Polyhedron(M, options->barvinok->MaxRays);
302 Matrix_Free(M);
303 M = Matrix_File_Read(in);
304 C = Constraints2Polyhedron(M, options->barvinok->MaxRays);
305 Matrix_Free(M);
306 param_name = Read_ParamNames(in, C->Dimension);
307
308 for (i = 0; i < nr_methods; ++i) {
309 struct tms st_cpu;
310 struct tms en_cpu;
311 options->barvinok->polynomial_approximation = methods[i].sign;
312 options->barvinok->approximation_method = methods[i].method;
313 if (methods[i].method == BV_APPROX_SCALE)
314 options->barvinok->scale_flags = methods[i].flags;
315 else if (methods[i].method == BV_APPROX_VOLUME)
316 options->barvinok->volume_triangulate = methods[i].flags;
317
318 times(&st_cpu);
319 EP[i] = barvinok_enumerate_with_options(A, C, options->barvinok);
320 times(&en_cpu);
321 result->ticks[i] = en_cpu.tms_utime - st_cpu.tms_utime;
322 /*
323 print_evalue(stdout, EP[i], param_name);
324 */
325 }
326 for (i = 0; i < nr_methods; ++i)
327 result->size[i] = evalue_size(EP[i])/4;
328 test(A, C, EP, result, options);
329 for (i = 0; i < nr_methods; ++i) {
330 free_evalue_refs(EP[i]);
331 free(EP[i]);
332 }
333
334 Free_ParamNames(param_name, C->Dimension);
335 Polyhedron_Free(A);
336 Polyhedron_Free(C);
337 }
338
339 int main(int argc, char **argv)
340 {
341 struct barvinok_options *bv_options = barvinok_options_new_with_defaults();
342 char path[PATH_MAX+1];
343 struct result_data all_result;
344 int n = 0;
345 static struct argp_child argp_children[] = {
346 { &barvinok_argp, 0, 0, 0 },
347 { &verify_argp, 0, "verification", BV_GRP_LAST+1 },
348 { 0 }
349 };
350 static struct argp argp = { argp_options, parse_opt, 0, 0, argp_children };
351 struct options options;
352
353 options.verify.barvinok = bv_options;
354 set_program_name(argv[0]);
355 argp_parse(&argp, argc, argv, 0, 0, &options);
356
357 if (options.verify.M == INT_MIN)
358 options.verify.M = 100;
359 if (options.verify.m == INT_MAX)
360 options.verify.m = -100;
361
362 result_data_init(&all_result);
363
364 while (fgets(path, sizeof(path), stdin)) {
365 struct result_data result;
366 FILE *in;
367 int i;
368
369 ++n;
370 result_data_init(&result);
371 fprintf(stderr, "%s", path);
372 *strchr(path, '\n') = '\0';
373 in = fopen(path, "r");
374 assert(in);
375 handle(in, &result, &options.verify);
376 fclose(in);
377
378 if (!options.quiet)
379 result_data_print(&result, 1);
380
381 value_addto(all_result.n, all_result.n, result.n);
382 for (i = 0; i < nr_methods; ++i) {
383 all_result.RE_sum[i] += result.RE_sum[i];
384 all_result.ticks[i] += result.ticks[i];
385 all_result.size[i] += result.size[i];
386 }
387
388 result_data_clear(&result);
389
390 if (!options.quiet) {
391 fprintf(stderr, "average\n");
392 result_data_print(&all_result, n);
393 }
394 }
395
396 result_data_clear(&all_result);
397
398 barvinok_options_free(bv_options);
399
400 return 0;
401 }
lattice point enumeration library