1 #include <assert.h>
2 #include <limits.h>
3 #include <math.h>
4 #include <string.h>
5 #include <isl/ctx.h>
6 #include <isl/stream.h>
7 #include <isl/options.h>
8 #include <isl/val.h>
9 #include <isl/point.h>
10 #include <isl/set.h>
11 #include <isl/polynomial.h>
12 #include <isl/printer.h>
13 #include <barvinok/options.h>
14 #include <barvinok/util.h>
15 #include "verify.h"
16 #include "config.h"
17 
18 #ifdef HAVE_SYS_TIMES_H
19 
20 #include <sys/times.h>
21 
22 typedef clock_t		my_clock_t;
23 
time_diff(struct tms * before,struct tms * after)24 static my_clock_t time_diff(struct tms *before, struct tms *after)
25 {
26 	return after->tms_utime - before->tms_utime;
27 }
28 
29 #else
30 
31 typedef int		my_clock_t;
32 
33 struct tms {};
times(struct tms * time)34 static void times(struct tms* time)
35 {
36 }
time_diff(struct tms * before,struct tms * after)37 static my_clock_t time_diff(struct tms *before, struct tms *after)
38 {
39 	return 0;
40 }
41 
42 #endif
43 
44 static struct {
45     int	    method;
46 } methods[] = {
47     { ISL_BOUND_BERNSTEIN },
48     { ISL_BOUND_RANGE },
49 };
50 
51 #define nr_methods (sizeof(methods)/sizeof(*methods))
52 
53 struct options {
54     struct verify_options    *verify;
55     int quiet;
56 };
57 
58 ISL_ARGS_START(struct options, options_args)
59 ISL_ARG_CHILD(struct options, verify, NULL, &verify_options_args, NULL)
60 ISL_ARG_BOOL(struct options, quiet, 'q', "quiet", 0, NULL)
61 ISL_ARGS_END
62 
63 ISL_ARG_DEF(options, struct options, options_args)
64 
65 struct result_data {
66     isl_val		    *n;
67     double		    RE_sum[nr_methods];
68 
69     my_clock_t		    ticks[nr_methods];
70     size_t		    size[nr_methods];
71 };
72 
result_data_init(isl_ctx * ctx,struct result_data * result)73 void result_data_init(isl_ctx *ctx, struct result_data *result)
74 {
75     int i;
76     for (i = 0; i < nr_methods; ++i) {
77 	result->RE_sum[i] = 0;
78 	result->ticks[i] = 0;
79 	result->size[i] = 0;
80     }
81     result->n = isl_val_zero(ctx);
82 }
83 
result_data_clear(struct result_data * result)84 void result_data_clear(struct result_data *result)
85 {
86 	isl_val_free(result->n);
87 }
88 
result_data_print(struct result_data * result,int n)89 void result_data_print(struct result_data *result, int n)
90 {
91     int i;
92 
93     fprintf(stderr, "%g", (double)result->ticks[0]/n);
94     for (i = 1; i < nr_methods; ++i)
95 	fprintf(stderr, ", %g", (double)result->ticks[i]/n);
96     fprintf(stderr, "\n");
97 
98     fprintf(stderr, "%zd/%d", result->size[0], n);
99     for (i = 1; i < nr_methods; ++i)
100 	fprintf(stderr, ", %zd/%d", result->size[i], n);
101     fprintf(stderr, "\n");
102 
103     fprintf(stderr, "%g\n", isl_val_get_d(result->n));
104     fprintf(stderr, "%g", result->RE_sum[0]/isl_val_get_d(result->n));
105     for (i = 1; i < nr_methods; ++i)
106 	fprintf(stderr, ", %g", result->RE_sum[i]/isl_val_get_d(result->n));
107     fprintf(stderr, "\n");
108 }
109 
110 struct verify_point_bound {
111 	struct verify_point_data vpd;
112 	struct result_data *result;
113 	isl_pw_qpolynomial *pwqp;
114 	isl_pw_qpolynomial_fold **pwf;
115 };
116 
verify_point(__isl_take isl_point * pnt,void * user)117 static isl_stat verify_point(__isl_take isl_point *pnt, void *user)
118 {
119 	struct verify_point_bound *vpb = (struct verify_point_bound *) user;
120 	const struct verify_options *options = vpb->vpd.options;
121 	int i;
122 	unsigned nparam;
123 	isl_val *max, *min, *exact, *approx, *t;
124 	isl_pw_qpolynomial *pwqp;
125 	isl_printer *p;
126 
127 	vpb->vpd.n--;
128 
129 	pwqp = isl_pw_qpolynomial_copy(vpb->pwqp);
130 
131 	nparam = isl_pw_qpolynomial_dim(pwqp, isl_dim_param);
132 	for (i = 0; i < nparam; ++i) {
133 		t = isl_point_get_coordinate_val(pnt, isl_dim_param, i);
134 		pwqp = isl_pw_qpolynomial_fix_val(pwqp, isl_dim_param, i, t);
135 	}
136 
137 	max = isl_pw_qpolynomial_max(isl_pw_qpolynomial_copy(pwqp));
138 	max = isl_val_floor(max);
139 
140 	min = isl_pw_qpolynomial_min(pwqp);
141 	min = isl_val_ceil(min);
142 
143 	exact = isl_val_sub(isl_val_copy(max), isl_val_copy(min));
144 	exact = isl_val_add_ui(exact, 1);
145 
146 	if (options->print_all) {
147 		p = isl_printer_to_file(isl_point_get_ctx(pnt), stderr);
148 		p = isl_printer_print_str(p, "max: ");
149 		p = isl_printer_print_val(p, max);
150 		p = isl_printer_print_str(p, ", min: ");
151 		p = isl_printer_print_val(p, min);
152 		p = isl_printer_print_str(p, ", range: ");
153 		p = isl_printer_print_val(p, exact);
154 	}
155 	isl_val_free(max);
156 	isl_val_free(min);
157 
158 	for (i = 0; i < nr_methods; ++i) {
159 		double error;
160 
161 		max = isl_pw_qpolynomial_fold_eval(
162 				isl_pw_qpolynomial_fold_copy(vpb->pwf[2 * i]),
163 				isl_point_copy(pnt));
164 		max = isl_val_floor(max);
165 
166 		min = isl_pw_qpolynomial_fold_eval(
167 				isl_pw_qpolynomial_fold_copy(vpb->pwf[2 * i + 1]),
168 				isl_point_copy(pnt));
169 		min = isl_val_ceil(min);
170 
171 		approx = isl_val_sub(max, min);
172 		approx = isl_val_add_ui(approx, 1);
173 		if (options->print_all) {
174 			p = isl_printer_print_str(p, ", ");
175 			p = isl_printer_print_val(p, approx);
176 		}
177 
178 		assert(isl_val_ge(approx, exact));
179 		approx = isl_val_sub(approx, isl_val_copy(exact));
180 
181 		error = fabs(isl_val_get_d(approx)) / isl_val_get_d(exact);
182 		if (options->print_all)
183 			fprintf(stderr, " (%g)", error);
184 		vpb->result->RE_sum[i] += error;
185 
186 		isl_val_free(approx);
187 	}
188 
189 	if (options->print_all) {
190 		p = isl_printer_end_line(p);
191 		isl_printer_free(p);
192 	} else if ((vpb->vpd.n % vpb->vpd.s) == 0) {
193 		printf("o");
194 		fflush(stdout);
195 	}
196 
197 	isl_val_free(exact);
198 	isl_point_free(pnt);
199 
200 	return vpb->vpd.n >= 1 ? isl_stat_ok : isl_stat_error;
201 }
202 
test(__isl_keep isl_pw_qpolynomial * pwqp,__isl_keep isl_pw_qpolynomial_fold ** pwf,struct result_data * result,struct verify_options * options)203 static void test(__isl_keep isl_pw_qpolynomial *pwqp,
204 	__isl_keep isl_pw_qpolynomial_fold **pwf, struct result_data *result,
205 	struct verify_options *options)
206 {
207 	struct verify_point_bound vpb = { { options }, result };
208 	isl_set *context;
209 	int r;
210 
211 	vpb.pwf = pwf;
212 	vpb.pwqp = pwqp;
213 	context = isl_pw_qpolynomial_domain(isl_pw_qpolynomial_copy(vpb.pwqp));
214 	context = isl_set_params(context);
215 	context = verify_context_set_bounds(context, options);
216 
217 	r = verify_point_data_init(&vpb.vpd, context);
218 	assert(r == 0);
219 	result->n = isl_val_set_si(result->n, vpb.vpd.n);
220 
221 	isl_set_foreach_point(context, verify_point, &vpb);
222 	assert(!vpb.vpd.error);
223 
224 	isl_set_free(context);
225 
226 	verify_point_data_fini(&vpb.vpd);
227 }
228 
handle(isl_ctx * ctx,FILE * in,struct result_data * result,struct verify_options * options)229 static void handle(isl_ctx *ctx, FILE *in, struct result_data *result,
230 	struct verify_options *options)
231 {
232     int i;
233     isl_pw_qpolynomial *pwqp, *upper, *lower;
234     isl_pw_qpolynomial_fold *pwf[2*nr_methods];
235     struct isl_stream *s;
236 
237     s = isl_stream_new_file(ctx, in);
238     pwqp = isl_stream_read_pw_qpolynomial(s);
239 
240     upper = isl_pw_qpolynomial_to_polynomial(isl_pw_qpolynomial_copy(pwqp), 1);
241     lower = isl_pw_qpolynomial_to_polynomial(isl_pw_qpolynomial_copy(pwqp), -1);
242 
243     for (i = 0; i < nr_methods; ++i) {
244 	int j;
245 	struct tms st_cpu;
246 	struct tms en_cpu;
247 
248 	times(&st_cpu);
249 	for (j = 0; j < 2; ++j) {
250 	    isl_pw_qpolynomial *poly = j == 0 ? upper : lower;
251 	    enum isl_fold type = j == 0 ? isl_fold_max : isl_fold_min;
252 	    isl_options_set_bound(ctx, methods[i].method);
253 	    poly = isl_pw_qpolynomial_copy(poly);
254 	    pwf[2*i+j] = isl_pw_qpolynomial_bound(poly, type, NULL);
255 	}
256 	times(&en_cpu);
257 	result->ticks[i] = time_diff(&en_cpu, &st_cpu);
258 	result->size[i] = isl_pw_qpolynomial_fold_size(pwf[2*i]);
259 	result->size[i] = isl_pw_qpolynomial_fold_size(pwf[2*i+1]);
260 	if (options->barvinok->verbose) {
261 	    isl_printer *p = isl_printer_to_file(ctx, stdout);
262 	    for (j = 0; j < 2; ++j) {
263 		p = isl_printer_print_pw_qpolynomial_fold(p, pwf[2*i+j]);
264 		p = isl_printer_end_line(p);
265 	    }
266 	    isl_printer_free(p);
267 	}
268     }
269     test(pwqp, pwf, result, options);
270 
271     for (i = 0; i < 2*nr_methods; ++i)
272 	isl_pw_qpolynomial_fold_free(pwf[i]);
273     isl_pw_qpolynomial_free(pwqp);
274     isl_pw_qpolynomial_free(lower);
275     isl_pw_qpolynomial_free(upper);
276 
277     isl_stream_free(s);
278 }
279 
main(int argc,char ** argv)280 int main(int argc, char **argv)
281 {
282     char path[PATH_MAX+1];
283     struct result_data all_result;
284     int n = 0;
285     struct options *options = options_new_with_defaults();
286     isl_ctx *ctx = isl_ctx_alloc();
287 
288     argc = options_parse(options, argc, argv, ISL_ARG_ALL);
289 
290     if (options->verify->M == INT_MIN)
291 	options->verify->M = 100;
292     if (options->verify->m == INT_MAX)
293 	options->verify->m = -100;
294 
295     result_data_init(ctx, &all_result);
296 
297     while (fgets(path, sizeof(path), stdin)) {
298 	struct result_data result;
299 	FILE *in;
300 	int i;
301 
302 	++n;
303 	result_data_init(ctx, &result);
304 	fprintf(stderr, "%s", path);
305 	*strchr(path, '\n') = '\0';
306 	in = fopen(path, "r");
307 	assert(in);
308 	handle(ctx, in, &result, options->verify);
309 	fclose(in);
310 
311 	if (!options->quiet)
312 	    result_data_print(&result, 1);
313 
314 	all_result.n = isl_val_add(all_result.n, isl_val_copy(result.n));
315 	for (i = 0; i < nr_methods; ++i) {
316 	    all_result.RE_sum[i] += result.RE_sum[i];
317 	    all_result.ticks[i] += result.ticks[i];
318 	    all_result.size[i] += result.size[i];
319 	}
320 
321 	result_data_clear(&result);
322 
323 	if (!options->quiet) {
324 	    fprintf(stderr, "average\n");
325 	    result_data_print(&all_result, n);
326 	}
327     }
328 
329     result_data_clear(&all_result);
330 
331     options_free(options);
332     isl_ctx_free(ctx);
333 
334     return 0;
335 }
336